{'fraction_non_alphanumeric': 0.03366719762486552, 'fraction_numerical': 0.024126629360598454, 'mean_word_length': 4.859844597566339, 'pattern_counts': {'":': 0, '<': 12, '<?xml version=': 0, '>': 6, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Basic Understanding of Condensed Phases of Matter via Packing Models

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 d-dimensional Euclidean space R d . Exclusionvolume 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. [1][2][3][4][5][6][7][8][9][10][11][12][13] Packings also serve as excellent models of the structure of multiphase heterogeneous materials, [14][15][16][17] colloids, 12,18,19 suspensions, 20,21 and granular media, 22,23 which enables predictions of their effective transport, mechanical and electromagnetic properties. 15 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, 24 "crowding" of macromolecules within living cells, 25 the packing of cells to form tissue, 15,[26][27][28] the fascinating spiral patterns seen in plant shoots and flowers (phyllotaxis), 29,30 and the competitive settlement of territories. 15,31 Understanding the symmetries and other mathematical properties of the densest sphere packings in various spaces and dimensions is a challenging area of longstanding interest in discrete geometry and number theory 32-34 as well as coding theory. 32,35,36 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 37 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. 38 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. 32 and 34 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, 39 which is regarded as one of the pillars of statistical mechanics. [40][41][42] 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 22,27,39,[43][44][45][46] and naturally arise in pure mathematical contexts. 32,47,48 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, 39,43 which can be regarded to be prototypical structural glasses, as well as the lowestdensity jammed packings. 49 Importantly, the simplified hard-sphere model embodies the primary structural attributes of dense manyparticle 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 systems 50 and highpressure crystal phases. Indeed, the equilibrium hardsphere model 9 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. 6,51 Moreover, the classic hard-sphere model provides a good description of certain classes of colloidal systems. 18,[52][53][54][55] 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 combustion 56 flow in packed beds, 57 and sintering of powders 58 ), 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. 59,60 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. 32,35 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 P is a collection of nonoverlapping solid objects or particles of general shape in d-dimensional Euclidean space R d . Packings can be defined in other spaces (e.g., hyperbolic spaces and compact spaces, such as the surface of a d-dimensional sphere), but our primary focus in this review is R d . 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 R d covered by the particles. A d-dimensional particle is centrally symmetric if it has a center C that bisects every chord through C connecting any two boundary points of the particle, i.e., the center is a point of inversion symmetry. Examples of centrally symmetric particles in R d 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 d-dimensional centrally symmetric particle for d ≥ 2 is said to possess d equivalent principal axes (orthogonal directions) associated with the moment of inertia tensor if those directions are two-fold rotational symmetry axes. Whereas a d-dimensional superball has d equivalent directions, a d-dimensional ellipsoid generally does not. The reader is referred to Ref. 39 for further discussion concerning particle symmetries.
A lattice Λ in R d is a subgroup consisting of the integer linear combinations of vectors that constitute a basis for R d . 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. 32 A lattice packing P L 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 R d . In a lattice packing, the space R d can be geometrically divided into identical regions F called fundamental cells, each of which has volume Vol(F ) and contains the centroid of just one particle of volume v 1 . Thus, the lattice packing fraction is are equivalent, i.e., D 4 ≡ D * 4 . The hypercubic lattice Z d ≡ Z d * and its dual lattice are equivalent for all d. A periodic packing of congruent particles is obtained by placing a fixed configuration of N particles (where N ≥ 1) 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 N 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 ρ = N/Vol(F ) is the number density, i.e., the number of particles per unit volume. For simplicity, consider a packing of N identical ddimensional spheres of diameter D centered at the positions r N ≡ {r 1 , r 2 , . . . , r N } in a region of volume V in d-dimensional Euclidean space R d . Ultimately, we will pass to the thermodynamic limit, i.e., N → ∞, V → ∞ such that the number density ρ = N/V is a fixed positive constant and its corresponding packing fraction is given by where ρ is the number density, is the volume of a d-dimensional sphere of radius R, and Γ(x) is the gamma function. For an individual sphere, the kissing or contact number Z is the number of spheres that may simultaneously touch this sphere. In a sphere packing, the mean kissing or contact number per particle, Z, is the average of Z over all spheres. For a lattice packings, Z = Z. For statistically homogeneous sphere packings in R d , the quantity ρ n g n (r n ) is proportional to the probability density associated with simultaneously finding n sphere centers at locations r n ≡ {r 1 , r 2 , . . . , r n } in R d ; see Ref. 9 and references therein. With this convention, each n-particle correlation function g n approaches unity when all particles become widely separated from one another for any system without longrange order. Statistical homogeneity implies that g n 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., g n = g n (r 12 , r 13 , . . . , r 1n ), where r ij = r j − r i . The pair correlation function g 2 (r) is of basic interest in this review. If the system is also rotationally invariant (statistically isotropic), then g 2 depends on the radial distance r ≡ |r| only, i.e., g 2 (r) = g 2 (r). The total correlation function is defined by h(r) ≡ g 2 (r) − 1. Importantly, we focus in this review on disordered packings in which h(r) decays to zero for large |r| sufficiently rapidly so that its volume integral over all space exists. 61 As usual, we define the nonnegative structure factor S(k) for a statistically homogeneous packing as whereh(k) is the Fourier transform of h(r) and k is the wave vector. The nonnegativity of S(k) for all k follows physically from the fact that it is proportional to the intensity of the scattering of incident radiation on a manyparticle system. 9 The structure factor S(k) provides a measure of the density fluctuations in particle configurations at a particular wave vector k. For any point configuration in which the minimal pair distance is some positive number, such as a sphere packing, g 2 (r = 0) = 0 and h(r = 0) = −1, the structure factor obeys the following sum rule: 62 For a single-component system in thermal equilibrium at number density ρ and absolute temperature T , the structure factor at the origin is directly related to the isothermal compressibility κ T via the relation 9 where k B is the Boltzmann constant. It is well known from Fourier transform theory that if a real-space radial function f (r) in R d decreases sufficiently rapidly to zero for large r such that all of its even moments exist, then its Fourier transformf (k), where k ≡ |k| is a wavenumber, is an even and analytic function at k = 0. Hence, S(k), defined by (5), has an expansion about k = 0 in any space dimension d of the general form where S 0 and S 2 are the d-dependent constants defined by and where is the surface area of d-dimensional sphere of radius R. This behavior is to be contrasted with that of maximally random jammed sphere packings, 43 which possess a structure factor that is nonanalytic at k = 0 (Ref. 63), as discussed in greater detail in Sec. III E. A hyperuniform 64 point configuration in R d is one in which the structure factor S(k) tends to zero as the wavenumber tends to zero, i.e., lim |k|→0 S(k) = 0, (12) implying that single scattering of incident radiation at infinite wavelengths is completely suppressed. Equivalently, a hyperuniform system is one in which the number variance σ 2 N (R) ≡ N (R) 2 − N (R) 2 of particles within a spherical observation window of radius R grows more slowly than the window volume, i.e., R d , in the large-R limit. This class of point configurations includes perfect crystals, many perfect quasicrystals and special disordered many-particle systems. [64][65][66] Note that structurefactor definition (5) and the hyperuniformity requirement (12) dictate that the following sum rule involving h(r) that a hyperuniform point process must obey: This sum rules implies that h(r) must exhibit negative correlations, i.e., anticorrelations, for some values of r.
The hyperuniformity concept was generalized to incorporate two-phase heterogeneous media (e.g., composites, porous media and dispersions). 65,67 Here the phase volume fraction fluctuates within a spherical window of radius R and hence can be characterized by the volumefraction variance σ 2 V (R). For typical disordered twophase media, the variance σ 2 V (R) for large R goes to zero like R −d . However, for hyperuniform disordered twophase media, σ 2 V (R) goes to zero asymptotically more rapidly than the inverse of the window volume, i.e., faster than R −d , which is equivalent to the following condition on the spectral densityχ V (k): 67 The spectral density is proportional to scattered intensity associated with "mass" (volume) content of the phases. 68

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; 69 see also Refs. 9, 70-72. Here we focus on packings of frictionless, congruent spheres of diameter D 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 A) or derived from one of an infinite number of nonequilibrium ensembles 15 (see Sec. III B for examples).

A. 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 p of a stable thermodynamic phase in R d at packing fraction φ and temperature T is simply related to the contact value of the pair correlation function, g 2 (D + ): 15 Away from jammed states, it has been proved that the mean nearest-neighbor distance between spheres, λ, is bounded from above by the pressure, 73 namely, Figure 1 schematically shows the 3D phase behavior in the φ-p 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 (φ ≈ 0.49). Increasing the density beyond the freezing point putatively results in an entropy-driven first-order phase transition 70,74 to a crystal branch that begins at the melting point (φ ≈ 0.55). 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. 70 Slow compression of the system along the crystal branch must end at one of the optimal (maximally dense) sphere packings with φ = π/ √ 18 = 0.74048 . . ., each of which is a jammed packing in which each particle contacts 12 others (see Sec. III D). This equilibrium state has an infinite pressure and is putatively entropically favored by the fcc packing. 75 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, 13,15,76 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) state 43 with φ ≈ 0.64. Accurate approximate formulas for the pressure along such metastable extensions up to the jamming points have been obtained. 15,77 This ideal amorphous state is described in greater detail in Sec. III E. 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. 78 For d ≥ 2, approximate formulas for g 2 (r) are known along liquid branches. 9 Approximate closures of the Ornstein-Zernike integral equation linking the direct correlation function c(r) to the total correlation function h(r), 79 such as the Percus-Yevick (PY) and hypernetted chain schemes, 9,80,81 provide reasonably accurate estimates of g 2 (r) for hard-spheres liquids. Because g 2 (r) decays to unity exponentially fast for liquid states, we can conclude that it must have a corresponding structure factor S(k) that is an even function and analytic at k = 0; see Eq. (8). Also, since the leading- 15. Three different isothermal densification paths by which a hardsphere 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 firstorder 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.
order term S 0 must be positive because the isothermal compressibility is positive [cf. (7)], hard-sphere liquids are not hyperuniform. Figure 2 shows g 2 (r) and the corresponding structure factor S(k) in three dimensions at φ = 0.35 as obtained from the PY approximation.

B. Nonequilibrium Disordered Packings
Here we briefly discuss three different nonequilibrium sphere packings: random sequential addition, "ghost" random sequential addition, and random close packing.

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 R d ; see Refs. 82-90. The RSA packing process in the first three space dimensions has been used to model a variety of different condensed phases, including protein adsorption, 85 polymer oxidation, 91 particles in cell membranes 84 and ion implantation in semiconductors. 92 In its simplest rendition, a RSA sphere packing is produced by randomly, irreversibly, and sequentially placing nonoverlapping objects into a large volume in R d that at some initial time is empty of spheres. If an attempt to add a sphere at some time t 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 t, obtaining an RSA configuration with a time-dependent packing fraction φ(t) but this value cannot exceed the maximal saturation packing fraction φ s = φ(∞) 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 C; moreover, the pair correlation function g 2 (r) decays to unity for large r superexponentially fast. 93 The latter property implies that the structure factor S(k) of RSA packings in R d must be an even and analytic function at k = 0, 89 as indicated in Eq. (8). It is notable that saturated RSA packings are not hyperuniform for any finite space dimension. 89,90 In the one-dimensional case, also known as the "carparking" problem, the saturation packing fraction was obtained analytically: φ s = 0.7475979202.... 82 However, for d ≥ 2, φ s in the thermodynamic limit has only been estimated via numerical simulations. The most precise numerical study to date 90 has yielded φ s = 0.5470735 ± 0.0000028 for d = 2 and φ s = 0.3841307 ± 0.0000021 for d = 3. Estimates of φ s in higher dimensions are discussed in Sec. IV. RSA packings have also been examined for spheres with a size distribution 94

"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 [0, 1]; see Ref. 104. When κ = 0, one recovers the standard RSA process and when κ = 1, 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 n-particle correlation functions g n for any n and d. The reader is referred to Ref. 104 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 t for spheres of unit diameter is given by is the volume of sphere of radius R. Thus, we see that as t → +∞, φ = 2 −d , which is appreciably smaller than the RSA saturation packing fraction φ s 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 d and all realizable densities, which has implications for high-dimensional packings, as discussed in Sec. IV.

Random Close Packing
The "random close packing" (RCP) notion was pioneered by Bernal 2,3 to model the structure of liquids, and has been one of the more persistent themes with a venerable history. [105][106][107][108][109][110][111][112] 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 welldefined 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) . . . 0.636"; see Ref. 106.
Torquato, Truskett and Debenedetti 43 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) 113 moleculardynamics growth algorithm to generate jammed packings, it was shown 43 that fastest particle growth rates generated the most disordered sphere packings with φ ≈ 0.64, but that by slowing the growth rates larger packing fractions could be continuously achieved up to the densest value φ max = 0.74048 . . . 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 43 replaced the notion of "close packing" with "jamming" categories (defined precisely in Sec. III C) 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 C 4 for details). This work pointed the way toward a quantitative means of characterizing all packings, namely, the geometric-structure approach.
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).

C. 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, 2,3,22,44,45,[115][116][117][118][119][120][121][122] 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, 43,49,114,[123][124][125][126][127] 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.

Jamming Categories
Three broad and mathematically precise "jamming" categories of sphere packings can be distinguished depending on the nature of their mechanical stability. 64,123 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 d + 1 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 volumenonincreasing 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, 128 and imply that there can be no "rattlers" (i.e., movable but caged particles) in the packing.

Polytope Picture and Pressure in Jamming Limit
A packing of N hard spheres of diameter D in a jammed framework in R d is specified by a N ddimensional configurational position vector R = r N ≡ {r 1 , . . . , r N }. Isostatic jammed packings possess the minimal number of contacts for a jamming category and boundary conditions. 133 The relative differences between isostatic collective and strict jammed packings diminish as N becomes large, and since the number of degrees of freedom is essentially equal to N d, an isostatic packing has a mean contact number per particle, Z, equal to 2d in the large-N 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. 127 Consider decreasing the density slightly in a sphere packing that is at least collectively jammed by reducing the particle diameter by ∆D, so that the packing fraction 1 is called the jamming gap or distance to jamming. There is a sufficiently small δ that does not destroy the jamming confinement property. For fixed N 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) P. 5,72 This polytope P is determined from the linearized impenetrability equations 114,129 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 p for the trapped but locally ergodic system by considering time averages as the system executes a tightly confined motion around the particular jammed configuration R J . The probability distribution P f (f ) of the time-averaged interparticle forces f has been rigorously linked to the contact value r = D of the pair correlation function g 2 (r). 133 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: 5,72,133 where T 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 p = +∞. 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 A). 15, 77 Kamien and Liu 134 assumed the same free-volume form to fit data for the pressure of "metastable" hard-sphere states.

Hard-Particle Jamming Algorithms
For many years, the Lubachevsky-Stillinger (LS) algorithm 113 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. 43,124,135 This is an event-driven or collisiondriven 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 136,137 to generate jammed packings of smoothly shaped nonspherical particles, including ellipsoids, 127,137 superdisks, 138 , and superballs. 139 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 140,141 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 142 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, 142-145 disordered hyperstatic packings 146 and the globally maximally dense inherent structures, 142,147-150 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. 143 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).

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 43,64,124,135 but also simple liquids, 13,135,151,152 glasses, 135,153 water, 154,155 disordered ground states, 156 random media, 157 and biological systems. 27 It is quite reasonable to consider "entropic" measures to characterize the randomness of packings. However, as pointed out by Kansal et al. 124 , 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 largesystem 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 E, including its significance for certain 2D hard-disk packings.
We know that a many-body system of N particles is completely characterized statistically by its N -body probability density function P (R; t) that is associated with finding the N -particle system with configuration R at some time t. Such complete information is virtually never available for large N 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 R; (2) it is subject typically to the normalization 0 ≤ ψ ≤ 1; and, (3) for any two configurations R A and R B , ψ(R A ) > ψ(R B ) implies that configuration R A is to be considered as more ordered than configuration R B . 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 158 and three dimensions, 159 translational order metrics, 43,135,156 and hyperuniformity order metrics. 64,65 These specific or-der 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: 39,124 (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. 15,43 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, 13,135,151,152 models of water, 154,155 disordered ground states of long-ranged isotropic pair potentials, 156 and models of amorphous polymers. 160

D. 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, Z, 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 49 (least dense jammed packings) as well as states between these and MRJ.
The simplest order map is the one that maps any hardsphere 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 43,49,114,124 from non-jammed packings in three dimensions. The boundaries of the jammed region delineate extremal structures (see Fig. 7). The densest sphere packings, 38 which lie along the locus B-B with φ max = π/ √ 18 ≈ 0.74, are strictly jammed. 114,123 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. 125 The most disordered subset of the stacking variants of the fcc packing is denoted by point B . In two dimensions, the strictly jammed triangular lattice is the unique densest packing 161 and so the line B-B in R 3 collapses to a single point B in R 2 .
The least dense strictly jammed packings are conjectured to be the "tunneled crystals" in R 3 with φ min = 2φ max /3 = 0.49365 . . ., corresponding to the locus of points A-A . 49 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 im-portance 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. 49 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 Sondhi 162 showed that their underlying topologies greatly simplify the determination of their antiferromagnetic properties. In R 2 , the "reinforced" kagomé packing with exactly four contacts per particle (isostatic value) is evidently the lowest density strictly jammed packing 114 with φ min = √ 3π/8 = 0.68017 . . . and so the line A-A in R 3 collapses to a single point A in R 2 .

E. Maximally Random Jammed States
Among all strictly jammed sphere packings in R d , the one that exhibits maximal disorder (minimizes some given order metric ψ) is of special interest. This is called the maximally random jammed (MRJ) state 43 ; 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 φ MRJ ≈ 0.64 (see Ref. 124), close to the traditionally advocated density of the RCP state, and with an isostatic mean contact number Z = 6 (see Ref. 133). 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, 124 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 2,3 , which, more recently, have been suggested to be the most probable jammed configurations within an ensemble. 44 The differences between these states are even starker in two dimensions, e.g., MRJ packings of identical circular disks in R 2 have been shown to be dramatically different from their RCP counterparts, including their respective densities, average contact numbers, and degree of order, 144 as detailed below.
The MRJ state under the strict-jamming constraint is a prototypical glass 49 in that it is maximally disordered (according to a variety of order metrics) without any long-range order (Bragg peaks) and perfectly rigid The pair correlation function g 2 (r) versus r/D − 1. 133 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 S(k) as a function of the dimensionless wavenumber kD/(2π) for a million-particle packing. 63 The inset shows that S(k) tends to zero (within numerical error) linearly in k hence an MRJ packing is effectively hyperuniform.
(i.e., the elastic moduli are unbounded. 39,125 ) Its pair correlation function can be decomposed into a Diracdelta-function contribution from particle contacts and a continuous-function contribution g cont 2 (r): where s 1 (R) is given by (11), δ(r) is a radial Dirac delta function and Z = 2d. The corresponding structure factor in the long-wavenumber limit is (18) For d = 3, g cont 2 (r) possesses the well-known feature of a split second peak, 163 with a prominent discontinuity at twice the sphere diameter, as shown in Fig. 8. Interestingly, an integrable power-law divergence (1/r α with α ≈ 0.4) exists for near contacts. 63,133 Moreover, an MRJ packing in R d has a structure factor S(k) that tends to zero linearly in |k| (within numerical error) in the limit |k| → 0 and hence is hyperuniform (see Fig. 8), belonging to the same class as Fermi-sphere point processes 164 and superfluid helium in its ground state. 165,166 This nonanalytic behavior at |k| = 0 implies that MRJ packings are characterized by negative "quasi-long-range" (QLR) pair correlations in which the total correlation function h(r) decays to zero asymptotically like −1/|r| d+1 ; see Refs. 62,63,145,[167][168][169]. The QLR hyperuniform behavior distinguishes the MRJ state strongly from that of the equilibrium hard-sphere fluid 170 , which possesses a structure factor that is analytic at k = 0 [cf. (8)], and thus its h(r) decays to zero exponentially fast for large r; see Ref. 39. 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. 168,169 This led to the identification of static nonequilibrium growing length scales in supercooled liquids on approaching the glass transition. 171 The following is a general conjecture due to Torquato and Stillinger 64 concerning the conditions under which certain strictly jammed packings are hyperuniform: Conjecture: Any strictly jammed saturated infinite packing of identical spheres in R d is hyperuniform.
To date, there is no rigorously known counterexample to this conjecture. Its justification and rigorously known examples are discussed in Refs. 62 and 145.
A recent numerical study of necessarily finite disordered packings calls into question the connection between jamming and perfect hyperuniformity. 172 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. [172][173][174] Current numerically generated disordered packings that are putatively jammed tend to contain a small concentration of rattlers; 44,63,142,143,175 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. 145 Although these packings appear to be jammed by conventional tests, rigorous linearprogramming jamming tests 114,129 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" 145 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 S(0) exactly down to its minimum value of zero if a true jammed critical state could be attained; typically, 63 S(0) ≈ 10 −4 .
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 threepart series, Klatt and Torquato 176 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. 176 In the second paper of this series, 177 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, 178 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. 179. In a separate study, Ziff and Torquato 180 determined the site and bond percolation thresholds of TJ generated MRJ packings to be p c = 0.3116(3) and p c = 0.2424(3), 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 A-A 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 0.6 < φ < 0.74048 . . . 43,124,146 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 142 was applied to yield disordered strictly jammed packings 146 with φ as low as 0.60, which are overconstrained with Z ≈ 6.4, and hence are more ordered than the MRJ state.
It is well known that lack of "frustration" 15,111 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 φ ≈ 0.88 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 pro-tocols, "entropic measures" of disorder would identify these as the most disordered, a misleading conclusion. Recently, Atkinson et al. 144 applied the TJ algorithm to generate relatively large MRJ disk packings with exactly isostatic (Z = 4) jammed backbones and a packing fraction (including rattlers) of φ MRJ ≈ 0.827. 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, 181 yielding the prediction φ MRJ = 0.834 in the monodisperse-disk limit, which is in very good agreement with the aforementioned recent numerical estimate. 144

F. 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). 182 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, 56 flow in packed beds, 57 sintering of powders, 58 colloids, 18 transport and mechanical properties of particulate composite materials. 15 Packings of biological cells in tissue are better modeled by assuming that the spheres have a size distribution. 27,28 Generally, we allow the spheres to possess a continuous distribution in radius R, which is characterized by a probability density function f (R). The average of any function w(R) is defined by w(R) = ∞ 0 w(R)f (R)dR. The overall packing fraction φ is then defined as where ρ is the total number density, v 1 (R) is given by (4). Two continuous probability densities that have been widely used are the Schulz 183 and log-normal 184 distributions. One can obtain corresponding results for spheres with M discrete different sizes from Eq. (19) by letting where ρ i and R i are number density and radius of type-i particles, respectively, and ρ is the total number density. Therefore, the overall volume fraction using (19) is given is the packing fraction of the ith component.

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 fluidfluid 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, 54,55,[185][186][187][188][189][190][191][192][193][194][195][196][197] which we only briefly touch upon here for both additive and nonadditive cases. In additive hardsphere 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 M components within the Percus-Yevick approximation. 185 Accurate approximate equations of state under liquid conditions have been found for both discrete 188,195 and continuous 189,198 size distributions with additive hard cores. Fundamental-measure theory can provide useful insights about the phase behavior of hard-sphere mixtures. 191,197 This theory predicts the existence of stable fluid-fluid coexistence for sufficiently large size ratios in additive binary hard-sphere systems, while numerical simulations 54 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 F 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. 186 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. 192,193 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 54,194 that belongs to the Ising universality class, 196 while those with negative nonadditivity have tendencies to form alloyed (mixed) fluid phases. 190 For certain binary mixtures of hard sphere fluids with nonadditive diameters, the pair correlation function has been determined in the Percus-Yevick approximation. 187,190 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" 27 systems at positive temperatures. 199,200

Maximally Random Jammed States
The study of dense disordered packings of 3D polydisperse additive spheres has received considerable attention. [201][202][203][204][205][206] 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. 207 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 x can be obtained whose packing fractions exceeds 0.64 and indeed can attain φ = 0.79 for α = 0.1 (the smallest value considered). Chaudhuri et al. 208 numerically generated amorphous 50-50 binary packings with packing fractions in the range 0.648 ≤ φ ≤ 0.662 for α = 1.4. It is notable Clusel et al. 209 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 202,206,210 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. 39,129 Using an the TJ linear programming packing algorithm, 142 Hopkins et al. 211 were recently able to generate high-fidelity strictly jammed, isostatic, disordered binary packings with an anomalously large range of packing fractions, 0.634 ≤ φ ≤ 0.829, with the size ratio restriction α ≥ 0.1. 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 x approach those of the corresponding densest known ordered packings. 147,148 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. 212 The LS algorithm has been used successfully to generate disordered strictly jammed packings of binary disks 213 with φ ≈ 0.84 and α −1 = 1.4. 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 phaseseparated crystal, it has been argued 213,214 that there is no "ideal glass transition". 215 The existence of an ideal glass transition remains a hotly debated topic of research.
Simulational 216-218 as well as experimental 219,220 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χ V (k), 216,217 which is defined by condition (14).

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 compounds 59 and solid rare-gas compounds 60 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 M different radii R 1 , R 2 , . . . , R M in R d . Specifically, the overall maximal packing fraction φ (M ) max of such a general mixture in R d [where φ is defined by (19) with (20)] is bounded from above and below in terms of the maximal packing fraction φ (1) max for a monodisperse packing in the infinite-volume limit: 15 The lower bound corresponds to the case when the M components completely demix, each at the density φ max . The upper bound corresponds to an ideal sequential packing process for arbitrary M in which one takes the limits R 1 /R 2 → 0, R 2 /R 3 → 0, · · · , R M −1 /R M → 0. 15 Specific nonsequential protocols (algorithmic or otherwise) that can generate structures that approach the upper bound (21) for arbitrary values of M are currently unknown and thus the development of such protocols is an open area of research. We see that in the limit M → ∞, the upper bound approaches unity, corresponding to space-filling polydisperse spheres with an infinitely wide separation in sizes 221 or a continuous size distribution with sizes ranging to the infinitesimally small. 15 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 φ max as a function of α and x. In practice, we have a sketchy understanding of the surface defined by φ max (α, x).
Among the 2D and 3D cases, we know most about the determination of the maximally dense binary packings in R 2 . Fejes Tóth 161 reported a number of candidate maximally dense packings for certain values of the radii ratio in the range α ≥ 0.154701 . . .. Maximally dense binary disk packings have been also studied to determine the stable crystal phase diagram of such alloys. 222 The determination of φ max 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. 223 All of these results can be compared to a relatively tight upper bound on φ max given by 224 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 R 3 in the physical sciences because of their relationship to binary crystal alloys found in molecular systems; see Refs. 148 and 225 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. 226,227 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. 147,148 have presented the most comprehensive determination to date of the "phase diagram" for the densest binary sphere packings via the TJ linear-programming algorithm. 142 In Ref. 148, 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. 59,60 Reference 148 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. 33, 39, 45, 61, 104, 149, 150, 167, 228-243. Physicists have studied high-dimensional packings to gain insight into liquid, crystal and glassy states of matter in lower dimensions. 45,167,229,230,232,235,236,239,240 Finding the densest packings in arbitrary dimension in Euclidean and compact spaces is a problem of long-standing interest in discrete geometry. 32,33,244,245 Remarkably, the optimal way of sending digital signals over noisy channels corresponds to the densest sphere packing in a highdimensional space. 32,35 These "error-correcting" codes underlie a variety of systems in digital communications and storage, including compact disks, cell phones and the Internet.

A. 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. 167,229,[239][240][241]246,247 This includes the evaluation of various virial coefficients across dimensions, 246,247 as well as the pressure along the liquid, metastable and crystal branches, 167,235,239-241 especially for d = 4, 5, 6 and 7. Using the LS packing algorithm, Skoge et al. 167 numerically estimated the freezing and melting packing fractions for the equilibrium hard-sphere fluid-solid transition, φ F 0.32 and φ M 0.39, respectively, for d = 4, and φ F 0.19 and φ M 0.24, respectively, for d = 5. 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.. 240 Finken, Schmidt and Löwen 248 used a variety of approximate theoretical methods to show that equilibrium hard spheres have a first-order freezing transition for dimensions as high as d = 50.
Any disordered packing in which the pair correlation function at contact, g 2 (D + ), is bounded (such as equilibrium hard spheres) induces a power-law decay in the structure factor S(k) in the limit k → ∞ for any dimension d given by There is a remarkable duality between the equilibrium hard-hypersphere (hypercube) fluid system in R d and the continuum percolation model of overlapping hyperspheres (hypercubes) in R d . 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. 249 This mapping becomes exact for d = 1 and in the large-d limit.

B. Nonequilibrium Disordered Packings Via Sequential Addition
In Sec. III B 2, we noted that the "ghost" RSA packing 104 is a disordered but unsaturated packing construction whose n-particle correlation functions are known exactly and rigorously achieves the infinite-time packing fraction φ = 2 −d for any d; 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 d = 6 (Ref. 89). A more efficient numerical procedure was devised to produce such packings for dimensions up through d = 8 (Ref. 90). The current best estimates of the maximal saturation packing fraction φ s for d = 4, 5, 6, 7 and 8 are 0.2600781±0.0000037, 0.1707761 ± 0.0000046, 0.109302 ± 0.000019, 0.068404 ± 0.000016, and 0.04230±0.00021, respectively. 90 These are lower than the corresponding MRJ packing fractions in those dimensions (see Sec. IV C). The quantity φ s apparently scales as d · 2 −d or possibly d · ln(d) · 2 −d for large d; see Refs. 89 and 90. While saturated RSA packings are nearly but not exactly hyperuniform, 89 as d increases, the degree of hyperuniformity increases and pair correlations markedly diminish, 90 consistent with the decorrelation principle, 61 which is described in Sec. IV E.

C. Maximally Random Jammed States
Using the LS algorithm, Skoge et al. 167 generated and characterized MRJ packings in four, five and six dimensions. In particular, they estimated the MRJ packing fractions, finding φ MRJ 0.46, 0.31 and 0.20 for d = 4, 5 and 6, respectively. To a good approximation, the MRJ packing fraction obeys the scaling form φ MRJ = c 1 2 −d + c 2 d · 2 −d , where c 1 = −2.72 and c 2 = 2.56, which appears to be consistent with highdimensional asymptotic limit, albeit with different coefficients. The dominant large-d density scaling d · 2 −d is supported by theoretical studies. 61,242,250,251 Skoge et al. 167 also determined the MRJ pair correlation function g 2 (r) and structure factor S(k) for these states and found that short-range ordering appreciably decreases with increasing dimension, consistent with the decorrelation principle. 61 This implies that, in the limit d → ∞, g 2 (r) tends to unity for all r outside the hard-core, except for a Dirac-delta function at contact due to the jamming constraint. 61 As for d = 3 (where φ MRJ 0.64), the MRJ packings were found to be isostatic, hyperuniform, and have a power-law divergence in g 2 (r) at contact; ; g 2 (r) ∼ 1/(r − D) α with α ≈ 0.4 as r tends to D + . Across dimensions, the cumulative number of neighbors was shown to equal the kissing (contact) number of the conjectured densest packing close to where g 2 (r) has its first minimum. Disordered jammed packings were also simulated and studied in dimensions 7-10; see Ref. 252.

D. Maximally Dense Sphere Packings
The sphere packing problem seeks to answer the following question: 32 Among all packings of congruent spheres in R d , what is the maximal packing fraction φ max and the corresponding arrangements of the spheres? Until 2017, the optimal solutions were known only for the first three space dimensions. 38 For d = 2 and d = 3, these are the triangular lattice (A 2 ) with φ max = π/ √ 12 = 0.906899 . . . and checkerboard (fcc) lattice (D 3 ) and its stacking variants with φ max = π/ √ 18 = 0.74048 . . ., respectively. We now know the optimal solutions in two other space dimensions; namely, the E 8 and Λ 24 lattices are the densest packings among all possible packings in R 8 and R 24 , respectively; see Refs. 244 and 245. For 4 ≤ d ≤ 9, the densest known packings are (Bravais) lattice packings. 32 The "checkerboard" lattice D d is believed to be optimal in R 4 and R 5 . Interestingly, the non-lattice (periodic) packing P 10c (with 40 spheres per fundamental cell) is the densest known packing in R 10 , which is the lowest dimension in which the best known packing is not a (Bravais) lattice. Table I lists the densest known or optimal sphere packings in R d for selected d. 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 64,65 as well as of two other well-known problems in discrete geometry: the covering and quantizer problems, 32,253 but such relationships may or may not exist for d ≥ 4, depending on the peculiarities of the dimensions involved. 254 The TJ linear-programming packing algorithm was adapted by Marcotte and Torquato 149 to determine the densest lattice packing (one particle per fundamental cell) in some high dimension. These authors applied it for 2 ≤ d ≤ 19 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. 255,256 The TJ algorithm was used to generate an ensemble of isostatic jammed hard-sphere lattices and study the associated pair statistic and force distributions. 150 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 d, lattice packings are most likely not the densest (see Fig. 9), but it becomes increasingly difficult to find explicit dense packing constructions as d increases. Indeed, the problem of finding the shortest lattice vector in a particular lattice packing grows super-exponentially with d and is in the class of NP-hard problems. 257 For large d, the best that one can do theoretically is Placing an additional sphere in Z 4 doubles the density of Z 4 and yields the four-dimensional checkerboard lattice packing D4, which is believed to be the optimal packing in R 4 . Right panel: A schematic "effective" distorted representation of the hypersphere within the hypercubic fundamental cell for large d, illustrating that the volume content of the hypersphere relative to the hypercube rapidly diminishes asymptotically. Indeed, the packing fraction of Z d is given by φ = π d/2 /(Γ(1 + d/2)2 d ). Indeed, the checkerboard lattice D d with packing fraction φ = π d/2 /(Γ(1 + d/2)2 (d+2)/2 ) becomes suboptimal in relatively low dimensions because it also becomes dominated by larger and larger holes as d increases.
to devise upper and lower bounds on φ max . 32 The nonconstructive lower bound of Minkowski 258 established the existence of reasonably dense lattice packings. He found that the maximal packing fraction φ L max among all lattice packings for d ≥ 2 satisfies where ζ(d) = ∞ k=1 k −d is the Riemann zeta function. Note that for large values of d, the asymptotic behavior of the Minkowski lower bound is controlled by 2 −d . Since 1905, many extensions and generalizations of the lower bound (24) have been derived 32,259-261 , but none of these investigations have been able to improve upon the dominant exponential term 2 −d ; they instead only improve on the latter by a factor linear in d.
It is trivial to prove that the packing fraction of a saturated packing of congruent spheres in R d satisfies 104 for all d. This "greedy" bound (25) has the same dominant exponential term as Minkowski's bound (24). Nontrivial upper bounds on φ max in R d for any d have been derived. 33,34,[262][263][264] The linear-programming (LP) upper bounds due to Cohn and Elkies 33 provides the basic framework for proving the best known upper bounds on φ max for dimensions in the range 4 ≤ d ≤ 36. Recently, these LP bounds have been used to prove that no packings in R 8 and R 24 have densities that can exceed those of the E 8 (Ref. 244) and L 24 (Ref. 245) lattices, respectively. Kabatiansky and Levenshtein 264 found the best asymptotic upper bound, which in the limit d → ∞ yields φ max ≤ 2 −0.5990d , proving that the maximal packing fraction tends to zero in the limit d → ∞. This rather counterintuitive high-dimensional property of sphere packings can be understood by recognizing that almost all of the volume of a d-dimensional sphere for large d is concentrated near the sphere surface.

E. Are Disordered Packings the Densest in High
Dimensions?
Since 1905, many extensions and generalizations of Minkowski's bound have been derived 32 , but none of them have improved upon the dominant exponential term 2 −d . The existence of the unjammed disordered ghost RSA packing 104 (Sec. III B 2) that rigorously achieves a packing fraction of 2 −d strongly suggests that Bravaislattice packings (which are almost surely unsaturated in sufficiently high d) are far from optimal for large d.
Torquato and Stillinger 61 employed a plausible conjecture that strongly supports the counterintuitive possibility that the densest sphere packings for sufficiently large d may be disordered or at least possess fundamental cells whose size and structural complexity increase with d. They did so using the so-called g 2 -invariant optimization procedure that maximizes φ associated with a radial "test" pair correlation function g 2 (r) to provide the putative exponential improvement on Minkowski's 100-year-old bound on φ max . Specifically, a g 2 -invariant process 250 is one in which the functional form of a "test" pair correlation g 2 (r) function remains invariant as density varies, for all r, over the range of packing fractions 0 ≤ φ ≤ φ * subject to the satisfaction of the nonnegativity of the structure factor S(k) and g 2 (r). When there exist sphere packings with a g 2 satisfying these conditions in the interval [0, φ * ], then one has the lower bound on the maximal packing fraction given by Torquato and Stillinger 61 conjectured that a test function g 2 (r) is a pair-correlation function of a translationally invariant disordered sphere packing in R d for 0 ≤ φ ≤ φ * for sufficiently large d if and only if the nonnegativity conditions on S(k) and g 2 (r). The decorrelation principle, 61 among other results, 39,237,243 provides justification for the conjecture. This principle states that unconstrained correlations in disordered sphere packings vanish asymptotically in high dimensions and that the g n for any n ≥ 3 can be inferred entirely (up to small errors) from a knowledge of ρ and g 2 . This is vividly exhibited by the exactly solvable ghost RSA packing process 104 as well as by computer simulations of high-dimensional MRJ 167 and RSA 234 packings. Interestingly, this optimization problem is the dual of the infinite-dimensional linear program (LP) devised by Cohn and Elkies 233 to obtain upper bounds on φ max .
Using a particular test pair correlation corresponding to a disordered sphere packing, Torquato and Stillinger 61 found a conjectural lower bound on φ max that is controlled by 2 −(0.77865...)d , thus providing the first putative exponential improvement on Minkowski's lower bound (24). Scardicchio, Stillinger and Torquato 237 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 2 −(0.77865...)d is much more general and robust than previously surmised.
Zachary and Torquato 265 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 Torquato 243 who showed that strictly jammed lattice sphere packings visibly decorrelate as d increases in relatively low dimensions.

F. Remarks on Packing Problems in Non-Euclidean Spaces
Particle packing problems in non-Euclidean (curved) spaces arises in a variety of fields, including physics, 266-268 biology, 29,269-272 communications theory, 32 and geometry. 32,36,48,273,274 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 S d−1 ⊂ R d . The reader is referred to the review by Torquato and Stillinger 39 for some discussion of negatively curved hyperbolic space H d .
Recall that the kissing number Z is the number of spheres of unit radius that simultaneously touch a unit sphere S d−1 . 32 The kissing number problem asks for the maximal kissing number Z max in R d . The determination of the maximal kissing number in R 3 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 Z max for d > 3 is only known in dimensions four, 275 eight, 276,277 and twenty four. 276,277 A packing of congruent spherical caps on the unit sphere S d−1 in R d yields a spherical code consisting of the centers of the caps. 32 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 36 and references therein for some developments in the mathematics literature. Interestingly, the jamming of spherical codes in a variety of dimensions has been investigated. 48 It is noteworthy that some optimal spherical codes 278 are related to the densest local packing of spheres around a central sphere. 279,280

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. 74,126,127,132,136,137,140,141,255, 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 19,140,304 are discussed in Sec.V D 5.  10. (color online). Superballs with different values of the deformation parameter p. We note that p = 0, 1/2, 1, and ∞ correspond to a 3D cross, regular octahedron, sphere and cube, respectively.

Ellipsoids
One simple generalization of the sphere is an ellipsoid, the family of which is a continuous deformation of a sphere. A d-dimensional ellipsoid in R d is a centrally symmetric body occupying the region where x i (i = 1, 2, . . . , d) are Cartesian coordinates and a i are the semi-axes of the ellipsoid. Thus, we see that an ellipsoid is an affine (linear) transformation of the sphere.

Superballs
A d-dimensional superball in R d is a centrally symmetric body occupying the region where x i (i = 1, . . . , d) are Cartesian coordinates and p ≥ 0 is the deformation parameter (not pressure, as denoted in Sec. III A), which controls the extent to which the particle shape has deformed from that of a d-dimensional sphere (p = 1). Thus, superballs constitute a large family of both convex (p ≥ 1/2) and concave (0 < p < 1/2) particles (see Fig. 10). A "superdisk",which is the designation in the 2D case, possesses square symmetry. As p moves away from unity, two families of superdisks with square symmetry can be obtained depending on whether p < 1 or p > 1. When p < 1/2, the superdisk is concave; see Ref. 138.

Spherocylinder
A d-dimensional spherocylinder in R d consists of a cylinder of length L and radius R capped at both ends by hemispheres of radius R , and therefore is a centrallysymmetric convex particle. Its volume V SC for d ≥ 2 is given by When L = 0, a d-dimensional spherocylinder reduces to a d-dimensional sphere of radius R. Figure 11 shows 3D examples. 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. 272.
An Archimedean solid is a highly symmetric, semiregular 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.

B. 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. 74,281,283,285,289,291,299,300,302,306,307,[310][311][312]314, and 317. Nonsphericalparticle 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, 281 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 hardnonspherical 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. 311,320 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. 317 Spherocylinders exhibit five different possible phases, depending on the packing fraction and aspect ratio: isotropic fluid, smectic, nematic, rotator, and solid crystal phases. 291 Apart from fluid and crystal phases, superballs can form rotator phases at intermediate densities. 299,310 Systems of tetragonal parallelepipeds can exhibit liquid-crystalline and cubatic phases, depending on the aspect ratio and packing fraction. 289 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. 302 Truncated cubes exhibit a rich diversity in crystal structures that depend sensitively on the amount of truncation. 312 A study of a wide class of polyhedra revealed that entropy maximization favors mutual alignment of particles along their facets. 307 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. 300 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. 314

C. Maximally Random Jammed States
The fact that M&M candies (spheroidal particles with aspect ratio α ≈ 1.9) were shown experimentally to randomly pack more densely than spheres (φ ≈ 0.66) motivated the development of a modified LS algorithm 136,137 to obtain frictionless MRJ-like packings with even higher densities for other aspect ratios. This included nearlyspherical ellipsoids with φ ≈ 0. 74,126,127 i.e., packing fractions approaching those of the densest 3D sphere packings. Note that these other densest MRJ packings are realizable experimentally. 132 Figure 14 shows separate plots of φ and mean contact number Z as a function of α as predicted by the more refined simulations of Donev et al.. 127 Each plot shows a cusp (i.e., non-differentiable) minimum at the sphere point, and φ versus aspect ratio α possesses a FIG. 14. Packing fraction φ versus aspect ratio α for MRJ packings of 10,000 ellipsoids as obtained in Ref. 127. The semi-axes here are 1, α β , α. The inset shows the mean contact number Z as a function of α. Neither the spheroid (oblate or prolate) or general ellipsoids cases attain their isostatic values of Z = 10 or Z = 12, respectively. 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 d f , which has been referred to as the isostatic conjecture. 321 This conjecture implies a mean contact number Z = 2d f , where d f = 2 for disks, d f = 3 for ellipses, d f = 3 for spheres, d f = 5 for spheroids, and d f = 6 for general ellipsoids. Since d f increases discontinuously with the introduction of rotational degrees of freedom as one makes the particles nonspherical, the isostatic conjecture predicts that Z should have a jump increase at aspect ratio α = 1 to a value of Z = 12 for a general ellipsoid. Such a discontinuity was not originally observed by Donev et al. 126 , rather, they found that jammed ellipsoid packings are hypostatic (or sub-isostatic), Z < 2d f , 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. 127,322 It has been rigorously shown that packings of nonspherical particles are generally not jammed to first order in the "jamming gap" δ [see Sec. III C 2], but are jammed to second order in δ due to curvature deviations from the sphere. 127 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. 127 The reader is referred to recent work on hypostatic jammed 2D packings of noncircular particles. 319 Jiao, Stillinger and Torquato 294 have computed the packing fraction φ MRJ of MRJ packings of binary superdisks in R 2 and identical superballs in R 3 . They found that φ MRJ increases dramatically and nonanalytically as one moves away from the circular-disk or sphere point (p = 1). 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. 323 Note that the geometric-structure approach was used to derive a highly accurate formula for the packing fraction φ MRJ of MRJ binary packings of convex superdisks 181 that is valid for almost all size ratios, relative concentrations and deformation parameter p ≥ 1/2. For the special limit of monodisperse circular disks, this formula predicts φ MRJ = 0.837, which is in very good agreement with the recently numerically discovered MRJ isostatic state 144 with φ MRJ = 0.827.
3D MRJ packings of the four non-tiling Platonic solids (tetrahedra, octahedra, dodecahedra, and icosahedra) were generated using the ASC optimization scheme. 301 The MRJ packing fractions for tetrahedra, octahedra, dodecahedra, and icosahedra are 0.763±0.005, 0.697±0.005, 0.716 ± 0.002, and 0.707 ± 0.002, 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 h(r) that decays to zero asymptotically with same power-law as MRJ spheres, i.e, like −1/r 4 . These results suggest that hyperuniform quasilong-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. 301 In addition, 3D MRJ packings of truncated tetrahedra with an average packing fraction of 0.770 were also generated. 314

D. 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 φ max of packings of nonspherical particles of general shape in R d 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 φ max of a packing of congruent nonspherical particles of volume v P in R d is bounded from above according to where v S is the volume of the largest sphere that can be inscribed in the nonspherical particle and φ S max is the maximal packing fraction of a packing of d-dimensional of identical spheres 140,141 (e.g., φ S max = π/ √ 12 for d = 2 and φ S max = π/ √ 18 for d = 3). 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 φ U max and φ L max , the packing fraction of the densest lattice packing. 140,141,258,286 As discussed below, these densest lattice packings are conjectured to be the densest among all packings. 140,141 Upper bound (30) will also be relatively tight for superballs (superdisks) for deformation parameters p near the sphere (circle) point (p = 1). Dostert et al. 318 recently obtained upper bounds on φ max of translative packings of superballs and 3D convex bodies with tetrahedral symmetry.
1. Ellipsoids   FIG. 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 α = 1 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 α = √ 3 and oblate ellipsoids with α = 1/ √ 3, as shown in the insets. For both α < 1/ √ 3 and α > √ 3, the packing fraction drops off precipitously holding the particle orientations fixed (blue lines). The presently maximal achievable packing fraction φ = 0.7707 . . . is highlighted with a thicker red line, and is constant for α ≤ 1/ √ 3 and α ≥ √ 3; see Ref. 288. The fact that MRJ-like packings of nearly-spherical ellipsoids exist with φ ≈ 0.74 (see Refs. 126 and 127) suggested that there exist ordered ellipsoid packings with appreciably higher densities. The densest known pack-ings of identical 3D ellipsoids were obtained analytically by Donev et al.; 288 see Fig. 15. These are exact constructions and represent a family of non-Bravais lattice packings of ellipsoids with a packing fraction that always exceeds the that of the corresponding densest Bravais lattice packing (φ = 0.74048 . . .) with a maximal packing fraction of φ = 0.7707 . . ., for a wide range of aspect ratios (α ≤ 1/ √ 3 and α ≥ √ 3). In these densest known packings, each ellipsoid has 14 contacting neighbors and there are two particles per fundamental cell.
While a convex "lens" fits snugly within an oblate spheroid of the same aspect ratio, the densest known lens packings are denser for general aspect ratios (except for a narrow range of intermediate values) than their spheroid counterparts, achieving the highest packing fraction of φ = π/4 = 0.7853 . . . in the "flat-lens" limit. 317

Superballs
Exact analytical constructions for candidate maximally dense packings of 2D superdisks were recently proposed for all convex and concave shapes. 138 These are achieved by two different families of Bravais lattice packings such that φ max is nonanalytic at the "circulardisk" point (p = 1) and increases significantly as p moves away from unity. The broken rotational symmetry of superdisks influences the packing characteristics in a non-trivial way that is distinctly different from ellipse packings. For ellipse packings, no improvement over the maximal circle packing fraction (φ max = π/ √ 12 = 0.906899 . . .) is possible, since the former is an affine transformation of the latter. 39 For superdisks, one can take advantage of the four-fold rotationally symmetric shape of the particle to obtain a substantial improvement on the maximal circle packing fraction. By contrast, one needs to use higher-dimensional counterparts of ellipses (d ≥ 3) in order to improve on φ max for spheres. Even for 3D ellipsoids, φ max increases smoothly as the aspect ratios of the semi-axes vary from unity, 288 and hence has no cusp at the sphere point. In fact, 3D ellipsoid packings have a cusp-like behavior at the sphere point only when they are randomly jammed. 126 Jiao, Stillinger and Torquato 139 showed that increasing the dimensionality of a superball from two to three dimensions imbues the optimal packings with structural characteristics that are richer than their 2D counterparts. They obtained analytical constructions for the densest known superball packings for all convex and concave cases, which are certain families of Bravais lattice packings in which each particle has 12 contacting neighbors. For superballs in the cubic regime (p > 1), the candidate optimal packings are achieved by two families of Bravais lattice packings (C 0 and C 1 lattices) possessing two-fold and three-fold rotational symmetry, respectively. For superballs in the octahedral regime (0.5 < p < 1), there are also two families of Bravais lattices (O 0 and O 1 lattices) obtainable from continuous deformations of the fcc lat- tice, which are apparently optimal in the vicinity of the sphere point and the octahedron point, respectively.
The exact maximal packing fraction φ max as a function of deformation parameter p is plotted in Fig. 16. As p increases from unity, the initial increase of φ max is linear in (p − 1) and subsequently φ max increases monotonically with p until it reaches unity as the particle shape tends to the cube. These characteristics stand in contrast to those of the densest known ellipsoid packings (see Fig.  15) whose packing fraction as a function of aspect ratios has zero initial slope and is bounded from above by a value of 0.7707 . . .. 288 As p decreases from unity, the initial increase of φ max is linear in (1 − p). Thus, φ max is a nonanalytic function of p at p = 1. 138 However, the behavior of φ max as the superball shape moves off the sphere point is distinctly different from that of optimal spheroid packings, for which φ max increases smoothly as the aspect ratios of the semi-axes vary from unity and hence has no cusp at the sphere point. 288 These distinctions between superball versus ellipsoid packings result from differences in which rotational symmetries in these two packings are broken. 139 For the small range 0.79248 < p < 1, Ni et al. 309 numerically found lattice packings that are very slightly denser than those of the theoretically predicted O 0 lattices of Ref. 139. All of the results, support the Torquato-Jiao conjecture that the densest packings of all convex superballs are their densest lattice packings. 141

Spherocylinders
The role of curvature in determining dense packings of smoothly-shaped particles is still not very well understood. While the densest packings of 3D ellipsoids are (non-Bravais) lattice packings, 288 the densest packings of superballs appear to lattice packings. 139,140 For 3D spherocylinders with L > 0, the optimal Bravais-lattice packing is very dense and might be one of the actual densest packings. This is because the local principal curvature of the cylindrical surface is zero along the spherocylinder FIG. 17. (color online). Packing fraction versus deformation parameter p for the lattice packings of concave superballs. 139 Insert: a concave superball with p = 0.1, which will becomes a 3D cross at the limit p → 0. axis (i.e., a "flat" direction), and thus, spherocylinders can have very dense lattice packings by an appropriate alignment of the spherocylinders along their axes. The packing fraction of the densest Bravais-lattice packing is given by where L is the length of the cylinder and R is the radius of the spherical caps. This lattice packing corresponds to stacking layers of aligned spherocylinders in the same manner as fcc spheres and hence there are a uncountablyinfinite number of non-Bravais-lattice packings of spherocylinders (in correspondence to the Barlow stackings of spheres 39 ) with the same packing densities (31). Thus, the set of dense nonlattice packings of spherocylinders is overwhelmingly larger than that of the lattice packing.
We will see that these dense packings are consistent with Conjecture 3 described in Sec. V D 5.

Polyhedra
About a decade ago, very little was known about the densest packings of polyhedral particles. The difficulty in obtaining dense packings of polyhedra is related to their complex rotational degrees of freedom and to the non-smooth nature of their shapes. It was the investigation of Conway and Torquato of dense packings of tetrahedra 290 that spurred the flurry of activity over several years to find the densest packings of tetrahedra, 140,141,293,[295][296][297][298]324 which in turn led to studies of the densest packings of other polyhedra and many other nonspherical convex and concave particles. 300,304,306,307,312,317 Torquato and Jiao 140,141 employed the ASC optimization scheme to determine dense packings of the non-tiling Platonic solids and of all of the Archimedean solids. For example, they were able to find the densest known packings of the octahedra, dodecahedra and icosahedra (three Unlike the densest tetrahedron packing, which must be a non-Bravais lattice packing, the densest packings of the other non-tiling Platonic solids found by the algorithm are their previously known densest (Bravais) lattice packings 258,286 ; see Fig. 18. These simulation results as well as other theoretical considerations led them to general organizing principles concerning the densest packings of a class of nonspherical particles, which are described in Sec. V D 5.
Conway and Torquato 290 showed that the maximally dense packing of regular tetrahedra cannot be a Bravais lattice. Among other non-Bravais lattice packings, they obtained a simple, uniform packing of "dimers" (two particles per fundamental cell) with packing fraction φ = 2/3. A uniform packing has a symmetry (in this case a point inversion symmetry) that takes one tetrahedron to another. A dimer is composed of a pair of regular tetrahedra that exactly share a common face. They also found a non-Bravais lattice (periodic) packings of regular tetrahedra with φ ≈ 0.72, which doubled that the density of the corresponding densest Bravaislattice packing (φ = 18/49 = 0.367 . . .), which was the record before 2006. This work spurred many studies that improved on this density. 140,141,293,[295][296][297][298]324 Kallus et al. 296 found a remarkably simple "uniform" packing of tetrahedra with high symmetry consisting of only four particles per fundamental cell (two "dimers") with packing fraction φ = 100/117 = 0.854700 . . .. Torquato and Jiao 324 subsequently presented an analytical formulation to construct a three-parameter family of dense uniform dimer packings of tetrahedra again with four particles per fundamental cell. Making an assumption about one of these parameters resulted in a two-parameter family, including those with a packing fraction as high as φ = 12250/14319 = 0.855506 . . .. Chen et al. 297 used the full three-parameter family to obtain the densest known dimer packings of tetrahedra with the very slightly higher packing fraction φ = 4000/4671 = 0.856347 . . . (see Fig. 19). Whether this is the optimal packing is an open question for reasons given by Torquato and Jiao. 295 The fact that the first nontrivial upper bound on the maximal packing fraction is infinitesimally smaller that unity 255 (φ max ≤ 1 − 2.6 × 10 −25 ) is a testament to the difficulty in accounting for the orientations of the tetrahedra in a dense packing.

Organizing Principles for Convex and Concave Particles
Torquato and Jiao 140,141 showed that substantial faceto-face contacts between any of the centrally symmetric Platonic and Archimedean solids allow for a higher packing fraction. They also demonstrated that central symmetry enables maximal face-to-face contacts when particles are aligned, which is consistent with the densest packing being the optimal lattice packing. The aforementioned simulation results, upper bound, and theoretical considerations led to three conjectures concerning the densest packings of polyhedra and other nonspherical particles in R 3 ; see Refs. 140, 141, and 295. Conjecture 1: The densest packings of the centrally symmetric Platonic and Archimedean solids are given by their corresponding optimal (Bravais) lattice packings.

Conjecture 2:
The densest packing of any convex, congruent polyhedron without central symmetry generally is not a (Bravais) lattice packing, i.e., the set of such polyhedra whose optimal packing is not a lattice is overwhelmingly larger than the set whose optimal packing is a lattice.

Conjecture 3:
The densest packings of congruent, centrally symmetric particles that do not possess three equivalent principle axes (e.g., ellipsoids) are generally not (Bravais) lattice packings.
Conjecture 1 is the analog of Kepler's sphere conjecture for the centrally symmetric Platonic and Archimedean solids. On the experimental side, it has been shown 305 that such silver polyhedral nanoparticles self-assemble into the conjectured densest lattice packings of such shapes. Torquato and Jiao 141 have also commented on the validity of Conjecture 1 to polytopes in four and higher dimensions. The arguments leading to Conjecture 1 also strongly suggest that the densest packings of superballs and other smoothly-shaped centrally symmetric convex particles having surfaces without "flat" directions are given by their corresponding optimal lattice packings. 141 Consistent with Conjecture 2, the densest known packing of the non-centrally symmetric truncated tetrahedron is a non-lattice packing with packing fraction φ = 207/208 = 0.995192 . . ., which is amazingly close to unity and strongly implies its optimality. 300 We note that non-Bravais-lattice packings of elliptical cylinders (i.e., cylinders with an elliptical basal face) that are denser than the corresponding optimal lattice packings for any aspect ratio greater than unity have been constructed, 284,325 which is consistent with Conjecture 3. A corollary to Conjecture 3 is that the densest packings of congruent, centrally symmetric particles that possess three equivalent principle axes (e.g., superballs) are generally Bravais lattices. 139 Subsequently, Torquato and Jiao 304 generalized the aforementioned three conjectures in order to guide one to ascertain the densest packings of other convex nonspherical particles as well as concave shapes. These generalized organizing principles are explicitly stated as the following four distinct propositions: Proposition 1: Dense packings of centrally symmetric convex, congruent particles with three equivalent axes are given by their corresponding densest lattice packings, providing a tight density lower bound that may be optimal. Proposition 2: Dense packings of convex, congruent polyhedra without central symmetry are composed of centrally symmetric compound units of the polyhedra with the inversion-symmetric points lying on the densest lattice associated with the compound units, providing a tight density lower bound that may be optimal. Proposition 3: Dense packings of centrally symmetric concave, congruent polyhedra are given by their corresponding densest lattice packings, providing a tight density lower bound that may be optimal. Proposition 4: Dense packings of concave, congruent polyhedra without central symmetry are composed of centrally symmetric compound units of the polyhedra with the inversion-symmetric points lying on the densest lattice associated with the compound units, providing a tight density lower bound that may be optimal.
Proposition 1 originally concerned polyhedra. It is generalized here to include any centrally symmetric convex, congruent particle with three equivalent axes to reflect the arguments put forth by Torquato and Jiao. 140,304 All of the aforementioned organizing principles in the form of conjectures and propositions have been tested in Ref. 304 against a comprehensive set of both convex and concave particle shapes, including but not limited to the Platonic and Archimedean solids, 19,140,300 Catalan solids, 326 prisms and antiprisms, 326 Johnson solids, 326 cylinders, 284,325 lenses 317 , truncated cubes 312 and various concave polyhedra. 326 These general organizing principles also enable one to construct analytically the dens-est known packings of certain convex nonspherical particles, including spherocylinders, and square pyramids and rhombic pyramids. 304 Moreover, it was shown how to apply these principles to infer the high-density equilibrium crystalline phases of hard convex and concave particles. 304 We note that the densest known 2D packings of a large family of 2D convex and concave particles (e.g., crosses, curved triangles, moon-like shapes) fully adhere to the aforementioned organizing principles. 308 Interestingly, the densest known packing of any identical 3D convex particle studied to date has a density that exceeds that of the optimal sphere packing value φ S max = π/ √ 18 = 0.7408 . . .. These results are consistent with a conjecture of Ulam, who proposed to Martin Gardner 327 , without any justification, that the optimal packing fraction for congruent sphere packings is smaller than that for any other convex body. There is currently no proof of Ulam's conjecture.

Dense Packings of Tori
We have seen that the preponderance of studies of dense packings of nonspherical shapes in R 3 have dealt with convex bodies that are simply connected and thus topologically equivalent to a sphere. However, much less is known about dense packings of multiply connected solid bodies. Gabbrielli et al. 313 investigated the packing behavior of congruent ring tori in R 3 , which are multiply connected non-convex bodies of genus one, as well as horn and spindle tori. Guided by the aforementioned organizing principles, they analytically constructed a family of dense periodic packings of individual tori and found that the horn tori as well as certain spindle and ring tori can achieve a packing density not only higher than that of spheres (i.e., π/ √ 18 = 0.7404 . . .) but also higher than the densest known ellipsoid packings (i.e., 0.7707 . . .). Moreover, they studied dense packings of clusters of pairlinked ring tori (i.e., Hopf links), which can possess much higher densities than corresponding packings consisting of unlinked tori; see Fig. 20 for a specific example.

VI. CHALLENGES AND OPEN QUESTIONS
Packing problems are fundamental and their solutions are often profound. This perspective provides only a glimpse into the richness of packing models and their capacity to capture the salient structural and physical properties of a wide class of equilibrium and nonequilibrium condensed phases of matter that arise across the physical, mathematical and biological sciences as well as technological applications.
Not surprisingly, there are many challenges and open questions. Is it possible to prove a first-order freezing transition in 3D hard-sphere systems? Is the FCC lattice provably the entropically favored state as the maximal density is approached along the stable crystal branch in this same system? Can one devise numerical algorithms that produced large disordered strictly jammed isostatic sphere packings that are rattler-free? Is the true MRJ state rattler-free in the thermodynamic limit? Can one prove the Torquato-Stillinger conjecture that links strictly jammed sphere packings to hyperuniformity? Can one identify incisive order metrics for packings of nonspherical particles as well as a wide class of manyparticle systems that arise in molecular, biological, cosmological, and ecological systems? What are the appropriate generalizations of the jamming categories for packings of nonspherical particles? Is it possible to prove that the densest packings of the centrally symmetric Platonic and Archimedean solids are given by their corresponding optimal lattice packings? Upon extending the geometricstructure approach to Euclidean dimensions greater than three, do periodic packings with arbitrarily large unit cells or even disordered jammed packings ever provide the highest attainable densities? Can one formulate a disordered sphere-packing model in R d that can be rigorously shown to have a packing fraction that exceeds φ = 1/2 d , the maximal value achievable by the ghost RSA packing? In view of the wide interest in packing problems across the sciences, it seems reasonable to expect that substantial conceptual advances are forthcoming. Weis with whom I have collaborated on topics described in this review article. I am very thankful to Jaeuk Kim, Duyu Chen, Zheng Ma, Yang Jiao and Ge Zhang for comments that greatly improved this article. The author's work on packing models over the years has been supported by various grants from the National Science Foundation.