Image Anal Stereol 2001;20:203-206 Original Research Paper SPATIAL STATISTICS FOR SIMULATED PACKINGS OF SPHERES Alexander Bezrukov1, Dietrich Stoyan1 and Monika Bargie³2 1Institut fu¨r Stochastik, TU Bergakademie Freiberg, 09596 Freiberg, Germany, 2Institute of Computer Science, AGH Krakow, al. Mickiewicza 30, 30-059 Krakow, Poland e-mail: stoyan@orion.hrz.tu-freiberg.de (Accepted August 24, 2001) ABSTRACT This paper reports on spatial-statistical analyses for simulated random packings of spheres with random diameters. The simulation methods are the force-biased algorithm and the Jodrey-Tory sedimentation algorithm. The sphere diameters are taken as constant or following a bimodal or lognormal distribution. Standard characteristics of spatial statistics are used to describe these packings statistically, namely volume fraction, pair correlation function of the system of sphere centres and spherical contact distribution function of the set-theoretic union of all spheres. Furthermore, the coordination numbers are analysed. Keywords: contact distribution, coordination number, pair correlation function, simulation, sphere packing, volume fraction. INTRODUCTION Random packings of spheres with constant and random diameters play an important role in many branches of physics and engineering. Simulated packings serve as models for real packings of particles, e.g. in the context of particle science, where, however, the assumption that the particles are spheres is often a simplification. Also many porous media can be represented as packed arrangements of spheres. In these and other cases, the characteristics of the simulated packings should be at least similar to those of the real structures investigated. This concerns primarily mean value characteristics such as porosity and mean coordination number and, additionally, functions such as coordination number distribution, pair correlation function and spherical contact distribution function. There is a wide variety of sphere packing algorithms. This paper uses the force-biased algorithm (Mos´cin´ski et al., 1989; Mos´cin´ski and Bargie³, 1991), and the Jodrey-Tory sedimentation algorithm (Jodry and Tory, 1979). The first algorithm starts with a very dense configuration of large spheres which can even overlap. It reduces stepwise the diameters and moves the spheres in order to reduce the degree of overlappings. The second algorithm simulates the successive packing of a container with spheres falling in the gravitational field. A sphere with three contacts with earlier packed spheres is considered as being in its final position. While the first algorithm is able to produce very dense nearly isotropic packings, the second one yields loose packings of volume fraction VV = 0.58 in the case of equal spheres, which show a light gradient in vertical direction. Other known packing algorithms use the discrete element method (Stroeven and Stroeven, 2000). For the physicist, engineer and statistician then the question arises which of these algorithms to choose. The decision should be based on statistical properties of the packings produced, i.e. as said above, the model characteristics should be as close as possible to the empirical counterparts and the diameter distributions should be identical. Therefore, this paper presents results of the two algorithms for some interesting cases of diameter distributions and packing densities. Fig. 1. A simulated sphere packing with a lognormal diameter distribution and volume fraction VV=0.70. 203 Bezrukov A et al: Spatial statistics Fig. 2. A simulated sphere packing with bimodal diameter distribution with diameters 1 and 10 and VV=0.80. 0.425 0.4 ^oC "^- ^ / / "" - >. -—- / I 0.375 ¦"- — _ _ — -- t L> " o L 0.35 • x o ¦ "'-"-";-. ~. /; 9, *-'/ o o o ¦ '¦"' o- • Sh X o ¦ o --L>--- ° 0.325 * X • 0.3 X • 0.275 X X • X X 0.2 0.4 0.6 0.8 1 relative volume of large particles Fig. 3. Porosity versus relative volume v of smaller spheres. The curves result from simulations with the sedimentation algorithm, the points are experimental data taken from Sohn and Moreland (1968). Size ratio 1 : k: o,——0.5; ,-------0.39; ,------0.29; x,..... 0.25. POROSITY First the most important parameter of sphere packings is considered, packing density or volume fraction VV of the space occupied by the spheres for two important cases, lognormal and two-point distributions. The case of a two-point distribution is of particular interest since in the literature just such packings were studied experimentally for real spheres (Sohn and Moreland, 1968; Yu and Standish, 1987). Here one diameter is assumed to be 1 and the other is k. The probabilities of occurrence of these values are p1 = p and pk = 1 — p. Fig. 3 shows porosity h = 1 — VV as a function of k and v for packings obtained by the sedimentation algorithm, where v = p p + k3(1-p) (1) It is interesting to compare these results with the empirical results, which are also shown in Fig. 3. It is not surprising that the empirical porosities are clearly smaller, since the sedimentation algorithm does not include any compression and rearragement step. For the force-biased algorithm a similar figure is not presented, since it is possible to find program parameters which yield just the empirical values given in Fig. 3. In the lognormal case two parameters play a role: m and s 2, where the mean is given by m = exp(m + s2/2). It can be observed that volume fraction VV depends only little on the mean-value parameter m but heavily on the variance parameter s; VV increases with s. PAIR CORRELATION FUNCTION AND SPHERICAL CONTACT DISTRIBUTION The classical descriptor of spatial variability with respect to a homogeneous system of spheres with mean number l of spheres per volume unit is the pair correlation function g(r); (see Stoyan et al, 1995, p. 129). We have determined g(r) statistically by means of the method described in Stoyan and Stoyan (2001) for many of the simulated packings. Fig. 4. Pair correlation function for a dense random packing of spheres with equal diameters and VV=0.65. For the case of packings of spheres of constant diameter obtained by the force-biased algorithm the result is similar to the function shown in Fig. 4.10 of (Stoyan et al, 1995) if VV = 0.64, while in the case of VV = 0.70 the higher order is clearly expressed by peaks much sharper as those shown in Fig. 4. For two-point distributions the number of peaks is smaller, and they result only from pairs of small spheres. 204 Image Anal Stereol 2001;20:203-206 1.2 COORDINATION NUMBER g(r) 0. 0.2- 0.2 0.4 0.6 0.8 Fig. 5. Pair correlation function for the case of a lognormal diameter distribution and VV=0.70. Finally, Fig. 5 shows the pair correlation function in a typical lognormal case. The given mixture of spheres of different diameters leads to the effect that this curve does not have any striking peak; it looks like the pair correlation function of a so-called softcore point process. HJr Fig. 6. Spherical contact distribution functions. —— equal spheres, VV=0.64;--------equal spheres, VV =0.70;.....lognormal diameters. The spherical contact distribution function Hs(r) is a descriptor for sphere systems of a quite different nature; it describes the size of the pore space between the spheres; Hs(r) is the distribution function of the random distance from a random test point in the pore space to the nearest sphere surface point (Stoyan et al., 1995, p. 206). Fig. 6 shows Hs{r) for equal spheres with VV = 0.64 and 0.70 and for the lognormal case with VV = 0.70. The different porosities and diameter distributions lead of course to different spherical contact distribution functions, but also the different arrangement of the spheres plays a role. A very important feature of sphere packings are the contacts between spheres. They determine, for example, the topological connectivity of the system of spheres or the transfer of forces in mechanically loaded systems of hard spheres. In the literature there are reports on the distribution of the number of contacts or coordination number c of a randomly chosen sphere and its mean c. It is clear that these characteristics depend on the diameter distribution as well as on porosity. For real packings of equal spheres, in Bernal and Mason (1960) the value of c = 6.4 is given for a ‘random close’ packing and c = 5.5 for a ‘loose’ random packing, with VV = 0.62 and 0.60, respectively. For packings obtained by the sedimentation algorithm values around c = 6 were observed. Fig. 7 shows the histogram of contact numbers c for the case of equal spheres and force-biased algorithm with VV = 0.64 and 0.70. Clearly, in the denser packing there are more contacts. The same figure also shows the empirical contact number distribution for the case of lognormal diameters as in Fig. 1. U.J J 0.3 0.25 0.2 J— 0.15 ! i j - J ~ ~- - 0.1 1 1 ~ J i 0.05 1 1 —, Z—j 4 6 8 10 coordination number 1.2 14 Fig. 7. Distribution of the number of contacts of the spheres in three packings. ——equal spheres, VV = 0.64;-------equal spheres, ¦ ¦ lognormal diameters, VV = 0.70. VV = 0.70; REFERENCES Bernal DJ, Mason J (1960). Coordination of randomly packed spheres. Nature 188:910-1. Jodrey WS, Tory EM (1979). Simulation of random packing of spheres. J Simulation 32:1-12. Mos´cinski J, Bargie³ M, Rycerz ZA, Jacobs PWM (1989). The force-biased algorithm for the irregular close packing of equal hard spheres. Mol Simulat 3:201-12. 205 Bezrukov A et al: Spatial statistics Mos´cin´ski J, Bargie³ M (1991). C-language program for the irregular close packing of hard spheres. Computer Phys Comm 64:183-92. Sohn HY, Moreland C (1968). The effect of particle size distribution on packing density. Can J Chem Eng 46:162-7. Stoyan D, Kendall WS, Mecke J (1995). Stochastic geometry and its applications. Chichester: Wiley. Stoyan D, Stoyan H (2001). Improving ratio estimators of second order point process characteristics. Scand J Statist 27:641-56. Stroeven P, Stroeven M (2000). Assessment of particle packing characteristics at interfaces by SPACE system. Image Anal Stereol 19:85-90. Yu AB, Standish N (1987). Porosity calculations of multi-component mixtures of spherical particles. Powder Technol 52:233-41. 206