aa r X i v : . [ s t a t . M E ] O c t Statistical Inference for Disordered Sphere Packings
Jeffrey D. PickaDepartment of Mathematics and Statistics, University of New Brunswick
Abstract
Sphere packings are essential to the development of physical models for powders, com-posite materials, and the atomic structure of the liquid state. There is a strong scientificneed to be able to assess the fit of packing models to data, but this is complicated by the lackof formal probabilistic models for packings. Without formal models, simulation algorithmsand collections of physical objects must be used as models. Identification of common as-pects of different realizations of the same packing process requires the use of new descriptivestatistics, many of which have yet to be developed. Model assessment will require the use oflarge samples of independent and identically distributed realizations, rather than the largesingle stationary realizations found in conventional spatial statistics. The development ofprocedures for model assessment will resemble the development of thermodynamic models,and will be based on much exploration and experimentation rather than on extensions ofestablished statistical methods.
A disordered sphere packing process can be loosely defined as the type of spatial stochasticprocess whose realizations describe the final positions of a collection of spherical objectswhich are tossed into a container and come to rest. Sphere packings are used in physicsand engineering to model the atomic structure of amorphous solids, the internal structure of1omposite materials, and the static and dynamic properties of powders. They are examplesof spatial stochastic processes which need to be fit to data in order to be properly used, yetno formal procedure has yet been developed for objectively assessing the quality of fittedmodels. These processes have a very different within-realization disorder to that found inprocesses normally studied in spatial statistics, and so very little study of these processes hasbeen undertaken by statisticians or probabilists. Most work on packings has been undertakenby physicists and engineers who tend to focus on the mean behaviour of the processes and toexclude consideration of their variability, although their work has highlighted a serious needfor fit assessment. By combining ideas and methods from spatial statistics and from physics,it is possible to develop a formal approach to inference which is useful in a scientific andengineering context. This approach to inference is also useful in the context of developingmethods for fitting spatial models to highly dependent spatial data which is not amenableto fitting by conventional statistical methods.It will be generally be assumed that the spheres in packings are of equal diameter( monodisperse ) rather than being of many different sizes ( polydisperse ). Spheres in R willbe termed discs, spheres in R will be termed spheres, and in higher dimensions they willbe termed hyperspheres. In the applied literature, planar packings of discs are sometimesreferred to as packings of rods are cylinders. This is misleading, since it assumes absoluterigidity for the physical rods and cylinders being modeled. A physical packing is a collec-tion of spherical objects assembled by tossing those objects into a container and allowinga force to form them into a jammed static arrangement. Packings will be assumed to bedisordered and not based on a point lattice unless otherwise stated. Ordered packings withpoint lattice structure are of use in the study of geometry in high-dimensional spaces (1; 2; 3)and in the development of efficient coding (4), but these applications have no relevance tothe modeling of disordered materials. Packings are also of great interest in more abstractmathematical contexts (e.g. the theory of discrete analytical functions (5)), but these uses2ave no applications to materials either. Sphere packings were first used in attempts to explain the patterns seen in X-ray diffractionstudies of liquids and glasses (6; 7). Bernal (8) and Scott (9) attempted to model macroscop-ically what an amorphous mass of densely arranged atoms might look like, using physicalpackings of monosized steel spheres as their model. More recent work (10) has suggestedthat this approach to modeling may not be feasible, since the packings of spheres are formedby the action of a different set of forces from those found at the atomic scale. If physicalmodels can be shown to explain the patterns seen, then this use of the physical model wouldbe analogous to the use of linear models in other sciences where the use of a simple modelcan be useful, even if it is partly incorrect.Sphere packings are more commonly used in the modeling of granular materials at rest.Granular materials are collections of solid particles surrounded by fluid or gas which formpackings at rest, but which flow like liquids or gases when a sufficiently large force is applied.When put into motion, the complex interaction between the grains makes it impossible tomodel granular flow by means of conventional methods from statistical physics (11; 12; 13).Packings are essential components of models for the static and dynamic behaviour of granularmaterials which are used the study of avalanches, of mudslides, of soil liquification anderosion, and of the flow of powders through pipes and hoppers.Sphere packings can be used in the study of composite materials. In the study of metalsintering (14), a packing can represent the initial state of a metal powder before it is com-pressed and heated to form a solid metal part. Polydisperse packings can be used as modelsfor some some colloids and for concrete. A colloid is a mixture of two immiscible liquids, inwhich one liquid may tend to form spherical inclusions within the other. For some compo-sitions, the spherical inclusions in a colloid may appear to form a packing (15). Concrete3s a composite of rocks of many different sizes held together by a matrix of cement (16).Sphere-packing-based models for concrete are highly idealized, but may be the only feasibleway to abstractly model such a complex material (17).Packings can be used as models for porous structures. A packed bed reactor is a largesteel vessel filled with solid catalyst-impregnated particles. Liquid reactants flow in at oneend of the vessel, and the product emerges from the opposite end. Design of these reactorsrequires being able to model both the reaction in the pore spaces around the particles andthe transfer of heat through the particles. Sphere packings are the simplest form of particlepacking that can be used to model the internal structure of these reactors (18).
A sphere packing process is a stochastic process whose realizations are arrangements ofspheres which closely resemble the arrangements of spherical objects or atoms produced byphysical processes in nature. Ideally, it would be possible to define the process mathemat-ically and then to analyze packing processes through mathematical arguments, as can bedone for Poisson processes and spatial Markov processes. Insurmountable difficulties posedin defining the packed state and in defining the process make it impossible to undertake ananalytical approach to inference.The first problem is to define what it means for spheres in a realization to be packed.Consider a configuration of spheres in R d . Construct the Delaunay triangulation of the spherecentres (19). An interior sphere is any sphere whose centre has no triangulation edges linkedto any of the synthetic points used to construct the triangulation. Boundary spheres arespheres which are not interior spheres. Identify all points of contact between spheres and theirneighbours, and define the contact network to be the subset of the triangulation which joinscentres of spheres in contact. If confining surfaces are present, then the contact network must4e augmented by edges connecting points of confining surface contact with the correspondingsphere centres. A sphere is jammed if it cannot be translated without overlapping anothersphere. A boundary sphere is jammed if it cannot be translated without overlapping anothersphere or penetrating a confining surface. Otherwise, it is a free boundary sphere .Given an configuration of spheres, then the simplest definition of a packing would requirethat all interior spheres be jammed and that the contact network consist of a single connectedcomponent. This definition is too strict, in that it is possible to have boundary spheres whichare not jammed (called rattlers ) inside of cages of jammed spheres which may also includeconfining surfaces. The simple definition is also not strict enough in many applications, inthat it may include arrangements of spheres which are unstable when the packing is underthe influence of forces. If forces are acting on the packing, then all of the spheres includingthose on the free boundary need to be in mechanically stable positions.Given these definitions, packings of monosized spheres can only be shown to exist bymathematical arguments if the spheres are arranged in a regular lattice structure, suchas a hexagonal or cubic close packing. There is no way to demonstrate mathematicallythe existence of a disordered packing of monosized spheres. Spherical physical objects canbe packed by mechanical processes, but the physical objects involved are neither exactlymonosized or exactly spherical. Given any list of sphere centre locations in a possible packing,the exact distance between the points can never be exactly known, and so it cannot be saidfor certain if the list represents a disordered packing of monosized spheres. It is possible thatthere is no such thing as a disordered packing of monosized spheres, which would complicatethe formal definition of a packing process. The same difficulties with existence apply topolydisperse packings, with some exceptions related to circle packings in the plane (5).If disordered packings could be shown to exist, then it would be necessary to formallydefine a stochastic process which had disordered packings as its realizations. One candidatefor such a model is to use the Gibbs process (20) to model the sphere centre positions, which5ould require that the packing be of infinite extent and be stationary. Gibbs processes areonly useful when the only interactions between spheres are pairwise (21). In a Gibbs modelfor a packing, there would be many higher-order interactions which could neither be ignorednor simply approximated.Instead of formally defining a Gibbs model for a packing process, it may be possibleto approximate it by explicitly listing the possible configurations for the spheres and thenassigning a probability that each configuration will appear in a realization of a packingprocess. The construction of such configurations would require arguments similar to thoseused in the identification of the largest fraction of a square (or equilateral triangle or regularhexagon) which can be covered by non-overlapping monosized discs. These calculationscan be done analytically for 9 or fewer discs, and computer-based proofs are needed for 10to 27 spheres (22). For larger configurations, approximate simulation-based arguments arerequired (23; 24) which makes this approach futile for any reasonably-sized packing.The lack of a useful formally defined model such as the Poisson process or the Gibbsprocess implies that likelihood-based inference is impossible for any physical packing process.Any inference or model assessment procedure must be based on replications of physicalprocesses or on independent realizations generated by computer programs.
Since mechanically stable packings of physical objects exist in nature, it is scientificallysensible to attempt to represent them by some type of abstract model which is simplerto deal with than the physical objects themselves. While formal mathematical modeling isimpossible, other physical systems and complex computer programs can be used as stochasticmodels instead. 6 .1 Physical Models for Packings
A physical model for a packing of objects is a different physical system which is easier toreplicate and study than the original packing. Examples include using packings of steelspheres to represent the positions of atoms in amorphous solids or to represent the internalstructure of catalyst pellets in a packed bed reactor. In the latter case, a scale model of thereal system could be prepared out of different materials from the reactor, so that internalimages of packed bed could be found. Physical models are also useful for investigatingpackings of soft objects and objects of irregular shape (25) when these objects are toodifficult to simulate.The physical model is stochastic because its realizations are generated by a process thatcan be modeled by a dynamical system which is sensitive to initial and boundary conditions.Pouring spheres from one container to another is an example of this type of physical process,since repetitions of the process under the same general experimental conditions will producedifferent disordered packed configurations. For any such process, there is a sample space ofpossible packed object configurations and a density which represents the relative likelihoodof each configuration. The goal of model selection in the case of the physical model isto match the sample space and distribution to that of the physical system which is to bemodeled. Both of these sample spaces and densities are inaccessible to theoretical argumentand difficult to approximate, and so any inference must be based on comparisons of thereplications.The first major use of physical models involved attempts to model the molecular structureof liquids (6). Scott (9) and Bernal (7; 8) both used packings of monosized steel spheres asmodels for atomic structure. Packings were made inside of thick rubber balloons bound withrubber bands in order to maximize packing density and reduce the effects of gravity. Bernalalso filled the balloons with black paint, which was then allowed to set. The agglomeratedmass was dissected carefully in order to identify which spheres were in contact, and physical7odels of the contact network were constructed. More elaborate measurement methods weredeveloped (26), but interest in this approach soon waned on account of the amount of timeand effort required to obtain useful results.Physical packings are also distinguished by the material that fills the voids around thespheres. A physical process that packs spheres in room temperature air is distinct from onethat packs spheres surrounded by warm liquid wax or paint. To make low-density packings,spheres have been packed in a fluid of equivalent density (15). It also possible to makesomething very much like a packing by combining two immiscible fluids, one of which formsnear-spherical inclusions in the other (27).Physical packing processes for discs can also be developed. These could be as simple asswirling a tea tray full of thick coins and then inclining the tray until the coins come torest. Processes of this kind produce realizations which are very close to hexagonal packings,and so more complicated methods are needed to obtain more disordered realizations. Non-overlapping discs can be placed at random locations on a slightly stretched sheet of rubber.The tension in the sheet is released, forcing the discs together into a smaller area. Thepositions of the discs are recorded, the sheet is retensioned, the discs are placed back on thesheet, and the tension is released again. By cycling through this procedure, a disc packingcan be generated (28).Interest in physical packings has been revived as new remote sensing methods have beendeveloped. X-ray tomography can be used to map the interiors of very large packings (29;30; 31; 32). Magnetic resonance imaging (33) may be less expensive, but it also produceslower-resolution images and cannot be used on very large packings. Confocal microscopy hasbeen used to study both sphere packings (34; 35) and colloidal packings (27).8 .1.1 What can be learned from physical packings
Any model of a physical process must obey the constraints on packed configurations whichhave been revealed by the study of packings of physical objects.The choice of experimental procedure which defines a physical packing process has agreat impact on the sample space of configurations associated with that procedure. Experi-ments with monosized spheres poured into containers showed that the density of the initialconfiguration of the pour could be made more dense by rodding the packing or by vibratingit up and down (9). There appeared to be a well-defined limit to the density achievable bythese methods, which resulted in the proposition that there was a well-defined physical statetermed a dense random packing of spheres. This state has been shown to not be well-defined(36), and it is possible to use mixing procedures in a Couette-Taylor cell (37) and othermethods (38) to produce denser near-ordered structure packed configurations. The implica-tion for modeling is that an arbitrary packing model cannot be assumed to be sampling fromthe same sample space of configurations as does the physical process that is being modeled.The apparent existence of a well-defined dense random packing opened the question asto whether or not a disordered packing could be denser than a hexagonal close packingof spheres. The mean volume fraction of dense random sphere packings was thought tobe between 0.63 and 0.64, but that did not rule out the existence of some configurationswhich might have higher volume fractions. For dimensions 2 ≤ d ≤
8, the densest latticepackings can be identified. For discs in the plane, the densest packing of any kind can beanalytically proven to be the hexagonal close packing (39). A computer-based proof with thesame conclusion has been constructed for spheres (40). In dimensions 10 and higher, somenon-lattice packings have been found which are denser than any known lattice packings (4).Polydisperse packings are far more disordered than monodisperse packings. One way ofdisrupting the near-order of physical packings of monosized discs is to randomly introducea small amount of size variation into the discs. Physical packings of a fixed sphere size9istribution can be very difficult to assemble, since the smaller spheres often segregate andaccumulate at the base of the specimen. These difficulties also occur in simulated polydis-perse sphere packings (41; 42).
A simulation model is a program which takes output from a random number generator andproduces an arrangement of spheres which is close to that of a packing. While these algo-rithms could in theory be represented by a formal probabilistic model, they are sufficientlycomplex that no such model could be stated or usefully analyzed.
Any simulated monodisperse sphere packing will be defined by a list of sphere centre loca-tions. With the exception of some cubic packings, most of these numbers will be irrational.In any simulation, sphere location coordinates can only be specified to a finite number ofdecimal places. It is not possible to determine if the list of points produced by any simulationprogram form a packing or not. It is generally assumed that any simulated packed configu-ration can be collapsed to a rigid jammed configuration, but it is not possible to determineexactly which spheres in the collapsed configuration are in contact and which ones are not.Examination of physical packings shows that they contain many pairs of spheres which areclose but not in contact (8; 43). If the simulation is intended to study a physical phenomenonwhich is very sensitive to the structure of the contact network, then this uncertainly couldcontribute additional uncertainty to any physical properties estimated from the simulatedpacking.If round-off errors are the only type of error present in the simulation, they can be min-imized through use of maximum precision variables in the code. It is also possible to adaptmethods used to deal with errors arising from detectors and from Fourier reconstruction in10-ray tomographic imaging of physical packings. In those reconstructions, the estimateddistance between the closest neighbours was found to be approximately Gaussian (29). Thisdistribution was used to estimate the mean number of neighbours in contact, and could beused to define a stochastic or deterministic method for deciding if two spheres are in contact.
Most physical packings are formed by a confining surface. This may be a plane orthogonalto the direction of the force driving the physical flow that produced the packing, or it maybe a completely closed surface that jams every sphere within the packing. Simulations oftenomit confining surfaces, either to simulate small parts of much larger packings or to avoiddifficulties in coding.It is possible to avoid using any boundaries. Spheres can be packed on the surface ofa hypersphere of one dimension higher (44), but this can introduce unwanted effects fromhypersphere surface curvature unless many spheres are used. Packings can also be assembledby simulating a force acting towards a central point. It is not clear that a packing simulatedby a central force would represent any packing found in nature.Confining surfaces can be avoided by imposing periodic boundary conditions . Considera packing in the plane and a rectangle W . If the discs are monosized with diameter δ ,then any disc whose centre lies within δ/ W has the part of the disc outsideof W appear on the opposite side of the rectangle. Since this can occur on any edge, thewindow W is toroidal and the plane is tiled by copies of W containing a disc packing. If theboundary effects vanish after some distance kδ from the periodic boundary, then the interiorof a realization may be indistinguishable from that of a stationary and isotropic packing.Efforts have been made to determine a useful value for k (45; 46), but these are often basedon a single statistic and may not reflect every effect of periodic boundaries on the structureof a realization. 11eriodic boundary conditions may be used within confining boundaries in order to savecomputational effort. Within cylindrical confining boundaries, it is possible to pack fourperiodic cells (47). In any small system, this symmetrical structure may have a strong effecton the distribution of statistics calculated from the packing. Most simulation models pack monodisperse spheres, since the code for these programs isrelatively easy to write. It is dangerous to assume that such models can represent packingsof more complicated objects, unless statistical methods can be used to show that simplermodels can represent the physical phenomenon of interest.Modeling packings of randomly shaped objects is very difficult. While a sphere requiresonly its centre and radius to completely specify it, a reasonable simulation of a small rockmight require a long list of numbers to specify the polynomial surfaces that are splicedtogether to produce a simplified description. When monodisperse spheres are packed, feasiblearrangements can be determined by comparing between-centre distances to the diameter ofthe spheres. For random shapes, very complex code is needed in order to efficiently determineif any overlaps have occurred. Until these problems can be overcome, the best results forrandom shapes will come from studying physical packings rather than simulated ones.It is possible to simulate the packings of simple non-spherical objects. Ellipses can bepacked in the plane and ellipsoids can packed in space because these shapes can be easilydescribed and overlaps can be easily determined. Ellipsoids can be shown to pack moredensely than spheres in simulations (48; 49; 50) and in physical experiments (51). Theshapes of ceramic rings have been approximated by assemblies of triangles (52) as part of asimulation of packed beds. 12 .3 Ballistic Algorithms
Ballistic algorithms are the simplest algorithms which can be used to produce configurationsof packed spheres or discs. They were the first developed and are relatively simple to code.The algorithms have no basis in physics, and packings generated by these algorithms areless dense than those generated by physical packing processes (53). Ballistic models mustbe used with caution in physical modeling, as there is no reason to believe that they havethe same sample spaces as do physical processes.The first ballistic algorithms were developed by Vold (54; 55) to model the formationof flocs. Spheres are dropped at random locations onto a surface. They stop falling eitherwhen they hit that surface or another sphere. If they hit another sphere, then they either arelocked in place with probability p , or are rolled down the packed spheres until they reach astable position. With p = 0, the algorithm can produce realizations of loose packings. Voldalso developed a central version of this algorithm (56) and one involving rods composed ofa set of k spheres in contact with their centres arranged along a line (57).The first high-density central packing simulators built up packings by placing spheres atthe closest possible position to the centre of a packing (58; 59). These algorithms producepackings with estimated volume fractions between 0.61 and 0.62. For central disc-packingalgorithms, near-ordered patterns can be avoided by seeding the realization with a non-triangular configuration of discs (60).Visscher and Bolsterli (61) developed a non-central algorithm which added periodicboundary conditions perpendicular to the base. A sphere is dropped at a random loca-tion, and falls vertically until it contacts another sphere. Then, it follows the shortest pathalong the surface of the already-placed spheres until it comes to rest in a gravitationallystable position. Each sphere drop is repeated k times, and the final position chosen is theposition that is closest to the base. The algorithm produces realizations with estimatedvolume fraction of 0.582 in R , which is less dense than a settled physical packing. In the13lanar version of this algorithm, the initial spheres dropped have a slightly different size inorder to prevent the formation of near-ordered realizations.The Visscher-Bolsterli algorithm has been extended in various ways. A central versionhas been developed (62), with the intention of producing packings that contain bridgingstructures. Non-periodic boundaries have been added to the algorithm, so that the spheresfill a box in space (63) or in the plane (64), or fill a cylinder (65; 66; 67). Using periodicboundaries, irregular shapes can be packed (52; 68).The Visscher-Bolsterli algorithm can be modified to pack spheres of different sizes. Whenthe largest and smallest spheres are greatly different in size, the smaller spheres tend toaccumulate at the base of the packing unless they are allowed to randomly stick at unstablepositions (41). More complicated deposition rules have been devised to avoid size segregationin the plane (69; 70).The density of physical packings can be increased by shaking them. The density ofballistic packings can be increased by subjecting the initial packing to rearrangements whichsuperficially imitate shaking. A simple rearrangement algorithm orders the initially packedspheres by height, then redeposits them at the same planar location using the Visscher-Bolsterli placement rules (71; 72). This algorithm was used to study segregation by sizeafter shaking. In a shaking algorithm for monosized spheres, each sphere in the packing isdisplaced upwards by a small Normally distributed perturbation, and is then subjected tomany small three-dimensional Normal perturbations which are allowed if no collisions occur.Once the number of collisions reaches a set threshold, the packing is collapsed from thebottom using VB deposition rules. This algorithm can increase the estimated mean volumefraction from 0.581 to 0.590, and still results in a loose packing (73). If the initial verticalperturbations are not sufficiently small, then the rearrangement can simplify the contactnetwork and decrease the volume fraction. Rearrangement methods can also be applied topackings of irregularly shaped objects in a container with hard boundaries (74).14 .4 Rearrangement Algorithms Algorithms which rearrange the points in point patterns were the first to achieve volumefractions similar to those seen in physical packings. These programs begin with a realizationof a Poisson or a regular point pattern. The points in the pattern are subjected to determin-istic or random translations which eventually produce a sphere configuration close to thatof a packing. These methods do not attempt to replicate the dynamic interactions betweenspheres which occur during the formation of a physical packing.The Jodrey-Tory algorithm (75; 76; 77) was the first simulation algorithm to producerealizations with estimated volume fractions similar to those for dense physical packings.It is initialized by a configuration of n points in the interior of a rectangular prism withperiodic boundaries. A sphere of unit radius is attached to each point. In the first stage ofthe process, the radius of each sphere shrinks by 0.0001 units. The distance between eachsphere and its nearest neighbour is found, and the closest pair defining overlapping spheresis identified. This pair are moved apart along a line joining the points until the spheres nolonger overlap. When the distance between the closest two spheres drops below a threshold,all spheres have their diameter increased by 0.0002 units, and the process repeats. After 2000cycles, a second routine of shrinking and translation removes all remaining overlaps. Oncethe initial configuration of points has been chosen, the algorithm is entirely deterministic.There are also versions of the algorithm which pack ellipsoids (78) and discs of two differentsizes (79).The Jodrey-Tory algorithm can be initialized with one of its own previous outcomes. Ifthis done many times, then the configuration crystallizes to become denser and less disordered(80). This phenomena is not seen in the simplest physical processes for generating densephysical packings, and it suggests that the Jodrey-Tory algorithm is not sampling from theset all feasible packed configurations in the same way that a physical packing process does.Other deterministic rearrangement procedures use more complicated rearrangement rules15hich may depend on all neighbours (81) or on fixing the contact network early in therearrangement (82). The latter strategy can be used to maximize disorder in planar discpackings (83).Rearrangement algorithms can also be inspired by models for the motion of gas molecules,but these models do not represent the physics of the condensation of a gas. Lubachevskyand Stillinger (84; 85) begin with points uniformly distributed within a region. Randomvelocities are assigned to each point, and each point begins to grow a disc at a constantrate. When discs collide, the collision is elastic and momentum and energy are conserved.After a few thousand iterations of disc expansion, the discs form a near-packing. There isa three-dimensional version of this algorithm (86) and a version which can pack ellipsoids(87).Rearrangement rules can also be based on the minimization of potential functions. Onecan begin with a uniformly distributed set of sphere centres and then rearrange these centresvia a conjugate gradient algorithm which minimizes a potential (88). The potential canbe one used in a thermodynamic model of a packing which is stable until a yield stress isexceeded.Rearrangement algorithms can also be written with non-periodic boundaries. They canpack spheres on the surface of a large sphere in R , avoiding any boundary effects (44).They can be used as the basis of programs to estimate the most efficient packing of a smallhard-boundaried region by a fixed number of spheres (23). When combined with elementsof ballistic algorithms, they can be used to pack cylinders (89). Ballistic and rearrangement algorithms have no basis in the physics of packing formation.Physical packings are generally formed by pouring physical objects into a container or intoa pile. In flow, the objects may be nearly packed or barely interacting. When objects strike16onfining surfaces or other already-packed objects, they expend energy through erosion,through fracturing, through friction, and through deforming themselves and the container.Exactly how this energy expenditure happens during physical packing formation is unknown,as these processes cannot be observed.It is possible to construct models for the formation of packings which incorporate idealizedversions of the unobservable energy expenditure processes. These models, known as DiscreteElement Method (DEM) models, were originally developed for modeling rock mechanics (90).DEM models have been proposed for both piling and container-filling processes for bothcohesive and cohesionless particles (91). They can be adapted to include processes whichhappen in viscous fluids and in environments where both viscous fluids and air are present(92). Since they are all based on unverifiable assumptions about interparticle behaviour, itcannot be assumed that any of them will capture the physics of packing formation.The first DEM packing model was that of Yen and Chaki (93), which simulated theformation of a packing in a rigid-walled box. The spheres are initially positioned above thebox in a configuration produced by a simple sequential inhibition process (94). As timeevolves, gravity acts upon each particle. When particles come into contact, they experiencea Herzian force associated with particle deformation and a tangential frictional force. A VanDer Waals force can be introduced to provide a cohesive force between the spheres. Whenfrictional forces are omitted, the algorithm produces realizations of estimated mean volumefraction 0.633 with higher volume fraction variance than expected. When frictional effectsare added, the mean volume fraction falls to 0.578.Since the initial model of Yen and Chaki, many improvements have been made to DEMmodels (92; 91). Rotational interactions have been added (95), as has the capacity forpacking polydisperse spheres (96) and complex shapes (97).17
Descriptive Statistics
Any one realization of a packing from a packing process consists of a long list of object loca-tions, each of which may be accompanied by size, shape, and orientation information. To fita model to a physical packing process, it is necessary to have statistics which can summa-rize what each realization from the same model have in common, and which can distinguishbetween realizations from different models. In physical and engineering applications, it isnecessary to have statistics which can characterize the physical response of interest and alsoto have predictive statistics which explain how the response depends upon packing structure.These statistics may be spatial averages of various kinds, but may also may be measures ofvariability of extreme behaviour.In physical examples, the search for descriptive statistics is analogous to the processwhich developed the original thermodynamic state variables used in describing ideal gases.Unlike the statistics in classical thermodynamics, descriptive statistics for packings must bethought of as random variables having distributional properties beyond their means whichmust be modeled. In the physics literature, it is customary to estimate a statistic from asingle large sample, and then to equate the estimate to its expected value via an implicitappeal to a law of large numbers. This may be reasonable in some contexts, but for anyone statistic there is no simple way to establish how large the packing must be before thebetween-realization variance of that statistic becomes negligible.To be useful for model fitting, descriptive statistics must be able to distinguish betweenrealizations from different models in spite of between-realization variability. Statistics forspatial processes often lack the power to make such distinctions, and so many statistics mayhave to be tried before a collection can be found that are sufficiently powerful. To findthese statistics requires having a library of statistics which can describe packings in as manydifferent ways as can be conceived of. Statistics are needed to describe what the expert eye18an see, and also what no human eye can identify. One of the most significant challenges indeveloping inference for packings is to develop statistics which can be used in these libraries.
Sphere packing processes are examples of random closed sets (98; 99). The random set Φ isformally defined by the collection of random indicator functions I Φ ( x ) = 1 if x ∈ Φ= 0 if x Φfor all x ∈ R d . Taking the expected value of this indicator function and its products allowsmoments M p to be defined at every point { x , . . . , x p } ∈ R d . M p ( x , . . . , x p ) = E [ I Φ ( x ) . . . I Φ ( x p )]= P r [ x ∈ Φ and x ∈ Φ and . . . x p ∈ Φ]These probabilities are defined on ensembles of many realizations, and do not necessarilydescribe any particular aspects the internal structure of a single realization. The k th momentis known as the k − point correlation function or the k − point probability function in statisticalphysics.If the random set is stationary, then reduced moments can be defined: m p ( x , . . . x p − ) = E [ I Φ (0) I Φ ( x − . . . I Φ ( x p − )]= P r [0 ∈ Φ and x ∈ Φ and . . . x p − ∈ Φ]The reduced first moment m is referred to as the volume fraction and m ( r ) as the covari-19nce.If Φ is also ergodic, then the reduced moments can be estimated from a single realization: c m = | φ ∩ A || A | c m ( x ) = | φ ∩ A ∩ φ x ∩ A x || A ∩ A x | where A is the window of observation, | A | is the volume of A , and A t = { x ∈ R d : x + t ∈ A } .These estimators are often themselves estimated by sampling at a regular or random arrange-ment of points within A . When the random set is non-stationary, c m estimates1 | A | Z A M ( t ) dt, whose interpretability depends upon the spatial variability of M ( t ).Both the packing and the closure of its complement can be considered to be random sets,so long as the boundary between them is well-behaved and not some sort of fractal. This isthe case for a sphere packing, and so moments m ck of the realization of the complement ofthe random set can also be defined. Mixed moments can also be defined: m ( x, y ) = P r [0 , x ∈ Φ , y ∈ Φ c ] = m ( x ) − m ( x, y ) . While mixed and complementary moment functions can be calculated from the moments ofthe random set, they may be of use if they can illustrate differences between processes moreaccurately than the moments can.The most commonly used descriptive statistic for packings is c m . When the randomclosed packing was believed to be a well-defined physical state, c m showed that the earliestsimulation algorithms were simulating looser packings. Later simulation algorithms werejudged by their ability to attain estimated volume fractions close to those of a dense random20acking (76). This statistic cannot be used alone for characterizing packing processes, sinceit reveals nothing about interactions between spheres. Local estimates of c m can be used todescribe differences in structure within non-stationary sphere packings (100; 52).Random set second moments m ( r ) are rarely estimated. Instead, point process statisticsare used to describe sphere interactions. For monodisperse sphere packings, there is a veryclose relationship between m ( r ) and the pair correlation function for the sphere centres(101).Moment statistics can also be applied to transformed random sets. If the spheres growat a constant rate until they fill space, estimates of m and the Euler-Poincar`e coefficientcan be plotted as a function of the degree of expansion (102; 103).The complement of the packing is its void structure. This structure can be describedby the spherical contact distribution function S ( r ), which is also known as the pore sizedistribution function (73; 82). Given an arbitrary point x in the void, S ( r ) is the probabilitythat the nearest point on a sphere to x lies within a distance r of x . It can be estimatedby finding the distance to the nearest point on a sphere from many locations which havebeen randomly or systematically sampled from within the void. It has been used to comparephysical packings with simulated packings and to investigate the applicability of theoreticalapproximations for S ( r ) which arise from statistical mechanics (104). The centres of each sphere in a packing constitute a realization of a point process. Thebasic point process statistics are based on the analogues of the first and second moment fora random measure, and on the Palm distribution of the point process (105).The intensity λ , defined for a stationary process to be the mean number of points perunit volume, is seldom estimated since λ is a constant multiple of m . Local estimates ofintensity can be used to quantify disorder (36) and to investigate the internal structure oflarge physical packings (29). 21he K − function K ( r ) is the second reduced moment function of a stationary ergodicpoint process. The quantity λK ( r ) is the expected number of points in a disc of radius r about a randomly chosen point, which is not counted in the expectation. The K − function isseldom used in the study of packings, and the pair correlation function g ( r ) is used instead.It is defined by g ( r ) = 1 db d r d − ∂K ( r ) ∂r , where b d is the dimension of the unit ball in R d . The pair correlation function may be alsobe referred to as the radial distribution function, although that name is also applied to thequantity RDF ( r ) = λdb d r d − g ( r ). The K -function has also been generalized to a functiondefined on the number of r − close triples (106) in a realization. Estimation of the K -functionfor stationary ergodic point processes requires careful treatment of window edge effects (107).If the point process is not stationary, estimators of g ( r ) for a stationary ergodic processare difficult to interpret. Estimates are made by taking each point on the realization andcreating a sequence of thick shells around it. For each shell, the number of other points inthe shell is counted. These counts are totaled for each shell type and rescaled based on theshell volume. A histogram of the rescaled totals is plotted. If the packed spheres have unitradius, a plot of this estimate has a maximum at r = 1, then has a pair of very sharp localmaxima at at r = 1 .
73 and r = 2 (43). This type of plot can be used to show the transitionfrom a packing to a looser structure in a DEM model as the Van Der Waals force increases(95). For planar packings, the estimated pair correlation function can be used to create apattern that would be seen if X-ray crystallographic methods were applied to a real materialwith the same structure as the realization (48). Those familiar with crystallography mayfind these patterns easier to interpret than estimates of g ( r ).Given an arbitrary location x ∈ R d , the empty space function H s ( r ) is the probabilitythat the nearest point in a realization to x lies within a sphere of radius r centered at x . For22 stationary point process, the nearest-neighbour function D ( r ) is the probability that thenearest neighbour to any point in the realization lies within a distance r of that point. Forthe Poisson process, these two probabilities are identical. For all other point processes, theycan be compared by the J-function (108), defined as J ( r ) = 1 − D ( r )1 − H s ( r ) . Nearest-neighbour functions can also be extended to the 2nd-nearest-neighbour, the 3rd-nearest neighbour, and so on (83). Estimation of all of these quantities require carefulattention to window edge effects.
The triangulation of the sphere centres in a packing forms a simple and natural descriptionof packing structure. Near neighbours can be clearly defined as spheres whose centres areconnected by a triangulation edge, and the contact network of the packing is a naturalsubgraph of the triangulation. The dual graph of the triangulation is a simplified descriptionof the void structure of the packing. Statistics based on triangulations and tessellations wereinitially applied to packings by physicists who were using physical packings as models forthe molecular structure of liquids and amorphous solids (8; 43).The Delaunay triangulation and the Voronoi tessellation are generally used as bases forstatistics. The tessellation is constructed by finding all points in R d which are equidistantbetween sphere centres. Euclidean distance is used in the construction, but other distancemeasures can be used to generate different tessellation structures (19). The triangulationis the dual of the tessellation, generated by joining pairs of sphere centres which define atessellation edge.The simplest statistics that can be extracted from triangulations and tessellations are lists23f characteristics for each cell. For both types of cells, the area, the perimeter, the largestand smallest angles, and the longest and shortest sides can be found. The coordinationnumber of each sphere is defined to be the number of triangulation edges which belong tothe contact network. The lists of observations are traditionally summarized by histograms orby summary statistics (mean, standard deviation, minimum, and maximum). In the physicsliterature, the histograms are often referred to as plots of the distribution of the statistic. Thisstatement is misleading, since the data is a collection of spatially dependent observations. Inthese cases, the histograms summarize aspects of the structure of an individual realizationrather than estimating a population parameter. For Voronoi cell areas, Gamma densitieshave been found to summarize the shape of the histograms (83; 109).Tessellation statistics have been used to investigate the differences between realizationsof different physical and simulated packings. Differences between packings are generallyspecified first in terms of differences in the estimated volume fraction of spheres c m . Plotscan be made of mean contact number versus c m (110; 29), or of the standard deviations ofVoronoi cell face area and volume against c m (111). Comparisons of histograms of contactnumbers from different algorithms are not powerful enough to clearly distinguish realizationsfrom different algorithms (53).Triangulation-based statistics can also be used to identify possible rattlers. In a physicalpacking, there are often a small number of spheres which are not immobilized. The packingsurrounding the rattlers consists entirely of jammed spheres forming arch structures, alsoknown as a bridging structures. These structures cause blockages when silos filled with bulksolids are emptied (112; 113; 114). The fraction of rattlers in a packing can be used as ameasure of packing efficiency (81). Identification of rattlers is complicated by uncertaintyabout whether or not neighbouring spheres are in contact.More elaborate statistics can be developed from the tessellation and the triangulation.The local density of a packing can be defined as the ratio of the sphere volume to the volume24f its Voronoi cell (29). The escape fraction statistic is the empirical cumulative distributionfunction of measurements of the largest sphere that could escape from each tessellation cellthrough the gaps between neighbouring spheres. It has been shown to distinguish betweenphysical packings of different volume fraction (29). A topological density has been definedwhich is based on a notion of topological distance. This distance is defined by choosinga sphere, and then identifying a sequence of shells radiating out from it. The first shellcontains only those spheres which touch the central sphere, the second contains all spheresin contact with the first shell but not contained within the first shell, and so on (115). Thedensity is defined as the leading coefficient of the quadratic fit of the number of spheres ineach shell to the shell number (116). The smallest values of the density are those of pointlattices. When physical packings were first proposed as models for the molecular structure of liquids,researchers sought to determine whether or not packings possessed some form of local order(6). This local order would take the form of small subunits with near-lattice structure,combined in some complicated way to produce the general disorder of the packing.The presence of an ordered structure in point patterns may not be obvious to the eye.Materials have been found in nature which are quasi-crystalline (117; 118). The arrangementof atoms in these materials would appear disordered to the eye, but x-ray diffraction revealsthat their structure can be modeled by a projection into R of a higher dimensional regularpoint lattice.Disorder can be defined as the absence of a regular point lattice structure for the spherecentres. The definition of disorder is not constructive, but instead merely identifies allconfigurations of spheres which are not point lattice structures. While it may be possibleto define statistics which identify the degree to which a particular configuration is distinct25rom a regular lattice, there will be many different ways of measuring these distinctions andnone of them will be clearly better than any other.It is necessary to distinguish between topological and geometric disorder. A packing istopologically ordered if the sphere centres can be continuously translated so as to transformits Delaunay triangulation into a point lattice without breaking any bonds (119). In a planarpacking, topological defects of the lattice structure can be easily identified and counted (120).A packing is geometrically disordered if any aspect of the arrangement of the spheresdiffers from that of a point lattice. If the spheres in a point lattice are reduced in size bya small amount and then randomly perturbed over a small distance, then topological ordercan be maintained in a geometrically disordered packing. Statistics which quantify geometricdisorder are based on local measurements which describe deviations from expected propertiesof point lattices.Statistics can be derived from the locations of contact points on individual spheres,expressed in spherical coordinates. The fourth- and sixth-order spherical harmonic functionscan be evaluated for the contact points on each sphere, and then averaged either over theindividual spheres or over many spheres. The fourth-order harmonics have non-zero averagesin the presence of local cubic lattice structure, while the sixth-order harmonics have non-zero averages in the presence of local hexagonal close-packed lattice structure. Averagesover single spheres have been used to study the structure of very large physical packings(29). Averages over many spheres were originally developed to study the emergence ofcrystallization in simulated liquids (121; 122), and improved averages were used to comparerealizations of simulated sphere packings (36). Averages over many spheres based on thesixth harmonic were used to describe the reduction in disorder observed in simulated spherepackings that had been cycled through the Jodrey-Tory packing simulation algorithm(123) .Two different many-sphere averages of the sixth harmonic were used to compare realizationsfrom three different rearrangement models (124). Neither average was powerful enough to26istinguish between realizations from different models.Statistics based on spherical harmonics cannot identify ordered structure that is foundwithin small clusters of neighbouring spheres. Statistics which reveal this type of orderedstructure can be constructed using the side lengths of Delaunay simplices, which are tetra-hedra formed by four Delaunay triangles which share common edges (125). In studies ofthe reduction of disorder over time in Jodrey-Tory packings, measures of tetrahedracity andquadroctahredracity were found to be more powerful at tracking changes than were were theaverages of spherical harmonics (80).Other measures of disorder can be found by estimating means and variances of localestimates of intensity (126; 127). Mathematical models for physical phenomena can be applied out of physical context to yieldnew statistics. If a composite of two electrically insulating materials can be simulated by asphere packing, then one or both of its phases can be assumed to be conducting. The bulkresistance of the material can be estimated if the internal structure is known. This resistancehas no physical meaning, but is based on the solution to a set of partial differential equationswhich use the packing structure as part of their boundary conditions.Statistics based on models for physical properties are traditionally used as response statis-tics in physical applications. When searching for statistics to explain the response, statisticsbased on physical models used out of context can be proposed as potential predictors of theresponse.
In frictionless flow, the flow experiences no internal resistance due to shear. Heat flow bymeans of conductance and the flow of electricity are examples. Models for frictionless flow27re relatively simple to construct.For electrical flow, the packing and its complement can be considered to be two materialswith differing electrical conductivity. On the boundary of the packing, two disjoint sets ofspheres can be considered to be connected to electrodes of infinite conductance. When aunit potential is placed between these electrodes, the potential can be calculated at all pointswithin the packing and its complement. From this, a bulk resistance for the composite can becalculated. This resistance defines a mean distance across the packing (128), whose form isdetermined both by the relative resistivity of the two phases and by the disordered structureof the packing.For a composite modeled by a disordered packing, the bulk resistance is difficult to cal-culate. There are no exact methods, and numerical methods require that the packing bediscretized very crudely. If equal resistance is assigned to all of the edges of the contactnetwork, then the potential at all vertices of the circuit can be found easily be means of theproperties of random walks through the network (129). The bulk resistance is a weighted av-erage of potentials at the electrode regions. If the network lies in the plane, the potentials canbe plotted and used as a diagnostic tool for finding problems with the potential-calculatingsoftware. A plot of the current along each edge, calculated from potential differences betweenthe defining vertices, is more useful than the potential as a descriptor of structure. Bulkresistance has been used to develop a test for the presence of anisotropy in sphere packingrealizations (130).Models for heat conduction can also be used. The spheres in a packing can be expandedin order to generate contact surfaces between neighbouring spheres. A plot of the bulk heatconductance of a packing as a function of the degree of spherical expansion has been usedto distinguish between packings generated by different models (131).28 .5.2 Statistics based on shearing flows
If a fluid flows through a fixed packing, it develops internal frictional losses which dependon the structure through which the fluid flows. If the packing itself is made to flow as apowder, then its flow is strongly affected by frictional losses arising from colliding spheres.These losses greatly complicate the modeling of the flows.If the packing is considered to be fixed, the flow of a fluid through its complement canbe modeled. All fluid flow is assumed to be laminar, since modeling turbulent flow througha complex structure is impractical. Major simplifications are required in order to be ableto calculate bulk flow properties. The complement can be represented by a piping networkwhose structure is determined from the Voronoi tessellation, and pipe resistances can beassigned on the basis of local void geometry (18). The flow through the void structurecan also be modeled by a lattice-Boltzmann model (132). The velocity profile and pressuregradient of the simulated flow can be used as statistics. If the diffusion of contaminantsthrough the flow is modeled (133), many different statistics related to contaminant flux andconcentration can be calculated.If the packing itself flows, this flow can take place in a vacuum, in air, in a liquid, orin both air and liquid. These flows can be modeled using DEM models which representcompaction processes that do not induce particle fracture (134; 135), flow in mixers anddrums (136; 137; 138), and flow during avalanches (139; 140). Flow in mutiaxial (141; 142;47), and shear test apparatus (143; 144) can also be simulated. In some cases, the simulationprogram would impose conditions on the structure of the sample. The best statistics wouldemerge from models for flows which are nearly packed. Flows of fluidized packings wouldbe less useful, since the packed structure would almost immediately be destroyed by thefluidization process.Use of DEM models is complicated by the complexity of modeling nearly packed states inmotion. In the applications literature, the unverifiable assumptions made in order to obtain29imulation results are not often clearly stated. Modeling some aspects of packing formationremain primitive. In physical powders, some energy is expended through deformation of thegrains. Modeling of this deformation can only be undertaken in the plane by means of finiteelement methods, and only for a small number of grains (145; 146). More idealized modelscan be used in model assessment. A plot of the number of sliding contacts as a function ofthe total number of contacts during simulated flow can distinguish between realizations ofdifferent packing models (147).
Given a packing of spheres, each sphere can be very slightly shrunken in place. The shrunkenspheres can be acted upon by a pre-specified potential, and the motion of the particles underthat potential can be tracked. The motion is confined to a region with periodic boundaries.Algorithms for simulating these motions have been developed as part of thermodynamicmodels for hard-sphere atoms (148; 149; 150), and also for non-spherical shapes (50). Oncethese algorithms have been run for some time, statistics representing pressure, temperature,and other thermodynamic quantities can be calculated.
New statistics can be developed through quantification of characteristic features of patternsseen in plots of packings. A simple plot of the discs in a realization can reveal the near-orderedarrangements of some planar packings. In space, cross-sections showing features of structurecan be made from tomographic data (29). Two-dimensional projections of tomographic datafrom physical sphere packings in cylinders reveal the ordered structure of the packing nearthe cylinder wall, in contrast to the disordered structure near the centre (32).In the planar case, extra features can be added to images of the packing. The Delaunaytriangulation can be plotted on top of the discs, and triangles can be coloured to identify30riples of discs in mutual contact (86). Edges of the triangulation which are associated withdefects can be coloured to illustrate the overall structure of those defects in large packings(64).
There is no evidence that statistics based on presently existing physical and statistical theorywill be powerful enough to distinguish between the realizations of a given pair of packingprocesses.Experimental physics is an important source of new statistics. Size segregation in agitatedgranular matter is a phenomenon that would not be predicted from theory or from simpleintuition. It occurs when physical packings of particles of differing sizes or densities areshaken (151; 152). An example of this is the Brazil nut paradox, in which the large andheavy Brazil nuts rise above smaller and lighter nuts when the container of nuts is shaken(42). Simulated sphere packings have been used to investigate the paradox (153; 154). Thesesimulation experiments could be used as sources of new statistics, particularly in cases wherethere is a size distribution among the spheres.
A statistic is smoothly defined if it is based on measurements taken at many locations inthe realization, and if the calculation of the statistic treats these measurements as if theyarose from a random sample. The observation locations can be taken on a fine grid orfrom a realization of a Poisson process, but their definition is completely independent of thestructure of the packing. To estimate a stationary 2-point correlation function at the vector r , a fine grid of sampling locations and its translate by r are imposed upon the realization.For each sampling location and its translate, a 1 is assigned if both points are in spheres,31lse a 0 is assigned. The estimator is the average of the measurements. While this statisticis considered to be a simple measure of interaction, its construction combines informationabout voids and neighbouring particles with no consideration of the overall structure of therealization.Other smoothly defined statistics can be constructed from the Delaunay triangulation orthe Voronoi tessellation. When a histogram is made of all triangle areas from a triangulation,all information about how the triangles fit together is lost. The estimation of the nearestneighbour function at a distance r is done by finding the length of the shortest triangulationedge from each vertex, and assigning a 1 if it is less than r in length and a 0 otherwise. Theestimator is again an average, and all other information about structure is ignored.Statistics based on physical properties differ from the smoothly averaging statistics inthat their construction depends on the disordered structure of the entire realization. Bulkelectrical resistance can be thought of as a distance across the network on which currentflows (128). This distance is based on a weighted average of all the possible paths that couldbe taken across the network, with the weights being chosen to be consistent with the lawsof electricity and the paths being defined by the disordered structure of the material. Thiswould also be true of any statistic based on path properties of a random walk, and alsoon the properties associated with mechanical deformation. The form of these statistics are structurally determined by the entire disordered structure of each realization.Statistics which are smoothly defined are ideal for use with the Poisson process, in whichthere is no interaction between the observed points. For a point process defined by the centresof packed spheres, smoothly defined statistics are often estimated on a length scale whichis much smaller than that of the individual sphere diameters. Information about individualsphere shape becomes confounded with interaction effects. This confounding makes thevalue of the statistic difficult to interpret in terms of features of the realization which canbe seen. Structurally determined statistics are determined by the disorder of the model.32hen physics-based statistics are used on packings, all of the information is captured on thesame length scale as the disorder. There is a vast literature associated with the electrical,granular-flow, and fluid-flow models which relates the models to easily-visualized physicalphenomena. Structurally-dependent statistics may be more useful than smoothly-definedstatistics at summarizing complex structural details that can be identified by eye. Finding the exact distribution of point process statistics is very difficult or impossible, exceptin the case of the Poisson process. If the realizations were very large, then it may be possibleto find asymptotic between-realization distributions for some statistics.The proofs of central limit theorems and laws of large numbers depend on the individualobservations being weakly dependent. If a point process statistic is a smoothly definedaverage of local measurements, then a central limit theorem can be found as long as thepair correlation function for the point process and some similar higher moments decay to 0sufficiently fast (155). This result is independent of the nature of the disorder within therealization. If the statistic is not a smoothly defined average, it is not clear that a centrallimit theorem could be proven.Laws of large numbers are of great importance in materials science and physics. If acomposite can be modeled by a packing with the spheres and the void having different re-sistivities, then sufficiently large specimens will have the same bulk resistance when a limitlaw for specimen resistance exists. In the case of ideal gases, equilibrium statistical thermo-dynamics (156; 157; 158) derives limit laws from probabilistic models of atomic behaviour.The limiting quantities, experienced as pressure, temperature, volume, and entropy, showno signs of variability arising from the dynamic chaos present at the atomic level. Attemptshave been made to create thermodynamic models for the bulk physical properties of granularmaterials (159; 160), but these attempts have been based on the assumption that any such33odel would take the same general form as classical thermodynamics. There is no reason tobelieve that limit laws for granular materials would resemble those for ideal gases, since themicroscale interactions of grains differ from those of atoms and molecules.Laws of large numbers are not applicable in many engineering problems. In the case ofgases, the limits are taken over collections of more than 10 atoms. For granular and com-posite materials, only thousands or millions of particles will be present, and the interactionsamong those particles will be much stronger than those found in gases. This suggests thatthe proper limiting behaviour is not that of convergence to a number, but rather conver-gence to a random variable. In the absence of theoretical arguments for constructing theselimits, limiting distributions for non-average-based statistics will need to be inferred fromlarge samples of model realizations.The mathematical intractability of probabilistic limit laws for disordered spatial patternshas resulted in the development of alternative approaches to finding constant large-samplelimits for statistics. If a composite is modeled by a lattice packing of spheres, then methodsfrom potential theory can be used to get very good approximations for the bulk electricalconductance (161; 162). When the packing is disordered, these methods cannot be usedon account of the lack of symmetry. Homogenization (163; 164) can be used to find thelimiting behaviour of statistics defined by differential equations which have spatially varyingboundary conditions. Often, these results are presented as bounds between which the limitingvalue of the statistic must fall. While the mathematical arguments behind homogenizationrepresent an elegant use of the calculus of variations, the underlying assumptions for thesearguments generally cannot be related to the microstructure of any real material. In onecase, a mathematical argument has been constructed to establish that a microstructure existswhich has the coherent potential approximation for bulk electrical conductance as its limit(165). 34 Inference
The central problem of inference for sphere packings is model assessment. No proper statis-tical methods yet exist for determining if a fitted model for a packing has anything beyonda superficial resemblance to the physical system that it is supposed to represent. This hasserious implications in physics and engineering, where complex simulation models basedon unverifiable assumptions are being used. The history of thermodynamics suggests howmethods for assessing fitted models can be developed.No algorithm for the simulation of a packing creates realizations in a manner whichfaithfully represents the physics of packing formation. The DEM models come the closest,but even these are based on assumptions about unobservable mechanisms of energy lossthrough friction and deformation. Parameters are few in number relative to the size andcomplexity of the realizations. Parameters often have no physical meaning, especially in thecase of rearrangement and ballistic models.Model parameters can be fit by methods similar to the method of moments, such as themethod of minimum contrast (166). This method is based on choosing a set of parametervalues which minimizes the difference between a statistic calculated from the model anda statistic calculated from the data. Often, the statistic chosen is the K − function. Thismethod ensures that there is some resemblance between the model output and the data, butit says nothing about how far that resemblance can be extended. Since the K − function lacksthe power to distinguish point patterns which can be distinguished by eye (167), there is noreason to believe that a fit based on it will result in a good fit for other statistics calculatedfrom the packing.A fitted model must be assessed to determine if the model is useful or if the fit is su-perficial. Since every realization from the model or from the data may be different, thisassessment must be carried out using descriptive statistics which can summarize common35spects of different realizations of the same process and which can identify evidence of impor-tant differences between the model and the data. Since there is no guide to which statisticsare useful, it is necessary to calculate many different statistics for each realization and thento compare their distributions by means of multivariate statistical methods. In physicalapplications, the need for many statistics will require the development of accurate methodsfor imaging the three-dimensional structure of static and flowing physical packings.Since no model directly and accurately represents the formation of a physical packing,it is likely that at least one statistic can be found which gives evidence that realizationsof a model differ from those of the data. This difference only matters if that statistic hasan effect on the physical property being modeled. In a typical application, the packing isbeing used as part of a model for some specific physical phenomenon such as the electricalconductance of a composite material. If the phenomenon of interest is summarized by a setof response statistics, these statistics can be used to fit the parameters of the model. Theonly descriptive statistics that are needed for model assessment are those which affect thedistribution of the response statistics. In general, there is no way to theoretically determineall of the statistics which can affect the response. Further modeling and testing will berequired to establish if statistics which identify a difference between the data and the modelalso affect the response.Some progress has been made towards the development of methods for model assess-ment in packings. Various methods have been developed for describing deviations fromcomplete spatial randomness in point patterns (168), but it is not clear if these methods willbe useful for distinguishing between different types of non-Poisson disorder. Residuals forpoint processes have been developed (169) which can be used to detect spatial trend andinterpoint interaction effects. Physicists have used descriptive statistics to quantify differ-ences between simulation models for packings (170; 171; 172) and between different typesof physical packings (29). Descriptive statistics have also been used to quantify the effects36f changing parameters in DEM models (173). These efforts have been used to suggest thatsome models may be physically wrong, but this condition is too strong. The identification ofdifferences between models is only the beginning of the development of objective and usefulmodel assessment procedures. Sphere packings form a very useful class of models for physical phenomena which cannot bemodeled formally. They form a class of statistical models consisting entirely of physical andcomputer simulation procedures which need to be fit to data.When packings are used as models, the fitted models must be further assessed to deter-mine whether or not the resemblance of the model to the physical system is superficial oruseful. If this is not done, there is a severe risk that the fitted model may be misleadingas a physical model. The development of assessment procedures requires taking the sameapproach as early physicists; a specific response must be defined, statistics which affect theresponse must be identified, and these statistics must be compared for the model and thedata. The phenomena being modeled are as well-understood as the properties of gases werein the 18 th century, and require a similar experimental approach to identify a useful way tosummarize everything of importance about them. Acknowledgments
I would like to thank Alan Karr of the National Institute of Statistical Sciences and SurendraShah of the Centre for Advanced Cement-Based Materials for providing the first circum-stances in which I found the need to fit sphere packing models to data. Bruce Ankenman,Takeru Igusa, Sanjay Jaiswal, Eric Slud, Dietrich Stoyan, Rolf Turner, and Viqar Husain also37rovided much useful advice. This work was funded by NISS, the University of Maryland,and NSERC.
References [1] C.A. Rogers.
Packing and Covering . Cambridge, 1964.[2] P.M. Gruber and C.G. Lekkerkerker.
Geometry of Numbers . North-Holland, 1987.[3] C. Zong.
Sphere Packings . Springer, New York, 1999.[4] J.H. Conway and N.J.A. Sloane.
Sphere Packings, Lattices, and Groups . Springer, 3rdedition, 1999.[5] K. Stephenson.
Introduction to Circle Packing . Cambridge, 2005.[6] J.D. Bernal. A geometric approach to the structure of liquids.
Nature , 183:141–147,1959.[7] J.D. Bernal. The structure of liquids.
Proc. R. Soc. A , 280(1382):299–322, 1964.[8] J.D. Bernal. Co-ordination of randomly packed spheres.
Nature , 188:910–911, Decem-ber 10 1960.[9] G.D. Scott. Packing of spheres.
Nature , 188:908–909, December 10 1960.[10] T. Aste and D. Weaire.
The Pursuit of Perfect Packing . Institute of Physics, 2000.[11] H.M. Jaeger, S.R. Nagel, and R.P. Behringer. Granular solids, liquids, and gases.
Rev.Mod. Phys. , 68:1259, 1996.[12] P.G. Degennes. Granular matter: A tentative view.
Rev. Mod. Phys. , 68:1259–1273,1996. 3813] P. Ball.
The Self-Made Tapestry . Oxford, 2001.[14]
ASM Handbook Vol. 7: Powder Metal technologies and Applications . ASM Interna-tional, 1998.[15] G.Y. Onoda and E.G. Liniger. Random loose packings of spheres and dilatancy onset.
Phys. Rev. Lett. , 64(22):2727–2730, 1990.[16] S. Mindess and J.F. Young.
Concrete . Prentice-Hall, 1981.[17] W. Shiu, F.V. Donz´e, and L. Daudeville. Influence of the reinforcement on penetra-tion and perforation of concrete targets: A discrete element analysis.
EngineeringComputation , 26(1):29–45, 2009.[18] C.F. Chu and K.M. Ng. Flow in packed tubes with a small tube to particle diameterratio.
AICHE J. , 35:148–158, 1989.[19] A. Okabe, B. Boots, K. Sugihara, and K. Chiu.
Spatial Tessellations: Concepts andApplications of Voronoi Diagrams . Wiley, 2nd edition, 2000.[20] D. Ruelle.
Statistical Mechanics: Rigorous Results . Benjamin, 1974.[21] J. Møller and R.P. Waagepetersen.
Statistical Inference and Simulation for SpatialPoint Processes . Chapman and Hall, 2004.[22] M. Locatelli and U. Raber. Packing circles in a square: a deterministic global opti-mization approach.
Disc. App. Math. , 122(1-3):139–166, 2002.[23] R.L. Graham and B.D. Lubachevsky. Dense packings of equal discs in an equilateraltriangle.
Electon. J. Combin. , 2:A1, 1995.[24] M.C. Mark´ot. Optimal packing of 28 equal circles in a unit square - the first reliablesolution.
Numer. Algorithms , 37:253–261, 2004.3925] T. Ikegami. Contacts and coordination numbers in a compact of polyhedral particles.
Journal of the American Ceramic Society , 79:148–152, 1996.[26] J.D. Bernal, I.A. Cherry, J.L. Finney, and K.R. Knight. An optical machine for mea-suring sphere coordinates in random packings.
J. Phys. E , 3(5):388–390, 1970.[27] J. Bruji´c, S.F. Edwards, I. Hopkinson, and H.A. Makse. Measuring the distributionof interdroplet forces in a compressed emulsion system.
Physica A , 327(3-4):210–212,2003.[28] T.I. Quickenden and G.K. Tan. Random packing in two dimensions and the structureof the monolayers.
J. Colloid Interf. Sci. , 48:382–392, 1974.[29] T. Aste, M. Saadatfar, and T.J. Senden. Geometrical structure of disordered spherepackings.
Physical Review E , 71(6):061302, 2005.[30] T. Aste. Variations around disordered close packing.
J.Phys. Condens. Matter ,17:S2361–S2390, 2005.[31] P. Richard, P. Philippe, F. Barbe, S. Bourl`es, X. Thibault, and Bideau. D. Analysisby X-ray microtomography of a granular packing undergoing compression.
Phys. Rev.E , 68(2):020301, 2003.[32] G.T. Seidler, G. Martinez, L.H. Seeley, K.H. Kim, E.A. Behne, S. Zaranek, B.D. Chap-man, S.M. Heald, and D.L. Brewe. Granule-by-granule reconstruction of a sandpilefrom X-ray microtomography data.
Phys. Rev. E , 62(6):8175–8181, 2000.[33] A.J. Sederman, P. Alexander, and L.F. Gladden. Structure of packed beds probed bymagnetic resonance imaging.
Powder Tech. , 117(3):255–269, 2001.[34] M.M. Kohonen, D. Geromichalos, M. Scheel, C. Schier, and S. Herminghaus. Oncapillary bridges in wet granular materials.
Physica A , 339:7–15, 2004.4035] M. Toiya, J. Hettinga, and W. Losert. 3d imaging of particle motion during penetrom-eter testing.
Granular Matter , 9(5):323–329, 2007.[36] S. Torquato. Is random close packing of spheres well-defined?
Phys. Rev. Lett. ,84(10):2064–2067, 2000.[37] D.S. Grebenkov, M.P. Ciamarra, M. Nicodemi, and A. Coniglio. Flow, ordering, andjamming of sheared granular suspensions.
Phys. Rev. Lett. , 100(7):078001, 2008.[38] O. Pouliquen, M. Nicolas, and P. D. Weidman. Crystallization of non-brownian spheresunder horizontal shaking.
Phys. Rev. Lett. , 79(19):3640–3643, 1997.[39] J.W.S. Cassells.
An Introduction to the Geometry of Numbers . Springer, Berlin, 1959.[40] T.C. Hales. The sphere packing problem.
J. Comput. Appl. Math. , 44:41–76, 1992.[41] A. Yang, C.T. Miller, and L.D. Turcoliver. Simulation of correlated and uncorrelatedpacking of random size spheres.
Phys. Rev. E , 53(2):1516–1524, 1996.[42] A. Rosato, K.J. Strandburg, F. Prinz, and R.H. Swendsen. Why Brazil nuts are on top:Size segregation of particulate matter by shaking.
Phys. Rev. Lett. , 58(10):1038–1040,1987.[43] J.L. Finney. Random packings and the structure of simple liquids I. The geometry ofrandom close packing.
Proc. R. Soc. A , 319:479–473, 1970.[44] J. Tobochnik and P.M. Chapin. Monte Carlo simulation of hard spheres near randomclosest packing using spherical boundary conditions.
Phil. Trans. R. Soc. A , 88:5824–5830, 1988.[45] K. Gotoh, W.S. Jodrey, and E.M. Tory. Variation on local packing density near thewall of randomly packed bed of equal spheres.
Powder Tech. , 20:257–260, 1978.4146] Taguchi. I., M. Kurashige, and K. Imai. Effects of cubic container’s wall or floor onrandom packing structures of spherical particles.
JSME Int. J. Series A , 49:265–272,2006.[47] L. Cui, C. O’Sullivan, and S. O’Neill. An analysis of the triaxial apparatus using amixed boundary three-dimensional discrete element model.
G´eotechnique , 57(10):831–844, 2007.[48] B.J. Buchalter and R.M. Bradley. Orientational order in random packings of ellipses.
Phys. Rev. A , 46(6):3046–3056, 1992.[49] S.R. Williams and A.P. Philipse. Random packings of spheres and spherocylinderssimulated by mechanical contraction.
Phys. Rev. E , 67(5):051301, 2003.[50] A. Donev, S. Torquato, and F.H. Stillinger. Neighbour list collision-driven moleculardynamics simulation for non-spherical hard particles. II. Applications to ellipses andellipsoids.
J. Comp. Phys. , 202(2):765–793, 2005.[51] W. Man, A. Donev, F.H. Stillinger, M.T. Sullivan, W.B. Russel, D. Heeger, S. Inati,S. Torquato, and P.M. Chaikin. Experiments on random packings of ellipsoids.
Phys.Rev. Lett. , 94:198001, 2005.[52] K. Nandakumar, Y. Shu, and K.T. Chuang. Predicting geometrical properties ofrandom packed beds from computer simulation.
AICHE J. , 45(11):2286–2296, 1999.[53] R. Jullien, A. Pavlovitch, and P. Meakin. Random packings of spheres built withsequential methods.
J. Phys. A , 25:4103–4113, 1992.[54] M.J. Vold. A numerical approach to the problem of sediment volume.
J. Colloid Sci. ,14:168–174, 1959. 4255] M.J. Vold. The sediment volume in dilute dispersions of spherical particles.
J. Phys.Chem. , 64:1616–1619, 1960.[56] M.J. Vold. Computer simulation of floc formation in a colloidal suspension.
J. ColloidSci. , 18:684–695, 1964.[57] M.J. Vold. Sediment volume and structure in dispersions of anisometric particles.
J.Phys. Chem. , 63:1608–1612, 1959.[58] C.H. Bennett. Serially deposited amorphous aggregates of hard spheres.
J. App. Phys. ,43(6):2727–2734, 1972.[59] D.J. Adams and A.J. Matheson. Computation of dense random packing of hardspheres.
J. Chem. Phys. , 56:1989–1994, 1972.[60] G.Y. Onoda and J. Toner. Deterministic, fractal defect structures in close packings ofhard discs.
Phys. Rev. Lett. , 57:1340–1343, 1986.[61] W.M. Visscher and M. Bolsterli. Random packing of equal and unequal spheres in twoand three dimensions.
Nature , 239:504–507, October 27 1972.[62] L. Rouill´e, J.-M. Missiaen, and G. Thomas. Collective random packing of disks in theplane under the influence of a weak central force.
J. Phys. Cond. Mat. , 2:3041–3051,1990.[63] N.D. Aparacio and A.C.F. Cocks. On the representation of random packings of spheresfor sintering simulations.
Acta Metall. Mater. , 43:3873–3884, 1995.[64] P. Meakin and R. Jullien. Periodic disc packings generated by random deposition innarrow channels.
Europhys. Lett. , 15:851–856, 1991.4365] W Soppe. Computer simulation of random packings of hard spheres.
Powder Tech. ,62:189–196, 1990.[66] G.E. Mueller. Numerical simulation of packed beds with monosized spheres in cylin-drical containers.
Powder Tech. , 92:179–183, 1997.[67] G.E. Mueller. Numerically packing spheres in cylinders.
Powder Tech. , 159(2):105–110,2005.[68] D. Coelho, J.F. Thovert, and P.M. Adler. Geometrical and transport properties ofrandom packings of spheres and aspherical particles.
Phys. Rev. E , 55(2):1959–1978,1997.[69] Y.T. Feng, K. Han, and D.R.J. Owen. Filling domains with disks: an advancing frontapproach.
Int. J. Numer. Methods Engrg. , 56(5):669–713, 2003.[70] K. Han, Y.T. Feng, and D.R.J. Owen. Sphere packing with a geometric compressionalgorithm.
Powder Tech. , 155(1):33–41, 2005.[71] R. Jullien, P. Meakin, and A. Pavlovitch. Particle size segregation by shaking in two-dimensional packings.
Europhys. Lett. , 22:523–528, 1993.[72] R. Jullien, A. Pavlovitch, and P. Meakin. Three-dimensional model for particle-sizesegregation by shaking.
Phys. Rev. Lett. , 69:640–643, 1992.[73] G.C. Barker and A. Mehta. Vibrated powders: structure, correlations, and dynamics.
Phys. Rev. A , 45:3435–3446, 1992.[74] X. Jia and R.A. Williams. A packing algorithm for particles of arbitrary shape.
PowderTech. , 120:175–186, 2001. 4475] W.S. Jodrey and E.M. Tory. Computer simulation of isotropic, homogeneous, denserandom packings of equal spheres.
Powder Tech. , 30:111–118, 1981.[76] W.S. Jodrey and E.M. Tory. Computer simulation of close random packing of equalspheres.
Phys. Rev. A , 32(4):2347–2351, 1985.[77] M. Bargiel and J. Moscinski. C-language program for the irregular close packing ofhard spheres.
Comp. Phys. Comm. , 64:183–192, 1991.[78] A. Bezrukov and D. Stoyan. Simulation and statistical analysis of random packings ofellipsoids.
Part. Part. Syst. Char. , 23(5):388–398, 2006.[79] A.S. Clarke and J.D. Wiley. Numerical simulation of the dense random packing of abinary mixture of hard spheres: amorphous metals.
Phys. Rev. B , 35:7350–7356, 1987.[80] K. Lochmann, A. Anikeenko, A. Elsner, N. Medvedev, and D. Stoyan. Statisticalverification of crystallization in hard sphere packings under densification.
Eur. Phys.J. B , 53(1):67–76, 2006.[81] Speedy. R.J. Random jammed packings of hard discs and spheres.
J. Phys. Cond.Mat. , 10:4185–4194, 1998.[82] A.Z. Zinchenko. Algorithm for random close packing of spheres with periodic boundaryconditions.
J. Comp. Phys. , 114(2):298–307, 1994.[83] E.L. Hinrichsen, J. Feder, and T. Jøssang. Random packing of discs in two dimensions.
Phys. Rev. A , 41(8):4199–4209, 1990.[84] B.D. Lubachevsky. How to simulate billiards and similar systems.
J. Comp. Phys. ,94(2):255–283, 1991. 4585] B.D. Lubachevsky and F.H. Stillinger. Geometric properties of random disc packings.
J. Stat. Phys. , 60(2):561–583, 1990.[86] B.D. Lubachevsky, F.H. Stillinger, and E.N. Pinson. Disks vs. spheres: Contrastingproperties of random packings.
J. Stat. Phys. , 64(2):510–524, 1991.[87] A. Donev, F.H. Stillinger, P.M. Chalkin, and S. Torquato. Unusually dense crystalpacking of ellipsoids.
Phys. Rev. Lett. , 92:255506, 2004.[88] C.S. O’Hern, S.A. Langer, A.J Liu, and S.R. Nagel. Force distributions near jammingand glass transitions.
Phys. Rev. Lett. , 86(1):111–114, 2001.[89] G.T. Nolan and P.E. Kavanaugh. Computer simulation of random packing of hardspheres.
Powder Tech. , 72:149–155, 1992.[90] P.A. Cundall and O.D.L. Strack. A discrete numerical model for granular assemblies.
Geotechnique , 29:47–65, 1979.[91] H.P. Zhu, Z.Y. Zhou, R.Y. Yang, and A.B. Yu. Discrete particle simulation of par-ticulate systems: A review of major applications and findings.
Chemical EngineeringScience , 63:5728–5770, 2008.[92] H.P. Zhu, Z.Y. Zhou, R.Y. Yang, and A.B. Yu. Discrete particle simulation of particu-late systems: Theoretical developments.
Chemical Engineering Science , 62:3378–3396,2007.[93] K.Z.Y. Yen and T.K. Chaki. A dynamic simulation of particle rearrangement in powderpackings with realistic interactions.
J. App. Phys. , 71(7):3164–3173, 1992.[94] P.J. Diggle, J. Besag, and J.T. Gleaves. Statistical analysis of spatial point patternsby means of distance methods.
Biometrics , 32(3):659–667, 1976.4695] R.Y. Yang, R.P. Zou, and A.B. Yu. Computer simulation of the packing of fine parti-cles.
Phys. Rev. E , 62(23):3900–3908, 2000.[96] F. Bertrand, T. Gange, E. Desauliniers, D. Vidal, and R.E. Hayes. Simulation of theconsolidation of paper coating structures: probabilistic versus deterministic models.
Computers and Chemical Engineering , 28(12):2595–2604, 2004.[97] H.G. Matuttis, S. Luding, and H.J. Herrmann. Discrete element simulations of densepackings and heaps made of spherical and non-spherical particles.
Powder Technology ,109(1-3):278–292, 2000.[98] G. Matheron.
Random Sets and Integral Geometry . Wiley, New York, 1975.[99] I. Molchanov.
Theory of Random Sets . Springer, 2005.[100] G.E. Mueller. Radial void fraction distributions in randomly packed fixed beds ofuniformly sized spheres in cylindrical constainers.
Powder Tech. , 72:269–275, 1992.[101] S. Torquato and G. Stell. Microstructure of two-phase random media.V. The n -pointmatrix probability functions for impenetrable spheres. J. Chem. Phys. , 82(2):980–987,January 1985.[102] P. Bhanu Prasad and J.P. Jernot. Topological description of the densification of agranular medium.
J. Microsc. , 163:211–220, 1991.[103] J.P. Jernot, P. Bhanu Prasad, and P. Demaleparde. Three-dimensional simulation offlow through a porous medium.
J. Microsc. , 167:9–21, 1192.[104] K. Gotoh, M. Nakagawa, M. Furuuchi, and A. Yoshiga. Pore size distribution inrandom assemblies of equal spheres.
J. Chem. Phys. , 85:3078–3080, 1986.47105] D.J. Daley and D. Vere-Jones.
Introduction To The Theory of Point Processes, Volume1 . Springer, 2nd edition, 2003.[106] K. Schladitz and A.J. Baddeley. A third-order point process characteristic.
Scand. J.Stat. , 27(4):657–671, 2000.[107] B.D. Ripley.
Statistical Inference For Spatial Processes . Wiley, 1988.[108] M.N.M. Van Lieshout and A.J. Baddeley. A non-parametric measure of spatial inter-action in point patterns.
Adv. Appl. Prob. , 28(2):337, 1996.[109] T. Aste and T. Di Matteo. Emergence of Gamma distributions in granular materialsand packing models.
Physical Review E , 77(2):012309, 2008.[110] M. Tassopoulos and D.E. Rosner. Microstructural descriptors characterizing granulardeposits.
AICHE J. , 38:15–25, 1992.[111] R. Jullien, P. Jund, D. Caprion, and D. Quitmann. Computer investigation of long-range correlations and local order in random packings of spheres.
Phys. Rev. E ,54(6):6035–6041, 1996.[112] R.T. Fowler and J.R. Glastonbury. The flow of granular solids through orifices.
Chem.Eng. Sci. , 10:150–156, 1959.[113] P.T. Stainforth. A flow-factor for no-arching at the transition between hopper andtrunk in tall cylindrical silos.
Powder Tech. , 9(1):53–55, 1974.[114] K. To, P.-Y. Lai, and H.K. Pak. Jamming of granular flow in a two-dimensional hopper.
Phys. Rev. Lett. , 86(1):71–74, 2001.[115] G.O. Brunner. The properties of coordination sequences and conclusions regarding thelowest possible density of zeolites.
Physica A , 29(1):41–45, 1979.48116] M. O’Keeffe. N-dimensional diamond, sodalite, and rare sphere packings.
Acta Crys-tallogr. Sect. A , 47(6):748–753, 1991.[117] D. Shechtman, I. Blech, D. Gratias, and J.W. Cahn. Metallic phase with long-rangeorientational order and no translational symmetry.
Phys. Rev. Lett. , 53:1951–1953,1984.[118] M. Senechal.
Quasicrystals and Geometry . Cambridge, 1995.[119] J.M. Ziman.
Models of Disorder . Cambridge, 1979.[120] P. Meakin and R. Jullien. Two-dimensional defect-free random packing.
Europhys.Lett. , 14(1):667–672, 1991.[121] P.J. Steinhardt, D.R. Nelson, and M. Ronchetti. Icosahedral bond orientational orderin supercooled liquids.
Phys. Rev. Lett. , 47(18):1297–1300, 1981.[122] P.J. Steinhardt, D.R. Nelson, and M. Ronchetti. Bond-orientational order in liquidsand gases.
Phys. Rev. B , 28(2):784–805, 1983.[123] P. Richard, L. Oger, J.-P. Troadec, and A. Gervois. Geometrical characterization ofhard-sphere systems.
Phys. Rev. E , 60(4):4551–4559, 1999.[124] A.R. Kansal, S. Torquato, and F.H. Stillinger. Diversity of order and densities injammed hard-particle packings.
Phys. Rev. E , 66(4):041109, 2002.[125] A.V. Anikeenko and N.N. Medvedev. Homogeneous crystallization of the Lennard-Jones liquid.
J. Struct. Chem. , 47(2):267–276, 2006.[126] S. Torquato, T.M. Truskett, and P.G. Debenedetti. Is random close packing of sphereswell-defined?
Phys. Rev. Lett. , 84(10):2064–2067, 2000.49127] S. Torquato and F.H. Stillinger. Local density fluctuations, hyperuniformity, and ordermetrics.
Phys. Rev. E , 68(4):041113, 2003.[128] D.J. Klein and M. Randic. Resistance distance.
J. Math. Chem. , 12:81–95, 1993.[129] P.G. Doyle and J.L. Snell.
Random Walks and Electrical Networks . MathematicalAssociation of America, 1984.[130] J.D. Picka and T.D. Stewart. A conduction-based statistic for description of orientationin anisotropic sphere packings. In
Powders and Grains 2005 Proceedings , pages 329–332. Balkema, 2005.[131] C. Argento and D. Bouvard. Modeling the effective thermal conductivity of a randompacking of spheres through densification.
Int. J. Heat Mass Tran. , 39(7):1334–1350,1996.[132] S. Chen, K. Diemer, G.D. Doolen, K. Eggert, C. Fu, S. Gutman, and B.J. Travis.Lattice gas automata for flow through porous media.
Physica D , 47:72–84, 1991.[133] R.S. Maier, D.M. Kroll, R.S. Bernard, S.E. Howington, J.F. Peters, and H.T. Davis.Hydrodynamic dispersion in confined packed beds.
Phys. Fluids , 15(12):3795–3815,2003.[134] O. Coube, A.C.F. Cocks, and C.-Y. Wu. Experimental and numerical study of diefilling, powder transfer and die compaction.
Powder Met. , 48(1):68–76, 2005.[135] Y. Chung and J. Ooi. Influence of discrete element model parameters on bulk behaviorof a granular solid under confined compression.
Part. Sci. and Tech. , 26(1):83–96, 2008.[136] B. Chaudhuri, F.J. Muzzio, and M.S. Tomassone. Modeling of heat transfer in granularflow in rotating vessels.
Chem. Eng. Sci. , 61(19):6348–6360, 2006.50137] P. Portillo, F. Muzzio, and M. Ierapetritou. Hybrid DEM-compartment modelingapproach for granular mixing.
AICHE J. , 53(1):119–128, 2007.[138] Sato Y., Nakamura H., and Watano S. Numerical analysis of agitation torque andparticle motion in a high shear mixer.
Powder Tech. , 186(2):130–136, 2008.[139] C.S. Chang. Discrete element method for slope stability analysis.
J. Geotech. Eng.ASCE , 118:1889–1905, 1992.[140] R. Brewster, G.S. Grest, J.W. Landry, and A.J. Levine. Plug flow and the breakdownof Bagnold scaling in cohesive granular flows.
Phys. Rev. E , 72(061301), 2005.[141] T.T. Ng and R. Dobry. Numerical simulations of monotonic and cyclic loading ofgranular soil.
J. Geotech. Eng. ASCE , 120(2):388–403, 1994.[142] C.T. David, R. Garcia-Rojo, H.J. Herrmann, and S. Luding. Powder flow testing with2d and 3d biaxial and triaxial simulations.
Part. Part. Syst. Charact. , 24(1):29–33,2007.[143] C. Thornton and L. Zhang. A numerical examination of shear banding and simpleshear non-coaxial flow rules.
Phil. Mag. , 86(21-22):3425–3452, 2006.[144] J. H¨artl and J.Y. Ooi. Experiments and simulations of direct shear tests: porosity,contact friction, and bulk friction.
Granular Matter , 10(4):263–271, 2008.[145] A. Zavaliangos. A numerical study of the development of tensile principal stressesduring die compaction.
Part. Sci. Tech. , 21(2):105–115, 2003.[146] A.T. Procopio and A. Zavaliangos. Simulation of multi-axial compaction of granularmedia from loose to high relative densities.
J. Mech. Phys. Solids , 53:1523–1551, 2005.51147] A. Mcnamara and H. Herrmann. Measurement of indeterminacy in packings of per-fectly rigid discs.
Phys. Rev. E , 70:061303, 2004.[148] T. Gruhn and P.A. Monson. Molecular dynamics simulations of hard sphere solidifi-cation at constant pressure.
Phys. Rev. E , 64(6):061703, 2001.[149] T. Gruhn and P.A. Monson. Isobaric molecular dynamics simulations of hard spheresystems.
Phys. Rev. E , 63(6):061106, 2001.[150] A. Donev, S. Torquato, and F.H. Stillinger. Neighbour list collision-driven moleculardynamics simulation for non-spherical hard particles. I. Algorithmic details.
J. Comp.Phys. , 202(2):737–764, 2005.[151] J.C. Williams. The segregation of particulate materials.
Powder Tech. , 15(2):245–251,1976.[152] J. Bridgwater. Fundamental powder mixing mechanisms.
Powder Tech. , 15(2):215–236,1976.[153] A. Rosato, K.J. Strandburg, F. Prinz, and R.H. Swendsen. Monte Carlo simulation ofparticulate matter segregation.
Powder Tech. , 49(1):59–69, 1986.[154] D.C. Hong, P.V. Quinn, and S. Luding. Reverse Brazil nut problem: Competitionbetween percolation and condensation.
Phys. Rev. Lett. , 86(1):3423–3426, 2001.[155] M.D. Penrose and J.E. Yukich. Limit theory for random sequential packing and depo-sition.
Ann. App. Prob. , 12(1):272–301, 2002.[156] L.D. Landau and E.M. Lifshitz.
Statistical Physics . Addison-Wesley, 1958.[157] S.G. Brush.
Statistical Physics and the Atomic Theory of Matter . Princeton, 1983.52158] R.A. Minlos.
Introduction to Mathematical Statistical Mechanics . American Mathe-matical Society, 2000.[159] S.F. Edwards and R.B.S. Oakeshott. Theory of powders.
Physica A , 157:1080–1090,1989.[160] C.C. Mounfield and S.F. Edwards. A model for the packing of irregularly shapedgrains.
Physica A , 210:301–316, 1994.[161] Lord Rayleigh. On the influence of obstacles arranged in rectangular order upon theproperties of a medium.
Phil. Mag. , 34:481–502, 1892.[162] W.T. Perrins, D.R. McKenzie, and R.C. McPhedran. Transport properties of regulararrays of cylinders.
Proc. R. Soc. A , 369(1737):207–225, 1979.[163] U. Hornung.
Homogenization and Porous Media . Springer, 1997.[164] G.W. Milton.
Theory of Composites . Cambridge, 2002.[165] G.W. Milton. The coherent potential approximation is a realizable effective mediumscheme.
Commun. Math. Phys. , 99(4):463–500, 1985.[166] P.J. Diggle and R.J. Gratton. Monte carlo methods of inference for implicit statisticalmodels.
J.R.S.S. B , 46:193–212, 1984.[167] A.J. Baddeley and B.W. Silverman. A cautionary example on the use of 2nd-ordermethods for analyzing point patterns.
Biometrics , 40(4):1089–1093, 1984.[168] N.A.C. Cressie.
Statistics For Spatial Data . Wiley, New York, 1993.[169] A. Baddeley, R. Turner, J. Møller, and M. Hazleton. Residual analysis for spatial pointprocesses.
J. R. Statist. Soc. B , 67(5):617–666, 2005.53170] P. Meakin and A.T. Skjeltorp. Application of experimental and numerical models tothe physics of multiparticle systems.
Adv. Phys. , 42(1):1–127, 1993.[171] R. Jullien and P. Meakin. Computer simulations of steepest descent ballistic deposition.
Colloids Surface A , 165:405–422, 2000.[172] A. Wouterse and A.P. Philipse. Geometrical cluster ensemble analysis of random spherepackings.
J. Chem. Phys. , 125:194709, 2006.[173] Z.P. Zhang, L.F. Liu, Y.D. Yuan, and A.B. Yu. A simulation study of the effects ofdynamic variables on the packing of spheres.