Jamming in hard sphere and disk packings

Hard-particle packings have provided a rich source of outstanding theoretical problems and served as useful starting points to model the structure of granular media, liquids, living cells, glasses, and random media. The nature of "jammed" hard-particle packings is a current subject of keen interest. Elsewhere, we introduced rigorous and efficient linear-programming algorithms to assess whether a hard-sphere packing is locally, collectively, or strictly jammed, as defined by Torquato and Stillinger [J. Phys. Chem. B 105, 11849 (2001)]. One algorithm applies to ideal packings in which particles form perfect contacts. Another algorithm treats the case of jamming in packings with significant interparticle gaps. We have applied these algorithms to test jamming categories of ordered lattices as well as random packings of circular disks and spheres under periodic boundary conditions. The random packings were produced computationally with a variety of packing generation algorithms, all of which should, in principle, produce at least collectively jammed packings. Our results highlight the importance of jamming categories in characterizing particle packings. One important and interesting conclusion is that the amorphous monodisperse sphere packings with density 0.64 were for practical purposes strictly jammed in three dimensions, but in two dimensions the monodisperse disk packings at previously reported "random close packed" densities of 0.83 were not even collectively jammed. On the other hand, amorphous bidisperse disk packings with density of 0.84 were virtually strictly jammed. This clearly demonstrates one cannot judge "stability" in packings based solely on local criteria. Numerous interactive visualization models are provided on the authors' webpage.© 2004 American Institute of Physics.

Contents

• I. INTRODUCTION
• II. BACKGROUND AND METHODOLOGY
• III. PACKING GENERATION ALGORITHMS
• IV. RESULTS
  • A. Periodic lattice packings
  • B. Periodic random packings
• V. DISCUSSION, FUTURE WORK, AND CONCLUSIONS
• ACKNOWLEDGMENTS
• REFERENCES
• FIGURES
• TABLES
• FOOTNOTES

Packings of hard particles interacting only with infinite repulsive pairwise forces on contact are applicable as models of complex many-body systems because repulsive interactions are the primary factor in determining their structure. Hard-particle packings are therefore widely used as simple models for granular materials,^1,2 glasses,³ liquids,⁴ and other random media,⁵ to mention a few examples. Furthermore, hard-particle packings, and especially hard-sphere packings, have inspired mathematicians and been the source of numerous challenging (many still open) theoretical problems.⁶

We focus our attention in this paper on the venerable idealized hard-sphere model, i.e., the only interparticle interaction is an infinite repulsion for overlapping particles, since this enables us to be precise about the important concept of "jamming." In particular, a hierarchical classification scheme for jammed packings into locally, collectively, and strictly jammed packings was proposed in Ref. 7. This classification is closely related to the concepts of "rigid" and "stable" packings found in the mathematics literature.^8,9 The idealized hard-sphere model is in a sense the "Ising model" for studying a variety of hard-particle physical systems, and the importance of understanding it in detail cannot be overstated. The term jamming is used in a different sense in the modeling of granular media, which includes effects such as friction, adhesion, particle deformability, etc., and, by definition, hard-sphere systems do not include these effects. It is also important to note that we do not discuss dynamical effects in hard-particle packings. In the present work, hard-sphere jamming is presented from a rigorous perspective that focuses on the geometry of the final packed states. We note that extensions of this work to packings of nonspherical particles (such as ellipses, ellipsoids, or spherocylinders) are possible and the subject of current and future research.

There are still many important and challenging questions open even for the simplest type of hard-particle packings, i.e., monodisperse packings of smooth perfectly impenetrable spheres. One important category of open problems pertains to the enumeration and classification of both ordered and disordered jammed circular disk and sphere packings for the various jamming categories described in the following. Since one cannot enumerate all possible packings even for a small number of particles, it is desirable to devise a small set of parameters that can characterize packings well. Two important scalar properties of packings are the density (packing fraction) and order metric . For any two states X and Y, _X>_Y implies that state X is to be considered as more ordered than state Y. Candidates for such an order metric include various translational and orientational order parameters,⁵ but the search for better order metrics is still very active.

Figure 1 from Ref. 10 shows a conjectured region of feasible hard-sphere packings in the – plane. It is clear that only a small subset of this feasible region will be occupied by jammed packings (for a given jamming category), as schematically indicated in Fig. 1. Several limit points in this region are particularly interesting.

Figure 1.

(1) Point A corresponds to the lowest-density jammed packing, and its location strongly depends on the jamming category used. It can be shown that there are zero-density locally jammed disk and sphere packings (see references and discussion in Ref. 11). However, for collectively and strictly jammed packings, it is not known what are the lowest possible densities.

(2) Point B corresponds to the most dense jammed packing. It has of course already been identified to be a triangular packing for disks and the FCC/HCP variant lattice for spheres. But much less is known about polydisperse packings,^11,12 or packings of nonspherical particles.

(3) MRJ point represents the maximally random jammed (MRJ) state,¹⁰ which has recently supplanted the ill-defined "random close packed" (RCP) state. The RCP state was widely believed to have a packing fraction 0.63–0.64 in three dimensions. The MRJ state is the most disordered jammed packing in a given jamming category (locally, collectively, or strictly jammed). The MRJ state is well-defined for a given jamming category and choice of order metric.

Numerical algorithms have long been the primary tool for studying random packings quantitatively. In a separate paper,¹³ we introduced two algorithms to assess whether a hard-sphere packing is locally, collectively, or strictly jammed.⁷ The first algorithm targets packings with perfect interparticle contacts, while the second allows for significant interparticle gaps. Both algorithms are based on linear programming and are applicable to both ordered and disordered disk and sphere packings of arbitrary polydispersity. In Ref. 13 we give a complete description of the algorithms. Here we demonstrate their applicability, usefulness, and efficiency in analyzing large disordered packings, as produced by various packing generation algorithms. Algorithms that generate large-scale hard-particle packings are very important, especially because experimental hard-particle configurations are difficult to obtain and are limited in applicability. Of particular interest are stochastic algorithms aimed at producing random (disordered) packings.

Through numerical investigations, we show here that several previously used packing algorithms generate collectively jammed packings under appropriate conditions. In particular, we study in detail monodisperse sphere as well as monodisperse and bidisperse disk packings produced by the Lubachevsky–Stillinger packing algorithm.¹⁴ We also tested a sample of monodisperse sphere and bidisperse disk packings produced by the algorithm described in Ref. 15, as well as monodisperse sphere packings produced by the Zinchenko packing algorithm,¹⁶ and observed similar behavior as for the Lubachevsky–Stillinger packings.

Our testing of these packings enables us to arrive at several important conclusions. First, we find that the amorphous monodisperse sphere packings (with covering fraction, or density, 0.64) and bidisperse disk packings (0.84) are practically strictly jammed (though not in the ideal sense). Second, we observe that large monodisperse disk packings are invariably highly crystalline (0.88) and are only collectively jammed. Previously reported¹⁷ low covering fractions for "random close packed" disks of 0.82–0.84 were not even found to be collectively jammed. This conclusion clearly demonstrates that the distinctions between the different jamming categories are important and one cannot judge "stability" in packings based solely on local criteria, as has been done extensively in the literature.^18,19,20 Preliminary investigations with an extension of the Lubachevsky–Stillinger algorithm indicate that it is possible to produce ideal strictly jammed packings, which is important in order to eliminate finite-size boundary effects, especially for monodisperse disk packings.

In Sec. II, we introduce important notation, definitions, and review basic concepts. In Sec. III, we describe various algorithms that are used to generate random packings. We will analyze the resultant packings. In Sec. IV we discuss the numerical implementation, and provide results for ordered periodic lattice packings and random packings. Finally, we conclude with a discussion of the results and future directions of investigation in Sec. V.

Here we briefly summarize some of the essential notation, problem statements, and methods, as described in detail in Ref. 13. We consider a sphere packing in Euclidean d-dimensional space ^d, characterized by the positions of the sphere centers R = (r₁,,r_N),

where the diameter of the ith sphere is D_i, and here we focus on monodisperse packings, i.e., packings where all the spheres are identical, D_i = D. Our perspective on jamming focuses on the set _R of configurations around a particular initial configuration R reachable via continuous displacements of the spheres R(t), subject to nonoverlapping constraints and certain boundary conditions. Here t is a time-like parameter, and we will often drop it for brevity, but it should be kept in mind that any change in configuration we consider must be reachable via a continuous deformation. If the extent of _R is small, in the sense that only small continuous displacements of the particles from their initial configurations are possible for all R_R, the packing is considered jammed. The natural length scale defining the meaning of "small" is the typical size of the particles, or the size of the interparticle gaps, depending on the context and the type of packing under consideration. In a jammed ideal packing, which has perfect interparticle contacts, the particles cannot at all be displaced continuously from their current configuration (modulo trivial rigid-body motions). By changing the boundary conditions, we get several different categories of jamming, namely local, collective, and strict jamming.⁷ We briefly review these definitions for the convenience of the reader in the following. We consider first ideal packings, and discuss interparticle gaps in more detail as an extension.

We specialize these jamming definitions for periodic sphere packings for concreteness, but packings in a concave hard-wall containers can also be considered. Periodic (repetitive) packings are characterized by a unit cell and a lattice = {₁,,_d}, where _i are linearly independent lattice vectors. We additionally allow the lattice to continuously change by (t) as the particles displace, where = ^T = ()^–1 is the symmetric macroscopic strain tensor.^21,22

Finite systems of spheres are characterized as follows:

(1) Locally jammed: Each particle in the system is locally trapped by its neighbors, i.e., it cannot be translated while fixing the positions of all other particles. Each sphere simply has to have at least d + 1 contacts with neighboring spheres, not all in the same d-dimensional hemisphere.

(2) Collectively jammed: Any locally jammed configuration where all finite subsets of particles are trapped by their neighbors. For periodic boundary conditions, collective jamming implies that there is no nonvanishing continuous periodic displacement of the particles R(t) that maintains impenetrability other than trivial uniform translations of the packing, while keeping the lattice fixed, (t) = 0.

(3) Strictly jammed: Any collectively jammed configuration that disallows all globally uniform volume-nonincreasing deformations of the system boundary. For periodic packings, the boundary is in fact the lattice, and strict jamming implies that there is no nonvanishing continuous periodic R(t) which maintains impenetrability other than trivial uniform translations of the packing, even if we allow a volume-nonincreasing continuous lattice deformation (t) (this translates to a strain tensor with a nonpositive trace).

Observe that these jamming categories are ordered hierarchically, with local being a prerequisite for collective and similarly collective being a prerequisite for strict jamming. It should be mentioned that jammed random particle packings produced experimentally or in simulations typically contain a small population of "rattlers," i.e., particles trapped in a cage of jammed neighbors but free to move within the cage. For present purposes we shall assume that these have been removed before considering the (possibly) jammed remainder (subpacking).

In Ref. 13, we presented a randomized linear programming (LP) algorithm to test whether a given packing is jammed or not, for each of the above-given jamming categories. The essential ingredient of this algorithm is to apply a randomly selected load (i.e., a force) on each particle (locally, or collectively) and then solve a linear program which takes into account a linearized version of the impenetrability constraints between neighboring particles to find whether (and how) the particles displace (and possible the lattice deforms) in order to support this applied load. If the particles do not displace then we apply the load of opposite sign and repeat the test. If the particles do not displace again, then the ideal packing under consideration is jammed.

Computer-generated packings, which we analyze, are never ideal and there are always small interparticle gaps between some particles, typically much less than a percent of the typical particle size D. One can safely consider such packings within the framework of ideal packings, with minor modifications to the algorithm, as described in detail in Ref. 13. However, the either-or character of the above-mentioned jamming criteria is often too restrictive or specialized when analyzing large disordered packings with possibly larger interparticle gaps, where particle displacements may be comparable to the typical particle size. Therefore, we investigate ways to study jamming in this practical sense for such nonideal packings. We focus here on trying to judge the extent of _R by trying to displace the spheres away from their current position by as much as possible. In Ref. 13 we describe an algorithm based on linear programming to do just this, and the basic idea is to repeatedly apply a random load on the particles, solve several linear programs, and displace the particles by as much as possible while still avoiding overlap, until the particles rearrange and form contacts that actually support the applied load. This is repeated for several random loads, in the hope of exploring _R along several directions. We can then actually quantitatively report the average/maximal displacement of the particles that was observed, and use this instead of a binary classification into packings which are jammed and not jammed. The numerous intricacies of the algorithm are discussed in detail in Ref. 13.

We have implemented these algorithms to test for jamming in sphere packings and here we apply them to monodisperse and bidisperse packings under periodic boundary conditions. We present some representative but nonexhaustive results for several periodic ordered lattice packings as well as random packings obtained via the Lubachevsky–Stillinger packing algorithm.¹⁴ We plot linear unjamming motions as suitably scaled "velocity" fields, showing the direction in which the particles can move (along a straight path with in this linear algorithm) without violating impenetrability. Numerous more illustrative and interactive Virtual Reality Modeling Language (VRML) animations can be viewed on our webpage.²³

We produced most packings using the Lubachevsky–Stillinger compression algorithm¹⁴ with periodic boundary conditions. This algorithm is essentially a hard-sphere molecular dynamics in which the spheres grow in size during the course of the simulation at a certain expansion rate, until a final state with diverging collision rate is reached.

We also obtained sample monodisperse sphere and bidisperse disk packings from the authors of Ref. 15. These packings are not of perfectly hard spheres, but rather soft spheres interacting via repulsive potentials when there is overlap between the cores of diameter D. They use energy minimization for harmonic and Hertzian potentials, descending to an energy minimum using the conjugate gradient algorithm from a random initial configuration (i.e., a rapid quench from T = to T = 0). The packings we analyzed were just above the "jamming threshold" density _c, meaning that there was only very small (less than 10^–5D) overlap between the outer cores. We therefore simply scaled the sizes of the particles by a factor very close to unity to obtain overlap-free hard-sphere packings. Since the jamming threshold densities found in Ref. 15 were very close to the final densities produced by the Lubachevsky–Stillinger algorithm (with reasonably large compression rates), we expected these packings to behave very similarly, and have confirmed this with computational tests. Therefore, here we focus on and present the results for the Lubachevsky–Stillinger packings.

We also had available disordered three-dimensional packings produced with the contact network building Zinchenko packing algorithm,¹⁶ and confirmed that they behaved like the packings produced by the other algorithms. Unfortunately, we do not know of a two-dimensional implementation of this algorithm, and it is important to develop one in the future and see whether it too produces near-triangular packings.

More detailed results will be given shortly, but we want to point out here that none of these algorithms produces truly strictly jammed packings a priori. Indeed, the packings that that we tested were never truly strictly jammed. This is not surprising because none of them incorporates deformations of the periodic lattice, but rather, they all use a fixed (typically cubical) unit cell. It is not hard to incorporate boundary deformations into these algorithms, and we are presently working on such extensions. In particular, the Lubachevsky–Stillinger algorithm can easily incorporate a deforming lattice in the spirit of Parinello–Rahman molecular dynamics.²⁴ We have in fact implemented such an extended Lubachevsky–Stillinger algorithm and used it to produce a priori strictly jammed packings. Details of this work in progress will be given in future papers, and here we will only analyze some of the final packings produced by the algorithm. In packing algorithms based on energy minimization, as in Ref. 15, one need only include the strain as part of the degrees of freedom in order to allow relaxation of the lattice and produce strictly jammed packings. The same is true of the Zinchenko packing algorithm.

On the other hand, all of these algorithms seem to produce collectively jammed packings in both two and three dimensions, excluding rattlers and allowing for appropriate numerical tolerances. This can be proved rigorously for the Zinchenko packing algorithm, and under certain additional assumptions also for the energy minimization algorithm. In principle, only locally jammed configurations are possible final states for the Lubachevsky–Stillinger algorithm since they give infinite collision rates, however, we believe that local configurations are unstable attractors for this algorithm and in fact under appropriate conditions all final states have a collectively jammed subpacking, excluding rattlers. We have recently devised a way to dynamically verify jamming during the packing process in the Lubachevsky–Stillinger algorithm for both packings of spheres and ellipsoids, however, details will be given in future publications.

We have developed an efficient numerical implementation of the randomized LP algorithm using the primal-dual interior-point algorithm in the LOQO optimization library.²⁵ Both FORTRAN 95 codes which directly invoke the LOQO library, and Algebraic Modeling Programming Language (AMPL) models have been developed, along with VRML visualization tools. Illustrations of results obtained using these implementations are given throughout this paper. We have applied the algorithms to test for the different jamming categories in practice and verified their utility and efficiency. Although reporting exhaustive results is not the primary aim of this work, in this section we present some relevant results for both ordered and disordered periodic packings. We have analyzed disordered packings produced by a variety of packing algorithms, namely the Lubachevsky–Stillinger packing algorithm,¹⁴ an energy minimization algorithm as presented in Ref. 15, as well as the Zinchenko packing algorithm.¹⁶

Table 1 in Ref. 7 gives a classification of some common simple lattice packings into jamming categories for hard-wall boundary conditions. Table I modifies this for periodic boundary conditions. The results in principle will depend on the choice of unit cell, so the terminology "lattice XXX is YYY jammed" is used loosely here. We illustrate some unjamming motions for lattice disk packings in Figs. 2 and 3.

Figure 2.

Figure 3.

Here we just point out for the curious that the triangular lattice is not the only strictly jammed ordered disk packing; two other examples are shown in Fig. 4.¹¹ It can be shown that one can remove at most one quarter of the disks from a triangular lattice packing and still maintain strict jamming. Using the Lubachevsky–Stillinger packing algorithm for small packings, we recently found a new family of strictly jamming packings obtained by reinforcing with triangular regions a particular tiling of the plane with three congruent pentagons. An example is shown in Fig. 4.

Figure 4.

We also tested a sample of periodic random packings in two and three dimensions. Both monodisperse and bidisperse packings were studied. The main reason for including bidisperse packings in this preliminary study is that monodisperse disk packings crystallize easily, forming large ordered almost-triangular domains (grains) with high packing fraction 0.88. This is because in two dimensions the locally densest configuration coincides with the globally densest triangular lattice, unlike in three dimensions, where the locally optimal (tetrahedral) configuration cannot tile three-dimensional space.⁵ It is only by introducing polydispersity that one can produce disk packings with no apparent (or little) short-range order (i.e., amorphous), as can be determined by, for example, bond-orientational order metrics,⁵ and in particular, the local Q₆ order metric. We used an equimolar mixture of disks with diameter ratio of 1.4 as done in Ref. 15. For amorphous monodisperse three-dimensional packings the typical packing fraction is around 0.64, and such a packing is shown in Fig. 5. For the aforementioned amorphous binary disk packings 0.84, and such a packing is illustrated in Fig. 6.

Figure 5.

Figure 6.

In a truly disordered (generic) packing, it is expected that the average number of interparticle contacts per particle (coordination number) will be Z = 2d (more precisely, twice the number of degrees of freedom per particle). Thus, it is expected that Z = 4 in two dimensions. However, collectively jammed monodisperse disks packings are rather dense (0.86–0.88) and crystalline and they have Z5.5 (This should be compared to Z = 6 for the triangular crystal). Disordered bidisperse disk packings do have Z4, and similarly in three dimensions monodisperse packings have Z6, consistent with an assumption of generic character. However, the exact number one gets depends rather sensitively on the criterion for assigning contacts and on whether rattling particles are excluded or not. Future work will give a more detailed and careful investigation of coordination number distribution in disordered packings. For this work, it is important to note that a large packing must have Z2d in order to be collectively or strictly jammed, and Zd + 1 to be locally jammed. We also point out that our algorithm to test for jamming does not depend sensitively on the criterion for selecting contacts.¹³

Although most of the packings we analyzed had small interparticle gaps and can also be studied within the framework of ideal packings and classified as jammed or not jammed, we instead consider them nonideal and explicitly deal with the interparticle gaps. We wish to stress that the results to follow are not averages over many packings with the same number of spheres/disks, but rather they are results for particular packings produced by the Lubachevsky–Stillinger algorithm. These packings seemed to be typical of the types of packings produced by the algorithm under a relatively wide range of expansion rates and packing sizes. We therefore believe that the numbers presented here serve well as a semiquantitative illustration of the behavior of random disk and sphere packings commonly used in many computational studies. The primary reason we do not give averaged results this is that detailed average results should be given only once it is determined what quantitative metric of jamming is physically appropriate (which is likely to be different for different types of packings and different applications), and results should also be correlated with more characteristics of the packings (i.e., not just the covering fraction) and to various relevant parameters of the algorithm used to generate the packing.²⁶

As a quantitative measure of jamming in these packings, we report the average particle displacement achieved during random loading. This choice is not ideal, and attaching a physical picture to the numbers is difficult. Furthermore, deciding when to terminate the Lubachevsky–Stillinger algorithm is nontrivial and we used the principle of allowing a certain number of binary collisions per particle and also limiting the total computational time, which results in larger packings not being as "well-packed" as smaller packings. Visualization of the resulting particle displacements is still the best way to analyze the results. For example, rattlers often contribute most to the average displacement for packings which might be "more jammed" if the rattlers are removed. Moreover, although an entry in Table IV below might say that the average displacement for a monodisperse disk packing was only 10% of the particle size, the character of the particle motion might be such that very significant rearrangements happen in the packing because grain boundaries move (see Fig. 8), and this has to be seen to be appreciated. We will share our VRML visualizations with interested readers, and many examples are provided on our webpage.²³

Figure 8.

Another statistic we report is the time (in seconds) spent by the AMPL implementation (with some FORTRAN) of the testing algorithm on a typical personal computer. (More precisely, calculations were performed on a 1666 MHz AMD Athlon PC running Linux.) Since most of the computational time is spent in LOQO, similar running times are typical of the FORTRAN codes as well. For each packing, we applied three different random loads (with opposite orientations), and for each load we successively solved three linear programs (so a total of 18 linear programs for each packing). The running times to follow should not be taken as a measure of the scaling of the LP solver computational effort with the number of spheres, but rather as typical runtimes for some representative packing sizes. This is because the computational effort depends nontrivially on many of the parameters in the algorithm, and on the exact implementation. We are currently developing more efficient and robust implementations of these algorithms, for both packings of disks/spheres and ellipses/ellipsoids.

Qualitatively different results were observed for the amorphous monodisperse sphere packings and binary disk packings, and the polycrystalline monodisperse disk packings.

For the amorphous packings, we give results in Table II for monodisperse packings in three dimensions and in Table III for bidisperse packings in two dimensions, with similar trends. In general, these packings were collectively jammed, in the sense that only small (average) displacements of the particles are possible. The small feasible displacements are mostly due to rattlers and/or early termination of the packing algorithm and we believe that any true final Lubachevsky–Stillinger packing with infinite collision rate will in fact have an ideal collectively jammed subpacking (similarly for the other packing algorithms). The packings were not strictly jammed for small system sizes, however, the magnitude of the feasible displacements decreased as the packings became larger, and therefore large amorphous packings were apparently collectively and strictly jammed. This can be understood by thinking of the distinction between collective and strict jamming as a boundary effect: As the packings become larger the boundary effects diminish. Therefore, even though none of the packing algorithms is meant to produce strictly jammed packings a priori, they do so for large amorphous packings.

Importantly, very different results were observed for monodisperse disk packings, which are invariantly nearly triangular (i.e., crystalline). We wish to point out that crystallization into a triangular lattice poses a convergence obstacle for the Lubachevsky–Stillinger algorithm since near triangular regions have very high collision rates even when the disks' diameters are not at their maximal value. Therefore it was only for monodisperse disk packings that some of the final packings were not collectively jammed (large particle rearrangements were possible near grain boundaries). Most packings were however collectively jammed just as for amorphous packings and we present results for these in Table IV. By using certain tricks in the Lubachevsky–Stillinger algorithm, such as collections of frozen particles or very large expansion rates, one can obtain apparently "jammed" amorphous monodisperse disk packings near a packing fraction 0.83. However, due to numerical instabilities or the presence of an artificial boundary of fixed disks, these packings were not collectively jammed, as illustrated in Fig. 7. One of the important observations is that none of the large Lubachevsky–Stillinger monodisperse disk packings were strictly jammed. In fact, typical grain boundaries are very fragile under shear, and so even for the large packings significant rearrangements of the grain boundaries are feasible, as illustrated in Fig. 8.

Figure 7.

It is also important to verify that any packing algorithm claimed to produce jammed packings can indeed produce jammed ideal packings, in the sense that all tolerances in the test for jamming can be tightened progressively as the numerical accuracy is increased and the convergence criteria in the packing algorithm are tightened. We demonstrate this for collective jamming in monodisperse sphere packings in Table V. The corresponding results for strict jamming, given in Table VI, illustrate that the (traditional) Lubachevsky–Stillinger packings do not have a strictly jammed ideal subpacking, but are practically strictly jammed for large system sizes. This is unlike monodisperse disk packings, which are far from being strictly jammed, as illustrated in Table VII.

Very recently, we have implemented an extension of the Lubachevsky–Stillinger packing algorithm in which the lattice deforms during the molecular dynamics run, as dictated by the collisional (contact) stress induced by the particle collection. Details of the algorithm and the packings it produces will be given elsewhere, but a short description can be found in Ref. 13. For relatively small numbers of particles, this algorithm typically produces truly strictly jammed packings, and for these packings is similar for both collective and strict jamming. The algorithm produces similar amorphous packings (in packing fraction and disorder) to the original Lubachevsky–Stillinger algorithm, however, for monodisperse disks it frequently terminates with completely crystal packings, and also produces complete triangular lattices with special types of defects, such as monovacancies and peculiar "dislocation cores." One such strictly jammed disk packing is shown in Fig. 9. Investigation of these strictly jammed disk packings as well as extensions of other packing algorithms to allow for deforming boundaries are being carried out at present. We note in passing that we also used the extended Lubachevsky–Stillinger algorithm to try the shrink-and-bump heuristic¹³ to test for strict jamming by also allowing the lattice to deform while the particles bump around. This seemed to detect disordered packings which are not strictly jammed, however, the test is significantly slower than the linear programming algorithm and is also very heuristic and much less reliable.

Figure 9.

Our results have important implications for the classification of random disk and sphere packings and suggest a number of interesting avenues of inquiry for future investigations. Random disk packings are less well understood than sphere packings. The tendency of disk packings to "crystallize" (to form ordered, locally dense domains) at sufficiently high densities is well established. For example, Quickenden and Tan experimentally estimated the packing fraction of the "random close packed" (RCP) state to be 0.83 and found that the packing fraction could be further increased until the maximum value of = 0.906 is achieved for the triangular lattice packing.²⁷ By contrast, random sphere packings at in the range 0.63–0.66 cannot be made more dense.

Our recent understanding of the ill-defined nature of random close packing and of jamming categories raises serious questions about previous two-dimensional studies, particularly the stability of such packings. Our present study suggests that disordered random disk packings are not collectively jammed at 0.83; at best they may be locally jammed. This brings into question the previous widespread belief that the two-dimensional analog of the RCP sphere-packing state has density about 0.82–0.83.¹⁷ Collectively jammed disk packings seem to have significantly higher densities 0.88 and consist of large triangular grains, but even at such high densities they are not strictly jammed. An interesting question is whether the grain size becomes small compared to the system size for large collectively jammed disk packings, or whether the appearance of grain boundaries is in fact a finite-size boundary effect. It may be that the preponderance of collectively and strictly jammed large disk packings are very crystalline, with a distribution of the local bond-orientational parameter Q₆ (see Ref. 5) highly peaked around some relatively large value. Furthermore, it is important to ascertain if the strong distinction between only collectively and strictly jammed disk packings persists in the limit of very large packings. Careful investigations of very large collectively and strictly jammed disk packings produced with a variety of packing algorithms are still required to answer these questions.

The old concept of the RCP state incorrectly did not account for the jamming category of the packing. Previous attempts to estimate the packing fraction of the "random loose" state¹⁸ are even more problematic, given that this term is even less well-defined than the RCP state. Furthermore, as our investigations of disk packings show, the "stability" of packings cannot be judged based solely on local criteria, as suggested in Ref. 20 for sphere packings, and using such local criteria in estimating mean coordination numbers or densities of packings^18,19 is at best an exercise in modeling locally jammed packings. The best way to categorize random disk packings is to determine the maximally random jammed (MRJ) state¹⁰ for each of the three jamming categories (local, collective, and strict). Such investigations will be carried out in the future, and we have some preliminary results and promising avenues of approach.

The identification of the MRJ state for strictly jammed disk packings is an intriguing open problem. On the one hand, we have shown that random packings exist with densities in the vicinity of the maximum possible value (= /(2)) that are not strictly jammed, and on the other hand, there is a conjectured achievable lower bound /8 corresponding to the "reinforced" Kagomé lattice (see Fig. 4). It may therefore be that the search for the MRJ state for strictly jammed disk packings should focus on randomly diluted triangular packings. For random sphere packings, an initial study undertaken in Ref. 26, using the LP algorithm described in this work, found that maximally disordered random packings around 0.63 were strictly jammed, suggesting a close relation between the conventionally accepted RCP state and the MRJ state for strictly jammed packings. Much less obvious is what the MRJ state for collectively jammed sphere packings is. Finally, a completely unexplored question concerns the identification of the MRJ state for locally jammed disk and sphere packings.

The jamming concepts and algorithms presented here can be extended to packings of nonspherical particles with certain nontrivial modifications, however, mathematical developments in this area are lacking. We are investigating such extensions and will report our findings in future work. Other important tasks include extending various packing generation algorithms to generate strictly jammed packings, as well as designing algorithms with guarantees of producing jammed packings. An even more challenging task is designing packing algorithms that can produce jammed packings with certain target properties, such as a certain density and degree of order. The algorithms to test for jamming, and more generally to explore the set of reachable configurations _R for hard-particle packings can be further improved. In particular, a carefully tuned implementation of linear solvers for three-dimensional packings is needed as a building block in implementations of various nonlinear programming algorithms related to packings. Development of such implementations and extensions is also under way. Finally, there are many open questions related to the enumeration and classification of random hard-particle packings that might be answered with the application of these tools.

In Ref. 13 and this work we have proposed, implemented, and tested a practical algorithm for verifying jamming categories in finite sphere packings based on linear programming. We demonstrated its simplicity and utility, and presented some representative results for ordered lattices and random packings. Interestingly, the large computer-generated monodisperse random packings that we tested were virtually strictly jammed in three dimensions, but not in two dimensions. Future extensions and applications of the proposed algorithms are awaiting exploration. Work is already under way to provide highly efficient implementations of various optimization algorithms for linear and nonlinear programming on large-scale (contact) networks.

We would like to thank Andrea Liu and Corey O'Hern for providing us with sample packings and helpful e-mail discussions. A.D. and S.T. were supported by the Petroleum Research Fund as administered by the American Chemical Society and by the MRSEC Grant at Princeton University, NSF DMR-0213706. R.C. was partially supported by NSF Grant No. DMS-0209595.

1. Granular Matter, edited by A. Mehta (Springer, New York, 1994). first citation in article
   
2. S. F. Edwards and D. V. Grinev, Chem. Eng. Sci. 56, 5451 (2001). [ISI] first citation in article
   
3. R. Zallen, The Physics of Amorphous Solids (Wiley, New York, 1983). first citation in article
   
4. J. P. Hansen and I. R. McDonald, Theory of Simple Liquids (Academic, New York, 1986). first citation in article
   
5. S. Torquato, Random Heterogeneous Materials: Microstructure and Macroscopic Properties (Springer, New York, 2002). first citation in article
   
6. T. Aste and D. Weaire, The Pursuit of Perfect Packing (IOP, Bristol, 2000). first citation in article
   
7. S. Torquato and F. H. Stillinger, J. Phys. Chem. B 105, 11849 (2001). [INSPEC] [ISI] first citation in article
   
8. R. Connelly, K. Bezdek, and A. Bezdek, Discrete Comput. Geom. 20, 111 (1998). [INSPEC] first citation in article
   
9. Rigidity Theory and Applications, edited by M. F. Thorpe and P. M. Duxbury, Fundamental Materials Research (Kluwer/Plenum, Dordrecht, 1999). first citation in article
   
10. S. Torquato, T. M. Truskett, and P. G. Debenedetti, Phys. Rev. Lett. 84, 2064 (2000). [ISI] [MEDLINE] first citation in article
    
11. F. H. Stillinger, H. Sakai, and S. Torquato, Phys. Rev. E 67, 031107 (2003). first citation in article
    
12. L. F. Toth, Regular Figures (Pergamon, New York, 1964). first citation in article
    
13. A. Donev, S. Torquato, F. H. Stillinger, and R. Connelly J. Comput. Phys. (in press). first citation in article
    
14. B. D. Lubachevsky and F. H. Stillinger, J. Stat. Phys. 60, 561 (1990); [INSPEC] [ISI]
    B. D. Lubachevsky, F. H. Stillinger, and E. N. Pinson, J. Stat. Phys. 64, 501 (1991), second part of Ref. 14. [INSPEC] [ISI] first citation in article
    
15. C. S. O'Hern, S. A. Langer, A. J. Liu, and S. R. Nagel, Phys. Rev. Lett. 88, 075507 (2002). first citation in article
    
16. A. Zinchenko, J. Comput. Phys. 114, 298 (1994). [INSPEC] [ISI] first citation in article
    
17. J. G. Berryman, Phys. Rev. A 27, 1053 (1983), see also references con-tained therein. first citation in article
    
18. W. Uhler and R. Schilling, J. Phys. C 18, L979 (1985). [INSPEC] first citation in article
    
19. D. E. G. Williams, J. Phys. C 18, L181 (1985). [INSPEC] first citation in article
    
20. E. A. J. F. Peters, M. Kollmann, T. M. A. O. M. Barenbrug, and A. P. Philipse, Phys. Rev. E 63, 021404 (2001). [MEDLINE] first citation in article
    
21. A. Donev and S. Torquato, J. Mech. Phys. Solids 51, 1459 (2003) (ScienceDirect). [INSPEC] first citation in article
    
22. S. Torquato, A. Donev, and F. H. Stillinger, Int. J. Solids Struct. 40, 7143 (2003). [ISI] first citation in article
    
23. A. Donev, http://atom.princeton.edu/donev/Packing, homepage for the sphere packing project, with useful supplementary materials. first citation in article
    
24. M. Parinello and A. Rahman, J. Appl. Phys. 52, 7182 (1981). first citation in article
    
25. WWW, the general purpose interior-point LOQO optimization library is not public-domain, but can be tried at http://orfe.princeton.edu/ loqo/. The public-domain PCx library implements interior point linear programming algorithms and can be found at http://www-fp.mcs.anl.gov/otc/Tools/PCx/. first citation in article
    
26. A. R. Kansal, S. Torquato, and F. H. Stillinger, Phys. Rev. E 66, 041109 (2002). first citation in article
    
27. T. J. Quickenden and G. K. Tan, J. Colloid Interface Sci. 48, 382 (1974). [INSPEC] [ISI] first citation in article
    

Full figure

Fig. 1. A highly schematic plot of the jammed subspace in the density-disorder plane. First citation in article

Full figure

Fig. 2. Simple collective mechanisms in the Kagomé and honeycomb lattices, respectively. These lattices are not collectively jammed with periodic boundary conditions, as the sample unjamming motions for the Kagomé (left) and for the honeycomb (right) packings shown here illustrate. The shaded disks represent periodic images. First citation in article

Full figure

Fig. 3. Shearing the honeycomb lattice. The honeycomb lattice is not strictly (or collectively) jammed, and an example of a lattice deformation, replicated on several unit cells to illustrate the shear character of the strain = ()^–1. Note that only three (original) spheres are involved in the actual calculation of this unjamming motion, the rest are image spheres. First citation in article

Full figure

Fig. 4. Examples of strictly jammed lattices in two dimensions (from Ref. 13). The 6/7th lattice (last packing in Ref. 11), top, is obtained by removing every seventh disk from the triangular lattice. The reinforced Kagomé lattice, middle, is obtained by adding an extra "row" and "column" of disks to the Kagomé lattice and thus has the same density in the thermodynamic limit, namely, it has every fourth disk removed from the triangular packing (see also Ref. 11). It can be proven that this is the lowest density strictly jammed subpacking of the triangular lattice. The pentagonal packing shown at the bottom with 10 disks in the unit cell is obtained from a particular tiling of the plane with three rotated congruent pentagons, and is just one member of a whole family of strictly jammed packings. First citation in article

Full figure

Fig. 5. Virtually strictly jammed sphere packing. This random packing of 500 spheres with density = 0.64 was produced using the (original) Lubachevsky–Stillinger algorithm and it is collectively jammed and practically strictly jammed. The (cubical) unit cell is also shown. First citation in article

Full figure

Fig. 6. Collectively jammed bidisperse disk packing. The algorithm to test for collective jamming in ideal packings was applied to this equimolar bidisperse disk packing of 250 disks (= 0.846) in order to identify a jammed subpacking of 232 disks, leaving 18 rattlers (colored black), which are not essential for jamming. The dotted disks represent periodic images. Note that the density would be significantly lowered if the rattling particles were removed. First citation in article

Full figure

Fig. 7. Locally jammed disk packing. A random packing (= 0.82) of 1000 disks that is not collectively jammed, and a representative periodic unjamming motion. More insightful animations can be found on our webpage (Ref. 23). First citation in article

Full figure

Fig. 8. Collectively jammed disk packing. A dense (= 0.89) random packing of 1000 disks that is collectively jammed but not strictly jammed, and a representative unjamming motion. One can see the grains gliding over each grain boundary due to the shear, bringing this packing closer to a triangular lattice. First citation in article

Full figure

Fig. 9. Strictly jammed disk packing. We show here 2×2 unit cells of a dense (= 0.88) random packing of 250 disks that is strictly jammed, modulo four rattling particles, shown in black. This packing was produced with the extended Lubachevsky–Stillinger algorithm which allows for deformations of the lattice during the compression. We also display the contact network of the packing. The striking feature of this and similar strictly jammed disk packings we have produced is the appearance of peculiar "dislocation cores" and the appearance of large perfectly triangular regions. First citation in article

Table I. Classification of some simple lattices into jamming categories for periodic boundary conditions. We give the packing (i.e., covering) fraction (to three decimal places), the coordination number Z, and the number of disks/spheres Ns per unit cell, as well as an assessment of whether the lattice is locally (L), collectively (C), or strictly (S) jammed (Y is jammed, N is not jammed). We chose the smallest unit cells for which an unjamming motion exists (illustrated on our webpage—Ref. 23), if there is one.
Lattice		Z	Ns	L	C	S
Honeycomb	0.605	3	4	Y	N	N
Kagomé	0.680	6	3	Y	N	N
Square	0.785	4	2	Y	N	N
Triangular	0.907	6	1	Y	Y	Y
Diamond	0.340	4	4	Y	N	N
Simple cubic	0.524	6	2	Y	N	N
Body-centered cubic	0.680	8	2	Y	N	N
Face-centered cubic	0.741	12	1	Y	Y	Y
Hexagonal close-packing	0.741	12	2	Y	Y	Y

First citation in article

Table II. Results for monodisperse sphere packings. The columns are as in Table III, and here we show the running times for both the testing for collective and strict jamming.
N		t (s) coll	t (s) strict	/D coll	/D strict
50	0.628	23	29	0.0012	0.12
100	0.644	53	76	0.00043	0.15
250	0.636	164	210	0.0021	0.031
500	0.641	480	597	0.0037	0.014
750	0.641	900	1017	0.0015	0.0035
1000	0.642	1822	1866	0.011	0.013

First citation in article

Table III. Results of the nonideal randomized LP algorithm for equimolar binary disk packings of diameter ratio 1.4. The first column shows the total number of particles N, the second the packing fraction, the third the running time for the AMPL model that tests for collective jamming, and the last two columns show the average particle displacement during collective (i.e., with a fixed lattice) and strict jamming (i.e., with a deforming lattice) testing. Notice that the displacements are significantly larger for the strict jamming test, especially for small packings.
N		t (s) coll	/Di coll	/Di strict
50	0.845	2.1	0.010	0.060
100	0.842	6.4	0.0034	0.011
250	0.846	21	0.0037	0.0053
500	0.847	72	0.0016	0.0067
750	0.849	88	0.0022	0.012
1000	0.849	130	0.0016	0.018
1500	0.848	247	0.0016	0.020
2500	0.849	248	0.0039	0.010

First citation in article

Table IV. Results for monodisperse disk packings. The columns are as in Table III. Notice the very large displacements during the test for strict jamming, even for large packings, as well as the high packing densities for larger packings.
N		t (s) coll	/D coll	/D strict
50	0.832	2.9	0.0022	0.39
100	0.863	8.9	5.4×10–8	0.18
250	0.886	21	0.0014	0.86
500	0.891	78	6.7×10–5	0.16
750	0.887	103	0.0040	0.26
1000	0.882	153	0.0017	0.23

First citation in article

Table V. The average particle displacement /D during the test for collective jamming is shown for a series of sphere packings produced by the (original) Lubachevsky–Stillinger algorithm. From top to bottom the packing size N increases, and from left to right the number of collisions per particle Ncoll (in thousands) increases (and thus the density also slowly increases). No special handling of rattlers was employed. It is easily observed that the packings uniformly become "more jammed" as the packing algorithm is run longer (though rattlers may continue to give a finite contribution to the observed displacements). Similar behavior is expected of any algorithm which in the limit of infinite numerical precision produces packings with a collectively jammed subpacking.
N/Ncoll(103)	1	5	10	25
50	0.041	0.015	0.0018	4.9×10–10
100	0.036	0.016	0.0011	0.000 14
250	0.050	0.023	0.0015	0.000 36
500	0.047	0.024	0.0028	0.001 4
750	0.046	0.019	0.0030	0.001 1
1000	0.052	0.020	0.0025	0.000 67

First citation in article

Table VI. The analog of Table V but for strict jamming. In this case it is seen that the average displacements do not converge uniformly toward zero, an indication that the packings do not have a strictly jammed ideal subpacking (similar results are observed for amorphous binary disk packings). However, the average displacements are quite small for large packings (this is even more pronounced for the binary disk packings).
N/Ncoll(103)	1	5	10	25
50	0.083	0.057	0.059	0.051
100	0.066	0.042	0.023	0.026
250	0.052	0.027	0.010	0.0097
500	0.056	0.024	0.012	0.010
750	0.048	0.027	0.014	0.014
1000	0.060	0.025	0.0040	0.0021

First citation in article

Table VII. Just as an illustration, shown are the results for a two-dimensional disk packing with 250 disks, corresponding to the results presented for amorphous sphere packings in Tables V and VI. It is seen that although the packing has an ideal collectively jammed subpacking, it is clearly far from being strictly jammed, as typical for monodisperse disk packings produced by the Lubachevsky–Stillinger packing algorithm with a fixed lattice.
Ncoll(103)	1	5	10	25
Collective	0.12	0.007	0.000 50	1.7×10–5
Strict	0.45	0.24	0.12	0.12

First citation in article

^aElectronic mail: torquato@electron.princeton.edu