Geometrical properties of a random packing of hard spheres

Physica A 176 (1991) 206-219 North- Holland

GEOMETRICAL PROPERTIES OF A RANDOM PACKING OF HARD SPHERES A. PAVLOVITCH Section de Recherches de M(tallurgie Physique, Centre d'Etudes Nucl~aires de Saclay, B.P. 2, 91191 Gif-sur-Yvette, France

R. JULLIEN Laboratoire de Physique des Solides, Universit( Paris-Sud, Centre d'Orsay, 91405 Orsay, France

P. MEAKIN Central Research and Development Deparonent, E. l. du Pont de Nemours and Company, Wilmington, DE 19880-0356, USA

Received 30 April 1991

Some geometrical properties of a random packing of identical hard spheres generated by a ballistic deposition model with complete restructuring are il,vestigated. Thc length distribution of chords in the space between spheres is numerically calculated and is shown to have an exponential form (up to chord lengths of about five diameters) as conjectured by Dixmier. The anisotropic properties of the packing are numerically investigated and are shown to modify the Dixmier relation between the packing fraction and the average coordination number.

1. Introduction Random packings of hard disks and spheres have been the subject of cnnsiderable attention for many years, mainly due to their role in understanding the structure of amorphous, porous and random materials such as glasses [1] and also due to their importance in powder technology [2]. In this paper we investigate some geometrical p~operties of a three-dimensional packing of hard spheres obtained by a ballistic deposition model with restructuring which has recently been studied using large scale computer simulations [3]. In particular, we compute the length distribution of chords F(c) in the empty space of tilis random packing, because there were interesting predictions by Dixmier [4] in 0378-4371/91/$03.50 © 1991 - Elsevier Science Publishers B.V. (North-Holland)

A. Pavlovitch et al. / Geometrical properties of random packb+g of hard spheres

207

1978 based on such a distribution. Dixmier's results can be briefly summarized in two formulas. The first one (formula 24.2 in ref. [4]): F(c=0)=

3

d(.) S

'

(1)

is an exact result which relates the c = 0 limit of the chord distribution (properly normalized to unity), F(c = 0), to the average coordination number (it) thfough the geometrical characteristics of the spheres, i.e. their diameter d and their surface S = wd 2. The second one (formula 24.18 in ref. [4]): 9

or

( n ) = 2 1 - or'

(2)

which relates ( n ) to the packing fraction or, is based on a conjecture asserting that the length distribution of chords should have the exponential form (formula 24.14 in ref. [4]) v(¢)=

1

-c] ,

(3)

where (c) is the average chord length. It is worth noticing that formula (1) was rederived in 1985 by Gardner [5] under a slightly different form, up to a multiplicative factor of 4, which was later corrected by Oger et al. [6]. The precise identification between Gardner's and Dixmier's formulas is made explicit in appendix A. In this paper we extend these previous works by studying in our case the role of the anisotropy and its influence on the Dixmier conjecture. After a brief recall of the construction method in section 2, we present a detailed analysis of the distribution of chords, from an analytical point of view, in sectioh 3 and numerical results in section 4.

2. Construction of three-dimensional random packing The packing used here is built according to the "steepest descent rule +', first introduced by Visscher and Bolsterlii [7] and reinvented by two of us as an example of off-lattice ballistic deposition of particles [3] and used in order to discuss various physical phenomena such as segregation of polydisperse spheres [8] and penetration of small particles in such a random medium [9]. Such simulations represent real processes in which particles are slowly deposited under a small gravity field via a low coneentratior, incident flux without

208

A. Pavlovitch et al. / Geometrical properties of random packing of hard spheres

shaking. We do not take into account many particle effects, such as avalanches or many particle restructuring (a particle does not move as soon as it has become part of the deposit). The packing is built by a sequential release of N spheres of diameter d on a basal plane of size Ld x Ld, with periodic boundary conditions at edges of the square. The addition of a new sphere is carried out in the following way: (i) The moving particle (MP) is released at a random location above the basal square, higher than the packing already constructed. (ii) The MP falls vertically until it contacts either the basal plane (in which case it stops) or another sphere of center C 1 in the packing. (iii) The MP rolls vertically on the sphere of center C1, until it either contacts the basal plane and stops, either reaches an equatorial ~ position (in which case it takes off and we go to (i), or contacts another sphere of center

c2. (iv) The MP rolls on both spheres of centers C~ and C 2, until it either contacts the basal plane and stops, either reaches a "taking off" position at which it looses its contact with one of the two spheres to roll on the other (in which case we go to (iii)), or contacts a third sphere of center C a. (v) If the position is stable under gravity the MP stops and another sphere is deposited, if it is unstable we go to (iv). A detailed description of the procedure is presented in ref. [9]. Fig. 1 shows a typical packing obtained with L = 32. Such 3D packings are of completely different nature than those obtained in 2D with an analogous procedure [10]: here randomness results from geometrical frustration between hcp and fcc structures, since one cannot tile space with regular tetrahedrons. Let us first recall how, in such packings, the average number of contacts per particle (or average coordination number), denoted ( n ) , is strictly equal to 6, as already shown [11]. Let us call c(N) the total number of contacts and ~ ( N ) the number of spheres with i contacts in a packing made with N spheres. As we add a new sphere, by construction we add 3 new contacts. Hence c(N + 1) = c(N) + 3. This gives

c(N) = 3(N - 2 ) ,

for N > 2

with c(1) = 0 and c(2) = 1.

However, we still have E if~(N) = 2c(N), and, in the limit N---> :z, the average number of contacts per sphere is • t An equatorial position corresponds to a situation where the two sphere centers are at the same height.

A. Pavlovitch et al. / Geometrical properties of random packing of hard spheres

209

Fig. 1. A typical example of random packing built with our model. Here 32 708 opkeres of diameter d have been released vertically on a square horizontal basal plane with an edge length ol 32d with periodic boundary conditions at the lateral edges of the square.

1

(n) : ~ ~ if~(N):6. Consequently any packing of hard spheres that is constructed so that each addition of a new sphere brings three new contacts (i.e. there is no degeneracy as in, for example, the fcc packing) will have on average 6 contacts per spheres, despite the fact that the density can change when the construction method is changed. Some preliminary results with an alternative procedure, similar to the one introduced by Bennett [1], in which the new particle is added at the lowest position insuring three contacts with previously deposited particles, leads to a more compact structure for which ( n ) is also exactly 6 [12].

3. Chord method

In this part we would like to summarize the derivation of Dixmier's results [4]. The whole reasoning is based on the chord method, which was introduced by Dirac in order to calculate escape probabilities of neutrons in an absorbing medium [13]. This method has been used in a number of works on neutron absorption and scattering [14-16].

A. Pavlovitch et al. / Geometrical properties of random packing of hard spheres

210

Consider the shape (S) shown in fig. 2. If we select a line (L) of random orientation which crosses (S), the segment of (L) that is contained in (S) (see fig. 2) will be called a chord. Selecting a large number of such random lines, we will obtain a set of chords of different lengths and orientations. The shape (S) can be characterized by F(c), the length distribution of chords: F(c)dc is the probability that a chord has a length lying between c and c + dc. If dS is a surface element of (S) and n the corresponding normal, then the number of chords in the direction g~ is proportional to g~. n. From this we obtain the expression

~ dS fn=~ dg2 ~ " n F(c)dc=

ydSSdOa'n

(4)

"

where R denotes the length of the chord in the direction g~. The denominator can be integrated, using dO = 2xr sin 0 dO, where 0 is the angle between the chord and a fixed direction in space:

f dS f dO g~ " n = 2arS f cos O sin O dO ='trS,

(5)

0

where S denotes the total area of (S). The average chord length (c) is defined as

(c)=

f 0

F(c) c d c = ~ j 1

f dS f

Rdaa.n=4~,

v

(6)

R = c

$

where V is the volume of (S). Our interest will be devoted to length distribution of chords in the empty space of the packing. We can obtain an expansion of F(c) for small c using the

Fig. 2. A chord of length R on a line (L) in the shape (S).

A. Pavlovitch et al. / Geometrical properties of random packing of hard spheres

211

fact that the dominant contribution to F(c = O) comes from contacting pairs of spheres as made explicit in appendix B. This leads to the first Dixmier formula (1), which can be rewritten as

e(c = 0 ) = 3 d

(7)

Before discussing the case of a random packing let us recall the results for F(c) in the case of a crystalline or quasi-crystalline packing. Explicit calculations were performed by Benoist [15i ifi the case of a two-dimensional hexagonal crystal of hard disks. Consider a regular hexagonal two-dimensional crystal of hard disks. We define ~b = 2a/s, where a is the disc radius and s the lattice constant (see fig. 3). It can be easily seen that for ~b > ½x/-3 the crystal will not contain corridors of infinite length and that for ~ < ~b < ~X/'3 such corridors will be present in three different directions. If ~b < ½new corridors will be added [15]. This has direct consequences for the length distribution of chords; when no corridors are present, the chord length has an upper limit, but when corridors are present, the length distribution of chords will decrease as a power of the chord length c, for large c [15]. In the case of a packing of hard spheres which represents a random medium, the behavior of the length distribution of chords is completely different. Dixmer [4] assumed that, for an isotropic packing, the length distribution of chords should have the form (3), which can be used to get an estimate of the length distribution of chords when only the mean chord length is known (for example from the density and the surface of the shape). I

s

I

©6 3 © Fig. 3. Sketch of a two-dimensional hexagonal array of hard disks. Corridors in two dilcctions are shown.

212

A. Paviovitch et al. / Geometrical properties of random packing of hard spheres

Now using formula (6) to express the average chord length in the voids of a random packing whose packing fraction is tr, we obtain (c)

=

2d 1 - t r 3 o"

(8)

'

and combining this with formulas (3) and (7), we get the second Dixmier formula (2), relating the packing fraction to the average coordination number. This second Dixmier tormula is not appropriate for an anisotropic random packing. In our case, using the value of the packing fraction o-= 0.5815 +-0.0005 obtained from large scale simulations [2], we get ( n ) = 6 . 2 5 . . . instead of ( n ) =6.

4. Numerical results

It is obvious that our particular packing, whose construction involves gravity, is anisotropic: the vertical direction is not equivalent to any horizontal direction. A clear illustration of the anisotropy is provided in fig. 4 where we have plotted the distribution of angles between the lines passing through the centers of spheres in contact and the vertical. N(0)sin 0 dO is proportional to the number of bonds whose angle with the vertical lies between 0 and 0 + d0. The average angle corresponds to a value between 9 and ½"rr.

0"00601

o ooo1 0,00404

z t

/

0.0030

o

lll

I

& LII

,f

oooo,j 0.0000|'~ 0

"

q

i

I

i

t

~

I

I

f

~0

20

30

40

50

60

70

80

90

0 Fig. 4. Distribution N(0) of angles 0 between the vertical and the lines joining centers of spheres in contact.

A. Pavlovitch et ai. / Geometrical properties of random packing of hard spheres

213

0"801V

,.;0

0.;0

R Fig. 5. The radial distribution function g(R) obtained from a packing of 32 768 spheres.

A classical way to describe a packing of hard spheres is to calculate the radial correlation function g(R), which is the distribution of distances R between particle centers in the packing. We have performed this numerical calculation with our packing (with 32 768 spheres and taking into account the periodic boundary conditions) and the result is reported in fig. 5. After a g-peak at R = d corresponding to spheres in contact, the curve exhibits two other relatively strong peaks at R = dV'-3 and at R = 2d corresponding to configurations shown in fig. 6. It is worth comparing these with the first peaks in compact 3D lattices which are located at d-V~, dx/3 and 2d in the fcc lattice and dXF2, d V ~ , dx/'3, dVr~ and 2d in the hcp lattice. The behavior of the radial correlation function between d and 2d in the random packing appears to be quite particular: g(R) behaves almost like a second order polynomial of R for d ~< R ~< 2d. This is probably related to the fact that our random structure

Fig. 6. Sketch of the typical neighboring distances in a closed packed structure.

214

A. Pavlovitch et al. I Geometrical properties of random packing of hard spheres

2-

0-

o_

": ~-.,5 "-::-.. -10-12-14-16(a) -1B

0

I

I

I

I

I

I

1

2

3

4

5

6 cld

2-

O-

U. 0

', :~.-~._:_:=.... -10-12 i

-1 (b)

-1 0

I

I

I

I

I

I

I

I

1

2

3

4

5

6

7

8

old 2-~

0-2LI.

-6-8-

-10m

--I 4-

-" "~..~ "~. -....

- 1 6 -

-18

-

(C)

t

I

I

1

2

3

"

"

..

I

I

I

I

I

4

5

6

7

8

old

Fig. 7. Logarithm of the length distribution of chords F(c) plotted as a function of c/d. (a), (b) and (c) correspond to 0 = O, 0 = ~ ~" and the average over O, respectively.

A. Pavlovitch et al, / Geometrical properties of random packing of hard spheres

215

results from the geometrical frustration between fcc and hcp lattices as it would have been more clearly seen if we have started from a compact triangular array of spheres as basal plane. We have calculated the length distribution of chords F(c) in our packing using a square base with edges of length L = 32d and containing 32768 spheres. The results are shown in fig. 7 using a logarithmic scale. It can be seen that the deviations from an exponential are almost imperceptible. Since our packing is anisotrcpic, we have also studied such distributions for chords of fixed direction in space making an angle e with the vertical. Each distribution is normalized such that its integral is equal to 1. It appears that each individual distribution is also almost exponential. The dependence on O is best seen in fig. 8 where we have plotted ( c ) (numerically calculated by integration of cF(c)) and 1/F(O) as functions of 0. The difference between these two curves is a quantitative measure of the deviation of the chord distribution from a pure exponential form. Moreover, averaging F(0) over 0, one gets 1/ (F(0)) = 0.4962, to be compared with the exact theoretical value ½ deduced from formula (7), while averaging directly ( c ) , one gets (c) = 0.4776 to be compared with the theoretical value ( c ) = 0.4798 deduced from formula (8) with o" = 0.5815. The small discrepancy, in each case, comes from the errors made by numerical integrations, but the significant difference between (c) and 1 / ( F ( 0 ) ) shows that the second Dixmier formula does not hold here.

0.5200.5100.5000.4900,480-

iiiiii J 0.450

0

^/

,

~

,

10

20

30

40

50

60

70

80

90

0

Fig. 8. Plot of mean chord length ( c ) as a function of O. (a) corresponds to a direct calculation from an integral over R. (b) corresponds to an indirect estimation from 1 IF(O).

216

A. Pavlovitch et al. , Geometrical properties of random packing of hard spheres

5. Conclusion In this paper we have analyzed some geometrical properties of a random packing of hard spheres built according to a ballistic model with steepest descent rules. We have shown that, although the length distribution of chords within the voids of the packing fits quite well an exponential form, the Dixmier formula relating the density and the mean number of contacts is not valid in our case. We have shown the important role played by the two characteristic lengths of the distribution of chords, namely (c) and l / F ( 0 ) , in quantifying the anisotropic character of such packing. In future works we intend to study other packings built with different rules leading to different packing fractions but with always ~he same mean number of contacts between particles (i.e. 6). It would be interesting to try to build a random packing within the same constraint ( ( n ) = 6) with no anisotropy, but this appears to be quite a difficult task.

Acknowledgements One of us (A.P.) would like to thank P. Benoist and G. Martin for stimulating discussions.

Appendix A In this appendix we show that the Gardner formula [5] modified by Oger et ai. [6], which reads n¢ =

3 Q(c = O) "a"

R

'

(A.1)

is equivalent to the first Dixmier formula (1). In Gardner's formula, n c is the number of contacts per unit volume, R is the sphere radius and Q(c) is the density of chords counted per unit length of a line selected at random. Obviously Q(c) is proportional to our F(c) and we can determine the proportionality factor knowing that the integral S~Q(c) dc should be equal to the number of chords (of any length) per unit length of the random line. Thus considering a cylinder of axis the random line, of radius R and length unity, we can show that the number of chords is equa! to the number of spheres whose centers are in the cylinder, giving

A. Pavlovitch et al. / Geometrical properties of random packing of hard spheres

f Q(c) dc = ,rRZp,

217

(A.2)

0

where p is the number of spheres per unit volume. It follows that F(c) = ~

1

Q(c).

(A.3)

Observing that n c = ½ (n)p and using (A. I) we recover (I).

Appendix B In this appendix we make explicit the low-c expansion of F(c) (see also refs. [5, 6]) and we show that only spheres in contact contribute to the zeroth order term. Consider two spheres of radii R and R' (R > R') separated by a gap 6 as shown in fig. 9. Using (4) and performing the azimuthal integral around the O'M axis, we can write

1: F(c, 6)dc- arS

0=o~

2~R2 d c°s 0 c°s/3 2"rrd c°s/3 "

(B.1)

0=0

R

o

0' Fig. 9. This figure shows two spheres nearly in calculate F(c) in the limit c = 0.

contact and

introduces the notations used to

218

A. Pavlovitch et al. / Geometrical properties of random packing of hard spheres

where F(c, 8) is the length distribution of chord for a given 6 and S is the surface of the two spheres which bounds a chord. Considering the limit of small c and using cos/3 = h/c, this integral reduces to e=ern

F(c, 6) = 4IrR2 f S

h2

~-~ 0 dO.

(B.2)

0=0

Using simple trigonometry we get the formulas (R' + h) sin 0' = R sin 0 ,

(B.3a)

( R ' + h) cos 0' + R c o s 0 = R + R' + 8 .

(B.3b)

Taking the limit of small O and 0' and eliminating 0' between the two equations~ we get (B.4) Then, inserting this relation in (B.2) and performing the integral over 0 (with 0 2 = 2 p ( c - 6 ) / R 2 and 1/O = I / R + 1/R') we get

F(c, 8 ) - 4 x r p 1 s3 7

(B.5)

To calculate F(c) for small c in a packing of equal spheres, one should integrate F(c, 6) over 6 taking into account the form of the radial distribution function near R = 2d: ~=¢

F(c)= f

[ ( n ) ~ ( 8 ) + A + O(8)lF(c, 6)d6 ,

(B.6)

8=0

where 6 denotes the Dirac distribution and A is a constant. After a straightforward integration one gets the first Dixmier formula (1) as the zeroth order leading term. It is worth noticing that only spheres in contact contribute here while nearly contacting spheres only ~,O.u.t, . . . . :'- .... utc to hiDher order terms.

References [1] J.D. Bernal, ~,oc. R. Soc. London A280 (1964) 299. C.H. Bennett, J. Appl. Phys. 43 (1972) 2727. M. Rubinstein and D.R. Nelson, Phys. Rev. B 2,7, (1982) 6254.

A. Pavlovitch et al. / Geometrical properties of random packing of hard spheres

[2] [3] [4] [5] [6] [7] [8] [9] [10] [11]

219

J.C. Williams, Powder Technol. 15 (1976) 245. R. Jullien and P. Meakin, Europhys. Lett. 4 (1987) 1385. M. Dixmier, J. Phys. (Paris) 39 (1978) 873. E. Gardner, J. Phys. (Paris) 46 (1985) 1655. L. Oger, M. Lic!ltenberg, A. Gervois and E. Guyon, J. Microsc. 156 (1989) 65. W.H. Visscher and H. Bosterlii, Nature 239 (1972) 504. R. Jullien and P. Meakin, Europhys. Lett. 6 (1988) 629; Nature 344 (1990) 425. P. Meakin and R. Jullien, J. Phys. (Paris) 51 (1990) 2673. R. Jullien and P. Meakia, Europhys. Lett. 14 (1991) 667. L. Oger, J.P. Troadec, D. Bideau, J.A. Dodds and M.J. Powell, Powder Teehnol. 46 (1986) 121. [12] R. Jullien an' P. Meakin, in preparation. [13] P.A.M. Dirac, ,)eclassified British Report M.S.D. 5 (1943) part I. [14] D.J. Behrens, Proc. R. Soc. A 62 (1949) 607. J. Lieberoth and A. Stojanovic, Nucl. Sci. Eng. 76 (1980) 336. [15] P. Benoist, Nucl. Sci. Eng. 86 (1984) 2240. [16] M.Y. Lin and S.K. Sinha, Math. Res. Soc. Symp. Proc. 195 (1990) 485.