Structure of sticky-hard-spheres random-aggregates: the viewpoint of contact coordination and tetrahedra
SStructure of sticky-hard-spheres random-aggregates: the viewpoint of contactcoordination and tetrahedra
Bl´etry, M., ∗ Russier, V., and Barb´e, E.
ICMPE-CNRS, 2-8 rue Henri Dunant, 94320 Thiais, France
Bl´etry, J
Honorary Professor, Bahia Blanca University, Argentina (Dated:)We study more than 10 random aggregates of 10 monodisperse sticky hard sphereseach, generated by various static algorithms. Their packing fraction varies from 0.370 upto 0.593. These aggregates are shown to be based on two types of disordered structures:random regular polytetrahedra and random aggregates, the former giving rise to δ peaks onpair distribution functions. Distortion of structural (Delaunay) tetrahedra is studied by twoparameters, which show some similarities and some differences in terms of overall tenden-cies. Isotropy of aggregates is characterized by the nematic order parameter. The overallstructure is then studied by distinguishing spheres in function of their contact coordinationnumber (CCN). Distributions of average CCN around spheres of a given CCN value showtrends that depend on packing fraction and building algorithms. The radial dependencyof the average CCN turns out to be dependent upon the CCN of the central sphere andshows discontinuities that resemble those of the pair distribution function. Moreover, it isshown that structural details appear when the CCN is used as pseudo chemical parame-ter, such as various angular distribution of bond angles, partial pair distribution functions,Ashcroft-Langreth and Bhatia-Thornton partial structure factors. These allow distinguish-ing aggregates with the same values of packing fraction or average tetrahedral distortionor even similar global pair distribution function, indicative of the great interest of payingattention to contact coordination numbers to study more precisely the structure of randomaggregates. Keywords: Sticky hard spheres ; random packings ; contact coordination number ; Delaunay tesse-lation ; tetrahedra distorsion ; distribution function of isocoordinated spheres ; Ashcroft-Langrethpartial structure factor ; Bhatia-Thornton partial structure factor ∗ Corresponding author [email protected] a r X i v : . [ c ond - m a t . d i s - nn ] S e p I. INTRODUCTION
Random aggregates of monodisperse spheres are of great interest to simulate various physicalsystems, such as amorphous solids, liquids, powders, etc. In particular, a fine characterization oftheir local or longer range structural properties is of interest to classify the various types of randomaggregates of spheres that can exist and study their properties.Many different approaches allow building random aggregates of spheres. Broadly, they maybe divided into two categories. The first one consists of algorithms for which all spheres areintroduced at once and then the system relaxes towards some more or less disordered structuralstate (eg molecular dynamics [1], Lubachevsky-Stillinger algorithm [2], Jodrey-Tory algorithm [3] orMonte-Carlo relaxation of chains of hard-sphere [4]). This family of algorithms produces aggregateswhose properties mimic what was found in many disordered systems, notably the dependency ofthe contact coordination number with packing fraction and they are also able to produce aggregateswith packing fraction equal or even superior to the random close packing (RCP) value ( γ RCP ≈ .
64) or equal to or lower than the random loose packing value (RLP, γ RLP ≈ .
555 see [5]).However, γ RLP presents a strong dependency on interparticles interactions [6] which suggests thatno geometrical property must be associated to the RLP transition, whereas the RCP state mighthave structural or geometrical constraints that set an upper bound to γ RCP , such as the proportionof spheres involved in quasi perfect tetrahedra, as observed by Anikeenko and Medvedev [7].The second large category of algorithms introduces spheres in the aggregate one by one andsets their position definitely at once (see, e.g., [8]). Aggregates built by such methods have twoimportant differences with those of the first category: first, their average contact coordinationnumber remains close to 6, second, to the best of our knowledge, it is impossible to produceaggregates with packing fraction higher than 0.6. Hence, the latter family of algorithms gives accessto random aggregates with somewhat different structural properties from the former. Moreover,using sticky hard spheres, contact neighbours are rigorously defined by Dirac δ functions and thisallows studying various properties related to the coordination number of spheres in a simple andnon ambiguous manner.In [9], a first study of several hundreds sequentially built aggregates was conducted. It wasfound that for the lowest packing fractions, δ peaks appear on the pair distribution function, whichcorresponds to the formation of a disordered polytetrahedral structure in the random packingaggregate.The present study extends the results obtained in [9] by looking at structural tetrahedra,isotropy, contact coordination numbers, partial distribution functions as well as partial structurefactors, corresponding to pairs of spheres with equal or different contact coordination numbers.It also introduces a new family of random aggregates formed by sphere added in regular buildingtetrahedra, which appears essential to understand the composite nature of other sequentially builtrandom aggregates. II. STUDIED AGGREGATES
More than 10000 aggregates, each one containing 10 spheres with radius r s = 1 (and, hence,diameter d = 2) were built and studied[10]. They fall into two broad families: random irregularpolytetrahedral aggregates (RIPA) and random regular polytetrahedral aggregates (RRPA). All ofthem were built by adding spheres one by one at their final position, tangentially to three existingones. A. Random irregular polytetrahedral aggregates
Most aggregates were of RIPA type. They were built using algorithms that have been presentedin details in [9] and are only summarized here. A seed of three spheres forming an equilateraltriangle is used. Each new sphere P is positioned tangentially to at least three already presentspheres (noted O , A and B ). The new sphere can be introduced in a hole whose size is maximized(MAX algorithms) or chosen randomly (RAN algorithms) in the vicinity of the local origin O ,itself chosen randomly in the aggregate. Moreover, it is possible to insert from 1 up to 9 spheresat once (according to the index N IN S , for
Number of sphere INSerted ) around a given origin O (algorithms MAX 1 to 9 and RAN 1 to 9). Finally, the neighbourhood explored to choose spheres A and B around O is a cube whose edge length a can be varied: it controls strongly the packingfraction of the aggregates. The larger the value of a , the higher the packing fraction. It is variedbetween about 3.4 and 8 as, for a < .
4, no aggregate can be generated and for a >
8, the maximumpacking fraction is reached and no evolution of the generated aggregates is noted for higher valuesof a (see figure 2 in [9]). More than 300 aggregates of 10 spheres were built, by varying a for eachfamily of these algorithms.An additional modification with regard to the aggregates studied in [9] was to choose the origin O as close as possible from the center (0,0,0) of the growing aggregate (RMIN algorithms), insteadof a purely random fashion. This change has the effect of increasing the maximum packing fraction,reached for (MAX-1, a > . a and N IN S ) with no significant changes concerning the structural resultspresented in [9]. This modification also allows for a more homogeneous growth of the aggregate.Only the results obtained for RMIN-MAX-1 will be used hereafter[11].Finally, as it turned out that the previous families of aggregates had some inhomogeneities oftheir packing fraction close to (0,0,0) (see hereafter), some aggregates were generated by usingas seed a set of N spheres taken in previously built aggregates, instead of an equilateral triangle(N-RMIN-MAX-1 algorithm). The positions of the spheres composing the seed are taken froman aggregate with the same value of other parameters ( a and N IN S ), far from the origin (0,0,0),which has the effect to remove the central area with higher packing fraction. Typically, the numberof spheres in the seed is between 30 and 600. However, this modification entails only slight changesof pair distribution functions or structure factors.
B. Random regular polytetrahedral aggregates
A last algorithm has been used, that produces aggregates with only regular building tetrahedra[12]. In this case, the newly inserted sphere P forms a regular tetrahedron with the three alreadycontacting spheres O , A and B , i.e. P O = P A = P B = OA = OB = AB = 2. For thesealgorithms, a has virtually no impact on packing fraction, but the number of inserted spheresaround a given local origin (NINS) does.Once again, the local origin O can be chosen at random, or as the closest one from the centerof the aggregate. When it is chosen randomly (RRPA), the maximum packing fraction (0.418) isreached for N IN S = 3 and a minimum of 0.408 is observed for
N IN S = 1 and beyond 4. When O is chosen as the closest possible origin from the aggregate center (RMIN-RRPA), the maximumpacking fraction is 0.456 and is reached for N IN S >
4. Similar aggregates, so-called saturatedpolytetrahedra, have been studied by Medvedev and Pilyugina [13]. They found a packing fractionof 0.435 for aggregates consisting of about 576 000 spheres. This value falls in between the maximumones obtained for aggregates choosing O randomly (lower bound) and those taking O as the closestsphere from the origin (higher bound).It should be noted that RIPA and RRPA distinguishes aggregates based on their buildingalgorithms, not their structure, which will be discussed in details below.The isotropy and the randomness of all aggregates has been systematically checked through thedistribution of i − j bonds and the nematic tensor formalism, and turns out to be satisfactory.More details on that latter point are provided in supplementary informations. III. PACKING FRACTIONA. Basic relations
The radius R m of a large spherical aggregate centered in (0,0,0) and made of N spheres centeredin (cid:126)R i is given to a good approximation by [12]: R ag = (cid:113) R (1)where R quad is the average quadratic radius of all spheres in the aggregate: R = 53 1 N N (cid:88) i =1 R i (2)However, finite aggregates are not fully spherical and exhibit local order oscillations. Therefore,their packing fraction varies as a function of the radius R SC of the sphere cut into the aggregatebulk and deserves special attention.The volume V s ( r ) shared by a sphere of radius r s = 1 whose center is at a distance r from theorigin (0 , ,
0) with another sphere S , of radius R SC , centered in (0 , ,
0) is [14]: V s ( r, R SC ) = 43 πr s r ≤ R SC − r s V s ( r, R SC ) = π ( R SC + r s − r ) ( r + 2 rr s − r s + 2 rR SC + 6 R SC r s − R )12 r R SC − r s < r < R SC + r s (3) V s ( r, R SC ) = 0 r ≥ R SC + r s These relationships can also be used for any sphere S centered in x , y and z by a mere change ofreference frame.Hence, the packing fraction of the sphere S can be directly determined for any radius R SC , as: γ ( R SC ) = (cid:80) Ni =1 V s ( r i , R SC ) πR (4)where i accounts for all spheres in the aggregate.Finally, it is possible to determine the packing fraction of shells of arbitrary thickness w = R o − R i , where R i is the inner radius and R o , the outer radius of the shell, simply by removing theportion of spheres outside of the shell, according to the volume complementary of relation 3. B. Packing fraction of spheres inscribed in the aggregate
The evolution of the packing fraction of spheres inscribed in the aggregate as a function of theirradius γ = f ( R SC /R ag ) (figure 1) shows that, for all aggregates, whatever their building algorithm,the packing fraction decreases slightly when R SC increases.A seed effect appears as for MAX-1 and RMIN-MAX-1 aggregates (3 spheres forming an equi-lateral triangle as seed), a first regime of fast decrease is observed for R SC < . R ag , i.e. for anumber of spheres below approximately 8000, then the packing fraction tends to plateau whereasfor N-RMIN-MAX-1 aggregates (seed consisting of spheres taken far from the origin in a previouslygenerated aggregate), this initial decrease is much faster. Nevertheless, the exact range of effect ofthe seed can only be asserted by the packing fraction of shells studied hereafter.This effect probably stems from the fact that contacting equilateral-triangle configurations areextremely rare in high packing fraction aggregates. As a matter of fact, this seed dependencydisappears for lower packing fraction aggregates, in which such configurations are rather frequent.On the other hand, N-RMIN-MAX-1 and RMIN-MAX-1 aggregates show exactly the samebehaviour for larger values of R SC . (a) (b) FIG. 1. γ ( R SC /R ag ) for several a) MAX-1 aggregates and b) RMIN-MAX-1 aggregates built with var-ious values of a and a single N-RMIN-MAX-1 aggregate (seed composed of 600 spheres), for the sake ofcomparison, with S centered in 0 , , For every aggregate, a second regime is observed when R SC → R ag , i.e. when S reaches thelimit of the aggregate: logically, the packing fraction decreases faster. For a perfectly sphericalaggregate, when R SC > R ag + r s , then the ”packing fraction” of S should decrease as R − . Whenthis decrease is at first more progressive, it shows that the aggregate has an imperfect shape andeither has protuberances on its surface or is not overall perfectly spherical.For the same value of a , RMIN-MAX-1 aggregates tend to have a higher packing fraction thanthe corresponding MAX-1 aggregates, as well as a sharper decrease of γ when R SC → R ag , whichshows that RMIN-MAX-1 aggregates have a more regular surface than MAX-1 aggregates. Forthe latter, the thickness of ”imperfect aggregate” is about 1 . d for a = 3 . γ = 0 . d for a = 1 .
79 ( γ = 0 . d/ a = 3 . γ = 0 . . d for a = 1 .
78 ( γ = 0 . C. Packing fraction of shells (a) (b)
FIG. 2. Packing fraction of shells with thickness w = 0 . r s on the densest aggregates built with threestrategies: MAX-1 (3 spheres seed), RMIN-MAX-1 (3 spheres seed, RMIN) and N-RMIN-MAX-1 with aspherical seed of radius 10 extracted from a previously generated aggregate with the same value for a and N IN S a) for relatively small r ( r/R ag ∈ [0; 0 . r ( r/R ag ∈ [0 .
5; 1]).
Figures 2.a and b represent the variation of packing fraction in shells with thickness w = 0 . r . Its average value decreasesfrom a higher value near the seed, to a smoother behaviour when r increases. N-RMIN-MAX-1aggregate presents virtually no effect of the seed: the packing fraction of shells reaches the averagebehaviour for very small values of r , which seems logical as for these aggregates the seed consistsof a set of spheres with the average structure. RMIN-MAX-1 and N-RMIN-MAX-1 converge for r ∈ [0 . R ag , . R ag ]: the effect of the initial equilateral-triangle seed of RMIN-MAX-1 aggregatesseems then to act on about 10 spheres in the whole aggregate, consisting of 10 spheres, i.e.significantly less than suggested above by the comparison, in figure 1.b, of γ = f ( R SC /R ag ) forthe two same aggregates. Figure 3 compares two aggregates built by RMIN-MAX-1 and N-RMIN-MAX-1 algorithms, where spheres are colored based on their CCN. It appears that the formeraggregate has a brighter atypical central area in the region of the seed, denoting unusual structuralproperties with respect to the rest of the aggregate, whereas the latter displays a seed area muchmore similar to the rest of the aggregate.Moreover, for aggregates with lower packing fraction (i.e. for aggregates built using lowervalues of a and in which regular polytetrahedra appear), the range of aggregate affected by theseed decreases and completely disappears for the lowest packing fraction aggregates. In that case,indeed, the structure contains a significant amount of equilateral triangles and the initial seedceases to be special in comparison with the rest of the structure. At large values of r (figure 2.b),oscillations can still be detected in the packing fraction of shells, however with a much smalleramplitude. A slight decrease of the average value is noticeable: the farther a shell is from thecenter of the aggregate, the lower its packing fraction, on average. The origin of this phenomenonis not obvious to us. For r > R ag , the packing fraction of shells falls rapidly to 0. (a) (b) FIG. 3. Slices of aggregates generated by two algorithms: a) RMIN-MAX-1 with a = 4, b) N-RMIN-MAX-1with a = 4 and a seed consisting of N = 400 spheres. Colors (gray scale) correspond to contact coordinationnumber: brighter spheres have a higher CCN. (These figures were generated with Ovito software [15].) IV. TETRAHEDRAL STRUCTURE
The structure of random packings of spheres is commonly assessed via the tetrahedra connectingsphere centers, forming the so-called Delaunay tesselation [16] (these tetrahedra are noted by thesubscript D in what follows). For the present study, Delaunay tessellations were built using thecgal library [17, 18].In [9], another type of tetrahedra was studied, called building tetrahedra (noted by the subscript BT hereafter). A building tetrahedron is formed by spheres O , A , B and P when adding the newsphere P tangentially to the three other ones, O , A and B . It should be noted that such tetrahedramay or may not belong to Delaunay’s tessellation. The distortion of building tetrahedra was shownto be a very significant structural parameter, allowing the correlation of various structural traitsof the aggregates. In this section, we focus on two distortion parameters of Delaunay tetrahedra. A. Distortion parameters
The first tetrahedral distortion parameter has been defined for the characterization of buildingtetrahedra [9] by relation: κ BT = 3 d + OA + OB + AB d (5)where O , A and B are the three sphere centers used to add the new sphere P and the term 3 d corresponds to the three necessary sphere contacts P O, P A, P B imposed to building tetrahedraby the algorithm. The maximum value of κ BT is 2 and is obtained for a centered equilateraltriangle with side d √ P at its barycenter, while κ BT minimum value, 1, corresponds to a regular tetrahedron.The definition of the parameter κ BT is immediately extended to Delaunay tetrahedra by relation: κ D = (cid:88) i (cid:88) j>i d ij d (6)where i and j are the vertices of the tetrahedron and d ij the vertices length. The smallest distancepossible between sphere centers is d ij = d = 2, hence the smallest possible value is obtained for aregular tetrahedron and is κ D = 1.The last distortion parameter to be studied hereafter is the longest edge length ( L max ) of theconsidered tetrahedron. The smallest possible value of L max , is d , which is found in the case of aregular tetrahedron. The behaviour of L max has already been studied, along with others, notably0by Anikeenko et al [19] on aggregates built using Jodrey-Tory (JT) algorithm [3] and Lubachevsky-Stillinger (LS) algorithm [2, 20]. Anikeenko et al [19] have found that L max , in spite of its simplicity,shows a great consistency when compared with two other parameters, namely the edge differencesand the procrustean distance. B. Distributions of distortion parameters (a) (b)
FIG. 4. Normalized distributions of a) κ D and b) L max for various packing fractions. Aggregates weregenerated with MAX-1 algorithm, with the exception of the RRPA. Error bars are smaller than point size. Globally, the distribution of distortion parameters suggest that two limiting aggregates exist,the densest one, produced by (RMIN)-MAX-1 algorithms (here for γ = 0 . κ D obtained for aggregates built using MAX-1 al-gorithm (RIPA) and one RRPA, and figure 4.b represents the distribution of L max for the same1 (a) (b) FIG. 5. Slices of aggregates generated by two algorithms: a) RMIN-MAX-1 with a = 4, b) RRPA. Colors(gray scale) correspond to contact coordination number: brighter spheres have a higher CCN. (These figureswere generated with Ovito software [15].) aggregates. Concerning MAX-1, these distributions present very similar behaviours: their maxi-mum decreases with packing fraction while their full width at half maximum (FWHM) increaseswhen packing fraction decreases. A bimodal component appears for the lowest packing fractionaggregates on both distributions, respectively centered around κ D ≈ L max ≈ .
5. These twovalues may be related in this way: assuming that tetrahedra with L max ≈ . . κ D value of approximately 1 . δ peaks appear for thelowest packing fraction, which is consistent with the existence of well defined distances observed onpair distribution function for the same aggregates, noting the existence of regular polytetrahedraand recurrent configurations of spheres in the structure (see [9]). The RRPA, on the other hand,appears as a limiting case. Indeed, the first peak observed for κ D and L max distributions in the caseof MAX aggregates completely disappears and is replaced by a series of δ peaks (some out of scale)and a slow evolution with respect to the second mode of the distributions of MAX aggregates.Moreover, the distributions of κ D present a discontinuity at κ D ≈ .
33, and the distributionsof L max have one at L max = 2 . L max = 4. Additional work is needed toanalyse precisely the configurations corresponding to these discontinuities. C. Average values of distortion parameters (a) (b)
FIG. 6. a) Dependency of ¯ κ D with packing fraction b) dependency of ¯ L max with packing fraction. Errorbars are smaller than point size. Although average distortion parameters are basic structural parameters, they do not deter-mine the packing fraction of the aggregate because they only involve the short range order ofspheres (through the tetrahedral description) but do not take into account longer range order. It istherefore worthwhile studying the variation of packing fraction with distortion parameters for allaggregates. Overall, on average, the higher the packing fraction, the less distorted the tetrahedra,which can appear as slightly counterintuitive as the lowest packing fraction aggregates have a highRP component in their structure and as RP means a regular tetrahedral basis.Figure 6.a presents the dependency of ¯ κ D with packing fraction, which decreases when γ in-creases and turns out to be linear for each family of algorithm. More specifically, the variousalgorithms are roughly distinguished as, for MAX- i aggregates, the average distortion increaseswith i for a given packing fraction. RMIN-MAX-1 stands on its own, with a slightly differentslope. The same is observed for RAN- i aggregates, which span a larger interval of ¯ κ D for a givenpacking fraction (approximately 3 times as wide as that of MAX aggregates) and a certain overlapis observed as MAX-1 to 4 are between RAN-4 and RAN-6 for γ < .
5, which is the maximumpacking fraction that RAN aggregates can reach.3The behaviour of ¯ L max = f ( γ ) is globally the same (figure 6.b) as it also decreases when packingfraction increases. However, its behaviour deviates more from linearity: in the case of MAX-1 andRMIN-MAX-1 aggregates a change of slope is observed between 0.47 and 0.5 and they behavevery differently from the rest of the other aggregates, with a higher value of ¯ L max than any otheraggregate for a given packing fraction. Furthermore, there is no overlap between RAN and MAXaggregates. RAN aggregates appear more dispersed than MAX aggregates and, for the latter, theymore or less converge on the same curve (that of MAX-9) with the noticeable exception of MAX-1aggregates.Hence, surprisingly, the average distortions measured by both indicators do not agree as, forexample, in the case of ¯ κ D , RMIN-MAX-1 aggregates appear as the least distorted and as the mostdistorted according to ¯ L max . D. Proportion of regular and quasi regular Delaunay tetrahedra
Using either distortion parameter, it is possible to evaluate the volume fraction of regularDelaunay tetrahedra (Φ V , with κ D = 1 or L max = 2), which globally decreases when packingfraction increases, as the RP component of aggregates decreases also.RRPA algorithms, in particular, give the highest Φ V , ranging from 0.101 for γ = 0 .
452 to 0.073for γ = 0 . V with packing fraction. The highest proportion is obtained for RAN-6 aggregates(Φ V = 0.141, γ = 0 . γ ∈ [0 .
45; 0 .
48] and outside of this interval, RMIN-MAX-1 outside of thisinterval. The proportion of perfect tetrahedra goes to 0 at the highest packing fractions.Anikeenko and Medvedev [7] have studied the volume fraction of quasi regular Delaunay tetra-hedra (PQRT, i.e. with L max < .
3) within aggregates generated with JT and LS algorithms (seefigure 7.b of the present article). They have found that this volume proportion increases with pack-ing fraction in their studied interval of packing fraction, i.e. roughly between 0.53 and 0.71. Figure7.b superimposes their results with the ones found in the present study. Interestingly, the curvesfor MAX-4 and MAX-3 match over a rather narrow packing fraction interval, around γ = 0 .
55, andall MAX- i and RMIN-MAX aggregates show a similar increase of the volume fraction of tetrahedrawith L max < . γ ∈ [0 .
56; 0 . V decreases when packing fraction increases and goes to 0 beyond a threshold that depends onthe algorithm but is roughly γ = 0 .
57, whereas PQRT goes to a minimum for a given packingfraction that depends on the algorithm (between γ = 0 .
55 and γ = 0 .
57) and then increases withpacking fraction beyond this value of γ . (a) (b) FIG. 7. a) Dependency of the volume fraction, Φ V , of perfect tetrahedra ( L max = 2) with packing fraction.b) Dependency of volume fraction, Φ V , of tetrahedra with L max < . V. EFFECT OF SPHERE COORDINATION
The contact coordination number (CCN) of each sphere is determined when the aggregate isbuilt: when a sphere is added in contact with another one, both their CCN are increased by 1.Hence, by the end of the building process, each sphere is associated to its CCN.
A. Contact coordination number
1. Partial distributions of contact coordination numbers
Let η ij be the number of contacts between spheres with CCN i and j respectively. The (normal-ized) distributions of η ij can be easily determined from the sphere positions for all values of i and5 j . Figure 8 introduces distributions of η ij for three aggregates built by RMIN-MAX-1 algorithms,from the highest to the lowest packing fraction (fig. 8.a to c) and a RRPA aggregate (with thehighest packing fraction among RRPA, fig. 8.d). (a) (b)(c) (d) FIG. 8. Normalized distributions of η ij (lines are mere guides for the eye) for aggregates produced byRMIN-MAX-1 algorithm with a) γ = 0 . γ = 0 .
538 and c) γ = 0 .
378 and by d) RRPA with γ = 0 . These distributions show a progressive shift of the maxima from high to low packing fraction(fig. 8.a to c). At high packing fraction, spheres with low CCN tend to be surrounded by sphereswith higher CCN thus reducing local fluctuations of the CCN. Then, as the packing fraction de-creases, all η ij curves more or less collapse, meaning that in this regime, all spheres–whatever theircontact coordination number η i –have the same η ij distribution, centered on the average contactcoordination number, hence a very similar environment in terms of contact neighbours. Finally,6for even lower packing fraction, an inversion is observed and high CCN spheres are preferentiallysurrounded by high CCN spheres, which corresponds to a contact segregation effect. At the sametime, the FWHM of the distributions widen as the packing fraction decreases and they becomeless symmetrical. For the lowest packing fraction (fig. 8.c) η ij distributions are very highly spread,with still a higher proportion of high CCN spheres in the vicinity of other high CCN spheres.The RRPA presents rather similar η ij distributions (fig. 8.d) as low packing fraction RIPAs,going through a minimum for η i . However, the order of the various η ij distributions appearsinverted in the case of the RRPA as, on the high j end, the spheres with the highest proportionof contact with high CCN spheres are spheres with lower CCN ( i ), i.e. η , > η , > ... > η , ,whereas for RIPA with low packing fraction, the opposite situation is observed, i.e. η , > η , >... > η , (fig. 8.c).
2. Evolution of (cid:104) η ij (cid:105) around i sphere (cid:104) η ij (cid:105) , the average value of η ij around spheres with CCN i , can easily be determined from thedistributions presented above. Figure 9 presents the evolution of (cid:104) η ij (cid:105) for all values of i as afunction of packing fraction for MAX-1 and RMIN-MAX-1 algorithms. The inversion suggested bythe shift of the maximum of distributions in the previous section appears clearly for RMIN-MAX-1aggregates (figure 9.b) as they form a crossover for γ ≈ .
52 which separates a low and a highpacking fraction regimes. In the low packing fraction regime, high coordination spheres tend tobe surrounded by spheres with higher CCN than the spheres surrounding low CCN sphere i.e. (cid:104) η ,j (cid:105) > (cid:104) η ,j (cid:105) > ... > (cid:104) η ,j (cid:105) , with the exception of the limit case (cid:104) η ,j (cid:105) , corresponding tothe segregation effect seen when discussing η ij distributions. In the high packing fraction regime,this situation is inverted: spheres with low CCN are surrounded – on average – by spheres withhigher CCN, thus reducing density fluctuations.In the case of aggregates generated by MAX-1 algorithm, this crossover is not captured butmight take place at higher values of packing fraction – unaccessible by this algorithm – as all curvesbegin to collapse for γ > .
55. In this case, only the low packing fraction regime is observed.
3. Radial dependency of the average contact coordination number
The radial dependency of the average contact coordination number (cid:104)
CCN (cid:105) of spheres within[ r ; r + dr ] from an i coordinated sphere has been determined for all aggregates and all values of7 (a) (b) FIG. 9. (cid:104) η ij (cid:105) for various values of i in the case of a) MAX-1 and b) RMIN-MAX-1 algorithms. i . Figure 10 introduces two examples obtained for the RMIN-MAX-1 aggregates with the highestand lowest packing fractions (figure 10.a and b). (a) (b) FIG. 10. Radial value of the average CCN around i coordinated spheres in RMIN-MAX-1 aggregates. a)RMIN-MAX-1, γ = 0 . γ = 0 . At high packing fraction (see figure 10.a), the average value of the CCN of spheres at a distance r < . r < . γ ≈ .
54, which incidentally matches thepacking fraction of the crossover of the various (cid:104) η ij (cid:105) curves in figure 9.b for the same aggregates.At even lower packing fraction, the inversion is complete as it is exemplified in figure 10.b: low8coordinated sphere are, on average, surrounded by quasi-first-neighbour spheres with low CCN andvice versa. This regime corresponds to a segregation effect in the range of quasi first neighboursinstead of contact neighbours. At even larger values of r ( r > r = r e , however r e increaseswhen packing fraction decreases (from r e ≈ γ = 0 .
593 up to r e ≈ . γ = 0 . γ decreases.Discontinuities are observed at r = d √ r = 2 d , which matches discontinuities of thevarious PDF of the same aggregates (cf. section V C 1 b below and [9]). The amplitude of thesediscontinuities increases with the CCN of the central sphere.At low packing fraction (fig 9.b), strong local increases of (cid:104) CCN (cid:105) are noticeable for valuesof r which matches δ -peaks on the PDF, hence corresponding to distances characteristics of thepresence of regular polytetrahedra in the aggregate (cf. section V C 1 b below). B. Bond angle distributions
Two spheres i and j form a bond when they are contact neighbours. The bond angles α arounda sphere i are defined as the angles formed between all possible vectors (cid:126)R ij between contactingneighbours. Hence, for a sphere i with n contact neighbours there are n ( n − /
1. Global bond angle distributions
Bond angle distributions α have been calculated. The smallest possible bond angle is betweenthree contacting spheres, i.e. α = π/
3, and the largest is, of course, π .Figure 11.a presents bond angle distributions for the most and least dense aggregates obtainedby MAX-1 algorithm. It shows that the distribution of the densest aggregate is mostly smooth,marked by two discontinuities: the first one, the minimum value, is α = π/
3: it corresponds tothe configuration of three spheres in contact with each other. The value here is so high that itis outside the scale of the figure. The second discontinuity occurs for α = 2 . ± . π/
3, which corresponds to the situation where four spheres form two coplanar equilateral trianglessharing a common side. In [21] Karayiannis et al obtained a very similar bond distribution in theirstructures of chains of joined monodisperse hard spheres in the MRJ (maximally random jammed)state. The main difference is the broader shape of the peaks in the distributions of ref [21], which9we interpret as a finite size broadening since the maximum number of spheres considered in [21] is54000 instead of 10 used here. (a) (b) FIG. 11. a) Bond angle distributions for high and low density aggregates (values for α = π/ δ -peak in the bond angle distribution with the irregularityindex of building tetrahedra ¯ κ BT . For the aggregate with the smallest density, many singularities appear, in the same way as δ -peaks appear on its pair distribution function. On PDF, these δ -peaks are characteristic of thepresence of regular polytetrahedra in the disordered structure (see [9]). The addition of anothersphere on top of three contacting ones does not introduce a new bond-angle. The addition of afifth sphere forms then a trigonal bi-pyramid, with α ≈ . δ -peak at r = d (cid:112) / P ( r = d (cid:112) /
2. Partial bond angle distributions
Bond angle distributions for specific contact coordination numbers are shown in figure 12.Their global behaviour is similar to global bond angle distributions, but peculiarities can be seen,in function of packing fraction and/or contact coordination number.It turns out that for aggregates with high packing fraction, the angular environment dependsstrongly on the CCN. Low CCN spheres have a stronger asymmetry between α = π/ α = 2 π/ (a) (b)(c) (d) FIG. 12. Partial bond angle distributions for various MAX-1 aggregates a) γ = 0 .
586 b) γ = 0 .
500 c) γ = 0 .
370 and d) one RRPA aggregate, γ = 0 . α = π/ show a depletion of lower bond angles and an excess of high angles and the distributions get moreeven when the CCN increases. When the packing fraction decreases, low and high CCN spherestend to have more similar angular environment.The appearance of regular polytetrahedra, associated with δ peaks, like in the case of the globalbond angle distribution, is logically correlated with a decrease of the continuum component of thedistribution. This continuum disappears completely for the RRPA (figure 12.d).1 C. Partial pair distribution functions
The structure of random sphere packings is usually characterized by the probability per unitvolume of finding a sphere center at a distance r from another sphere center, P ( r ) × N/V , where N is the number of spheres in the aggregate of volume V and the pair distribution function (PDF) P ( r ) is normalized to 1 when r → ∞ .Distinguishing spheres by their contact coordination number allows a much more detailed struc-tural description by partial pair distribution functions (PPDF), which are defined as the probability P ij ( r ) of finding a sphere with contact coordination number j at a distance r from another spherewith contact coordination number i , normalized to 1 at large r .
1. Principle for partial PDF P ii ( r ) , P ij ( r ) and P i ( r ) In practice, coordination numbers range from 3 to 12, the maximum CCN in 3D space. For aspherical aggregate with radius R , P ij ( r ) write: P ij ( r ) = V N i N j ∆ N ij (2 − δ ij ) S ( r )∆ r (7)where N i and N j are the number of spheres with CCN i and j , respectively, in the volume V =4 πR / N ij is the number of sphere pairs of CCN i and j respectively, lyingin the interval [ r : r + ∆ r ]; δ ij is the usual Kronecker symbol; S ( r ) is the spherical shape factor ofthe aggregate [22]: S ( r ) = π r (2 R − r ) (4 R + r ) (8) P ii ( r ) PPDF describe the arrangement of i coordinated spheres, while P ij ( r ) PPDF with i (cid:54) = j describe the mutual arrangement or ”chemical order” between i and j coordinated spheres.In the case of sticky hard spheres with diameter d , the peak of contacting neighbours in P ij ( r )is represented by [23]: P ij ( d ) = ¯ η ij V πN j d δ ( r − d ) i . e . numerically η ij = 3 N j R d σP ( d ) (9)where ¯ η ij is the average number of j coordinated spheres contacting an i coordinated sphere and σ = 0 .
01, is the length step used in P ( r ) calculations. Numerical values of P ij ( d ) fall out of rangeof the P ij ( r ) figures presented below.2One can then define the probability P i ( r ), to find a sphere with any coordination number at adistance r from a sphere with coordination number i . It writes: P i ( r ) = (cid:88) j =3 C j P ij ( r ) (10)where C j = N j /N is the concentration of j coordinated spheres ( N = (cid:80) j =3 N j ).These P i ( r ) characterize the global arrangement of spheres around a sphere with coordination i . Finally, the global PDF writes: P ( r ) = (cid:88) i =3 C i P i ( r ) (11) a. Random regular polytetrahedral aggregates Figure 13 presents P ii ( r ) curves (with i = 4 , γ = 0 . FIG. 13. Global and partial pair distribution functions ( P ( r ) and P ii ( r ) for i = 4 , Concerning the global PDF, a striking difference from what was observed for RIPA (see [9]),is the disappearance of the topological discontinuities at r = √ d and r = 2 d . The continuousstructure of the PDF observed in RIPAs almost disappears and is replaced by a set of polyhedral δ peaks ([13, 23]), which are due to a large (virtually infinite) regular polytetrahedron. The positionsof these peaks, which are due to precise configurations of spheres with well defined distances, areidentical to those studied by Medvedev et al [13]. The continuous structure grows almost linearlyas a function of r and goes to 1 for r ≈ . d = 7.3The partial PDFs ( P ii ( r )) show a systematic tendency, as their continuous regime in the regionof the quasi first neighbours (QFN) starts from 0 for spheres with coordination number i equalor superior to 8 and increases like the global P ( r ). The lower coordination numbers ( i <
8) leadto P ii ( r ) starting from higher values in the QFN area. This can be qualitatively understood: themore contacting first neighbours a sphere has, the less quasi first neighbours it can accept. b. Random irregular polytetrahedral aggregates The global PDF of RIPA was studied in [9],with the exception of RMIN aggregates. However, the latter bring no qualitative differences tothese results. P ii ( r ) The PPDF P ii ( r ) together with the corresponding PDF are presented in figures 14.aand b for MAX-1 aggregates with the two most extreme packing fraction, γ = 0 .
586 and γ = 0 . r close to d ) decreases when the coordination number ( i ) increases, like in the case ofRRPA and for the same reasons: the most contacting neighbours a sphere has, the less quasi firstneighbours it can accept. (a) (b) FIG. 14. Sample P ( r ) and P ii ( r ) obtained for aggregates generated by algorithms a) MAX-1 ( γ = 0 . γ = 0 . On the other hand, the comparison between figure 14.a and b, shows that the number of quasifirst neighbours decreases with packing fraction. At low packing fractions, they form plateauswhose level increases when the coordination number decreases, whereas they have a very distinctbehaviour for high packing fractions: low coordination number spheres possess a lot of quasi firstneighbours while high coordination number ones have a very limited number of QFN.4Both topological discontinuities (at r = √ d and r = 2 d ) increase when the coordination numberincreases and when packing fraction increases.Finally, (see supplementary informations for figures) ”iso-packing-fraction” and ”iso-¯ κ D ” ag-gregates obtained by different building algorithms show significant differences in their global andpartial PDF, confirming that parameter γ and ¯ κ D are greatly insufficient for a full structural de-scription of random aggregates. Aggregates sharing similar values of these two parameters ( γ and¯ κ D ) still present significant differences in structural properties. Hence, even when used together, γ and ¯ κ D are not satisfying predictors of structural properties of disordered systems. Besides, it hasbeen impossible to find aggregates with similar ¯ L max and γ to compare them in a similar fashionleaving the question open for a possible combination of these two parameters as good predictor ofrandom aggregates global structure. Conversely, it turns out that PPDFs allow a more sensitivedistinction between aggregates built by different algorithms than the corresponding global PDFand are thus interesting structural descriptors. P ij ( r ) with i (cid:54) = j A sampling of P ij ( r ) curves is shown in supplementary information. Theyshow that QFN are favoured by higher packing fractions and lower coordination numbers i and j .The implicit ”chemical” ordering in these curves will be studied in more details hereafter. P i ( r ) P i ( r ) (for i = 4 , i and packing fraction. On the one hand, for each packing fraction: • quasi first neighbours increase when the coordination number decreases; • topological peaks at r = d √ d are more intense for high coordination numbers;On the other hand, the main differences between low and high packing fraction aggregates arethat: • polytetrahedral δ peaks are noticeable on aggregates with low packing fractions; • the transition from the high to the low packing fraction is mostly due to a decrease of quasifirst neighbours around low coordination spheres; • in the case of the densest aggregate, P i oscillations are increasingly shifted toward smaller r values as the coordination number i increases, while such a shift is not observed for thelowest packing fraction aggregate.5 (a) (b) FIG. 15. P i ( r ) for two packing fractions, algorithm MAX-1. D. Radial evolution of local packing fraction from pair distribution functions
1. Principle
Knowing pair distribution functions, it becomes possible to study the variation of the localpacking fraction around an average sphere as a function of r . Using relation 3, the packing fractionof a sphere S with arbitrary diameter R SC > d situated within the aggregate writes: γ ( R SC ) = 4 πr s / ρ (cid:82) R SC + r s P ( r )4 πr V s ( r, R SC ) dr πR / ρ = N/V , the number of spheres per unit volume in the aggregate; the 4 πr s / r = 0; V s ( r ) is defined by equation 3; P ( r ) can be a total or partial pairdistribution function. It also comes γ ( r ≤ r s ) = 1.It is also possible to remove the contribution of the central sphere and its contact first neighboursto packing fraction, which writes: γ CN ( R SC ) = 4 πr s / iV s ( r = d, R SC )4 πR / i is the number of contacting neighbours of the central sphere. It then comes: γ W CN ( R SC ) = γ ( R SC ) − γ CN ( R SC ) (14) γ W CN is the contribution to packing fraction of quasi-first and further neighbours around anaverage sphere.6
2. Local packing fraction around average sphere
Using the global PDF, the formalism just introduced gives access to the variation of the localpacking fraction as a function of r around an average sphere. Figure 16.a represents γ ( R SC )and figure 16.b represents γ ( R SC ) /γ , which goes to 1 as R SC goes to ∞ . All aggregates behavedifferently depending on their packing fraction. For the densest ones, the packing fraction fallsbelow the average value, before converging more rapidly at large r than the least dense ones. Onthe other hand, in the case of the least dense ones, packing fraction remains above the average valueand converges more slowly. Oscillations of local packing fractions are damped at about r = 5 d forthe lowest packing fraction and around r = 3 d for the highest packing fraction.
3. Local packing fraction around i coordinated spheres Using P i ( r ), it is possible to determine the local packing fraction as a function of r for anaverage sphere with coordination number i , γ i . Figures 16.c and d present the sets of correspondingcurves for the two MAX-1 aggregates with the lowest and highest packing fractions. The short-range local-density grows with the coordination number and the packing fraction of the aggregate.Significant oscillations are observed for the aggregate with the highest packing fraction. Theiramplitude is damped for the aggregate with the lowest packing fraction. Figures 16.e and f showthe evolution of packing fraction without the contribution of contact neighbours ( γ W CN ) and thecentral sphere, hence only that of quasi-first and further neighbours. As it could be expected,the lower number of quasi first neighbours (QFN) around spheres with high coordination numbersleads to a significantly lower packing fraction due to QFN than in the case of spheres with lowcontact coordination number.In the case of the aggregate with the highest packing fraction (fig. 16.e), the various curvesconverge for R SC ≈ . d and the highest difference is observed between i = 10 and i = 3 for R SC ≈ . γ W CN, − γ W CN, ≈ .
28 = 0 . γ ag . Between i = 10 and i = 6, the maximum isfound at R SC ≈ .
57, with γ W CN, − γ W CN, ≈ .
16 = 0 . γ ag .There is a much lower difference in contribution of the QFN between spheres with differentcontact coordination number in the case of the aggregate with the lowest packing fraction (interms of absolute as well as relative value). However, the maximum of the difference between γ W CN, and γ W CN, is obtained for approximately the same value R SC ≈ .
55 and γ W CN, − γ W CN, ≈ .
06 = 0 . γ ag . Between i = 10 and i = 6, the maximum is found at R SC ≈ .
51, with7 (a) (b)(c) (d)(e) (f)
FIG. 16. a) Local packing fraction b) local packing fraction divided by average PF for various MAX-1aggregates. Radial evolution of the packing fraction around i coordinated spheres for MAX-1 aggregateswith c) γ = 0 .
586 and d) γ = 0 . i coordinated sphereswithout contact neighbours (relation 14) for MAX-1 aggregates with e) γ = 0 .
586 (the bottom lines representthe difference between the highest and lowest bound, i.e. γ W CN, ( R SC ) − γ W CN, ( R SC ) and intermediate,i.e. γ W CN, ( R SC ) − γ W CN, ( R SC )) and f) γ = 0 . γ W CN, − γ W CN, ≈ .
06 = 0 . γ ag , which is remarkably similar to the former. In this case, figure16.e shows that γ W CN ( R SC ) for coordination numbers below 7 behave very similarly, suggestingthat these types of spheres have more or less the same environment in terms of number of QFN,irrespective of their contact coordination number. E. Global and partial structure factors
The previous structural description achieved by the PPDF analysis of random packings leadsnaturally to the study of the corresponding global and partial structure factors which can becompared with the results of diffraction experiments on disordered materials.
1. Global and Ashcroft Langreth partial structure factors
Let us first introduce the global structure factor: S ( Q ) = 1 + NV (cid:90) ∞ (cid:16) P ( r ) − (cid:17) sin( Qr ) Qr πr dr (15)and the partial Ashcroft Langreth (AL) structure factors, defined by: S ij ( Q ) = δ ij + (cid:112) N i N j V (cid:90) ∞ (cid:16) P ij ( r ) − (cid:17) sin( Qr ) Qr πr dr (16)where δ ij is the Kronecker symbol.The diagonal terms, S ii , represent the structure factors of the partial aggregates formed byspheres with contact coordination i and will be studied first. The non diagonal terms S ij with i (cid:54) = j describe the mutual or ”chemical” arrangement between i and j coordinated spheres and willbe studied later by the more suited Bhatia-Thornton formalism. S ( Q ) and S ij ( Q ) values for Q (cid:28) /R cannot be calculated from relations 15 and 16 owing tothe finite radius R of the aggregates [9]. However S ( Q = 0) values are determined by the statisticaldensity fluctuations according to relation [24]: S ij (0) = N i N j − N i N j (cid:113) N i N j (17)where upper bars represent averages.Fluctuations were directly derived from the positions of the sphere centers. To that purpose, alarge cube with edge 50 d centered on the aggregates origin has been subdivided into 1000 sub-cubeswith edge 5 d , each of them containing N i and N j sphere centers with respective coordination i and9 j . Average values of N i , N j and N i N j were taken over the 1000 sub-cubes to get the statisticalfluctuations involved in relation 17. a. Random regular polytetrahedral packing – global structure factor Figure 17 introduces theglobal structure factor of two RRPA: 17.a (small Q ) and b (large Q or asymptotic regime). Itturns out that the structure factors of RRPA exhibit a prepeak, whose position depends on theused algorithm or packing fraction. The position and intensity of this pre-peak corresponds tothe folding of the polytetrahedra which controls the space correlation of density fluctuations. Itaccounts for interferences between high density and low density zones of the single polytetrahedron.At high Q , S ( Q ) appears aperiodic due to the infinite series of δ ( r p ) peaks in P ( r ) correspondingto the distances between vertices in the infinite RRPA and is only weakly affected by small changesin packing fraction (depending on the used algorithm). (a) (b) FIG. 17. a) Structure factor of two RRPA aggregates and experimental results obtained on amorphous Ge(data obtained in [25]) at small Q b) large Q behaviour of the structure factor of the RRPA aggregates andone (high packing fraction) RA aggregate. Comparison with the structure of tetravalent elements
The polytetrahedral struc-ture factor can be usefully compared with the experimentally measured structure factors of puretetravalent elements (C [26], Ge [25], Si [27]) or with the tetrahedral structure of oxygen atoms inglassy water [28, 29]. As a matter of fact the structure of these elements is based on body centeredregular tetrahedra (BCT) connected by their vertices [30], while RRPA’s are made from simpleregular tetrahedra connected by their faces. As a consequence, a ”pseudo superstructure” peak forsuch hard-sphere based centered-tetrahedra should correspond to Q s /Q = d /d s = √ / . d is the distance between contacting atoms and d s is the edge length of the regular centeredtetrahedron. Experimentally, Etherington finds Q s /Q = 0 .
572 [25] for the structure factor ofamorphous Ge, Laaziri et al. find 0 .
583 for the same ratio in a-Si [27] and Gaskell et al. find0 .
540 for a-C [26]. The discrepancy from 0.612 for these three values could be due to the errorintroduced by hard sphere approximation of tetrahedral bond angle between first neighbours ofsofter potentials, which does not allow reproducing perfectly the behaviour of covalent materials.However, this superstructure peak is lacking from the structure factor of RRPA’s but anothersmall angle peak at significantly lower value is found around Q s /Q ≈ .
25, varying with packingfraction. The latter can be attributed to the pseudo periodic correlation of density fluctuationspresent in the RRPA model. These correlations happen over larger distances than those of theBCT, hence the smaller Q values of the corresponding pre-peak. It corresponds roughly to the endof significant oscillations of pair distribution functions.The reason why this longer correlation-distances pre-peak does not appear for tetravalent el-ements has probably to do with the fact that tetrahedra forming their structures are connectedby vertices, allowing for many more degrees of freedom and preventing longer distances correla-tions than in the case of RRPA aggregates where tetrahedra are connected by their faces, whichis geometrically more constraining. On the other hand, the second peak of the structure factorof the RRPA splits into two subpeaks in a completely analogous manner with what is found inamorphous tetravalent elements.Finally, for RRPA aggregates, characteristic distances are strictly defined, giving rise to δ peakson the PDF resulting in aperiodic large Q behaviour, also absent from the structure factor oftetravalent elements. In contrast, the lack of such well defined distances, because of higher degreesof freedom in the orientations of tetrahedra relatively to one another in the case of tetravalentelement, results in a single δ peak on their PDF, which accounts for their periodic behaviour inthe large Q regime [23]. b. Random irregular polytetrahedral aggregates The AL partial structure factors of random ir-regular polytetrahedral aggregates are presented in figure 18.a to c, along with the partial structurefactors of the RRPA (18.d). These structure factors all have in common that the global amplitudeof their oscillations increases from i = 3 up to i = 6 and then decreases as i increases beyond6. Moreover, the position of the first peak of S ii goes to small Q when packing fraction increases(whatever i and γ ). At large Q , there is no phase shift, whatever i and γ , and S ii oscillates likesin( Qd ) /Qd (not shown here).Some peculiarities are also observed. A shoulder appears on the left of the first peak of S Q = 3 . S and S an intense pre-peak around Q = 0 . S and decreases when the packing fraction increases for S . This pre-peak is due to correlated density fluctuations between low or high density regionsseparated by an average distance of 2 π/ . r s ≈ r s lying beyond the main P ii oscillations; theyare observed whatever the packing fraction, even if there is no noticeable pre-peak on the globalstructure factor. In the case of S , S (0) decreases when γ increases and barely any pre-peakmay be observed at all. The first peak remains symmetrical; its position depends little on packingfraction. The first peak of S becomes strongly asymmetrical for the highest packing fraction.Finally, a splitting of the second peak is observed on the AL structure factors with high coordi-nation numbers and low packing fraction. This splitting has been noticed above for the RRPA andit is consistent with the existence of a RP structure embedded in a ”continuum” random struc-ture (respectively RP and FR structural components), as was already concluded in [9] for theseaggregates. (a) (b) FIG. 18. Ashcroft-Langreth partial structure factor for a) highest packing fraction RIPA, b) RRPA.
2. Bhatia Thornton partial structure factor
The partial structure factors introduced by Bhatia and Thornton [31] are associated with densityand concentration correlations in binary alloys. They have been extended to alloys consisting ofmore than two components [32] and can then be used by considering random aggregates as alloysof spheres with various coordination numbers. a. Formalism
The partial structure factor corresponding to:2 • density-density correlation function, S NN , writes: S NN ( Q ) = 1 + ρ (cid:90) ∞ (cid:16) P ( r ) − (cid:17) sin( Qr ) Qr πr dr (18) • concentration-concentration correlation function, S C i C j , writes: S C i C j ( Q ) = 1 + ρ (cid:90) ∞ (cid:16) − P ( r ) + P i ( r ) + P j ( r ) − P ij ( r ) (cid:17) sin( Qr ) Qr πr dr (19) • density-concentration, S NC i , writes [12, 23]: S NC i ( Q ) = ρ (cid:90) ∞ (cid:16) P i ( r ) − P ( r ) (cid:17) sin( Qr ) Qr πr dr (20)For the sake of simplicity, we only consider hereafter sphere mixtures made of two components,namely i and j coordinated spheres. The corresponding relations are provided in supplementaryinformations. b. Results Between 15 (for high packing fraction) and 45 sets of Bhatia-Thornton structurefactors corresponding to different coordination pairs i − j were calculated for 500 aggregates. Inthe case of high packing fraction aggregates where ”extreme” coordination numbers (i.e. i = 3 and i >
9) are scarce, BT sets could only be calculated for coordination numbers lying between 4 and8, hence the lower number of BT structure factors calculated for them.A sample of Bhatia Thornton partial structure factors is provided in supplementary informationsfor three MAX-1 aggregates, with maximum, minimum and intermediate packing fraction alongwith those of one RRPA aggregate.Their main characteristics are the following: S N i N j The S N i N j structure factor is related to the overall structure of the two coordinations i and j considered. It looks like the average AL structure factor for small coordination differences( i − j = 1 or 2 at most), while it resembles more to the global structure factor of the aggregate for”large” differences in coordination numbers ( i − j ≥ S NC i structure factors exhibit significant oscillations whose intensity increases with the coor-dination difference i − j . Furthermore, for low coordination differences these oscillations decreasewhen the coordination numbers increase. All these oscillations come from the fact that the av-erage environment of sphere i and j characterized by P i and P j differ from the average globalenvironments P ( r ) (see relation 20). This behaviour differs completely from the case of binary”substitution” alloys (with equal atomic diameters) for which S NC does not oscillate at all because3 P i ( r ) ≈ P ( r ), independently of the chemical order of the two alloy components (which only affects S CC [23]). S C i C j Most interesting is the study of S C i C j , which characterizes the ”chemical” order between i and j coordinated spheres through its dependence on P ii − P ij (see equation 19). For smallcoordination differences (at most 2) and all packing fractions, S C i C j oscillates weakly around 1 andthere is no ”chemical” order effect between i and j coordinated spheres.For large coordination differences ( i − j ≥
4) a ”pre-deep” is observed on S CC around Q ≈ . Q .This pre-deep indicates a segregation effect [33], between spheres with large CCN differences. Itvaries non uniformly with packing fraction (see figure 19.a) and its maximum amplitude lies inthe packing fraction interval γ ∈ [0 .
5; 0 . (cid:104) CCN (cid:105) seen in sectionV A 3. As a matter of fact, the maximum segregation measured by this pre-deep occurs in a rangeof packing fraction where there is virtually no contact segregation in the case of RMIN-MAX-1aggregates (by comparing fig 9.b and fig 9.a), suggesting that there can exist contact and longer-range segregation effects between spheres of various coordination numbers and that these two typesof segregation are not necessarily concomitant. (a) (b)
FIG. 19. a) Amplitude of S CC − minimum ( Q ∈ [1 . − . Q maxima inBhatia-Thornton S CC − structure factors ( Q <
Furthermore, for these high coordination differences a new S CC prepeak appears at very small Q values, around 0 . Q , which seems to vanish in the aggregates with the highest packing fraction(fig 19.b), with the notable exception of MAX-4. Accordingly, this S CC prepeak is also observed4in the RRPA. This peak should be associated with the existence of segregation effects extendingto larger r range, i.e. to the formation of isocoordinated aggregates in an ”average matrix”. Thissuggests the existence of regions of low or high coordination in the aggregate matrix. VI. DISCUSSION AND CONCLUSION
The so called random packings of adhesive hard spheres cover a wide range of packing fractions,from a lower limit of 0.15 at the three dimensional percolation threshold [34] with contact coordi-nation number 2 [23], up to a maximum value of about 0 .
636 in the close packed random packing(RCP) with coordination number approximately 8.1, as determined experimentally by Scott [35].In this paper we focused on the detailed structural characterization of such packings, using awide variety of static aggregates-building algorithms.The randomness or isotropy of the 10 aggregates (containing 10 spheres each one) was checkedfirst by studying the angular distribution of pairs of spheres separated by a given distance. Theeigenvalues of the corresponding nematic order tensor were shown to be always less than 2 . − and confirmed that aggregates are fully isotropic.Accurate methods for the determination of the aggregate packing fraction were then introduced.It was thus shown that seed effects at the origin of the aggregates do not extend beyond thefifth neighbour range, i.e. are limited to the first 8 000 spheres, if the seed consists of spherearrangements that are rare in the aggregate. This seed effect can be totally removed by using asseed a sample with a structure similar to the built aggregate.On the other hand it has been shown that the effect of the imperfect spherical surface of theaggregates depends on the building algorithm and packing fraction, and is limited to a thicknessof approximately 6 r s in the ”worst case”. Finally, the accuracy on the (average) packing fractionincreases with the sphere number of the aggregate which must reach 10 in order to get a 10 − accuracy.The detailed structural analysis of random packings could then be undertaken.As a first step, these structures have been tackled via the Delaunay tessellation of tetrahedraconnecting sphere centers. The distributions of these tetrahedra have been characterized by twodistortion parameters, L max and κ D . Their distributions show a bimodal character that varieswith packing fraction and appear related, they both include discontinuities whose origin remainsunclear. Their average values decrease with increasing packing fraction. However, they providecomplementary characterization of the aggregates produced by the various studied algorithms,5suggesting that beyond their similarities, these distributions also present subtle differences thatwill be the object of a future study.Special attention was paid to the populations of regular tetrahedra (formed by 4 mutuallycontacting spheres, L max = 2) and quasi regular Delaunay tetrahedra ( L max < .
3) which wereshown to behave quite differently. As a matter of fact, the volume fraction of regular tetrahedrahas been shown to decrease with increasing packing fraction and reaches a maximum value of 0.165for RPPA aggregates (only built with regular construction tetrahedra). This raises a fundamentalquestion: is there a maximum geometrically defined (i.e. irrespective of their building mechanism)proportion of regular tetrahedra in random aggregates? Conversely, the proportion of quasi regularDelaunay tetrahedra goes to a minimum around γ = 0 .
56 and then increases with increasing packingfraction.New structural characterization methods could then be introduced by taking advantage of theunequivocal definition of sphere contacts and, hence, contact coordination numbers.First, partial characterization was carried for short distances, i.e. contacting neighbours. Thedistributions of pairs of spheres with respective CCN i and j ( η ij ) were first studied. Their FWHMincreases when packing fraction decreases and their average values ( (cid:104) η ij (cid:105) ) show distinct behaviourswith respect to γ for low and high values of i respectively: a non uniform decrease with γ isobtained only in the case of low i values. The evolution of these distributions shows that contactsegregation can exist between spheres of various CCN, whereas the evolution of (cid:104) CCN (cid:105) for spheresof various coordination shows that another segregation may exist, on the basis of CCN, over alarger range. Partial distribution of bond angles for spheres with different CCN were also studied.It was shown that high CCN spheres tend to have a smoother distribution, whereas low CCN onespresent a depletion of low angle bonds and have a higher proportion of high angle bonds. Thesedifferences are particularly distinct in high packing fraction aggregates.This structural description was then extended to all distances by introducing the partial pairdistribution functions P ij ( r ) i.e. the probability of finding a sphere with CCN i at a distance r from another sphere with CCN j (normalized to 1 for large r) and local packing fractions around i coordinated spheres. The different shapes and discontinuities of the P ij ( r ) curves were given andthis detailed analysis allowed a clear cut distinction between different random packing structureswhich cannot be distinguished by their packing fraction and/or κ D parameters. The questionremains of whether these aggregates could be characterized by characterizing them through thecombination of packing fraction γ and distortion parameter L max .Distinguishing (with Bernal [36, 37]) contacting from quasi-contacting spheres, it was shown6that high CCN spheres have few quasi contacting neighbours (QCN) while low CCN spheres have ahigher number of QCN (as could be expected) and that this local effect extends to larger distances.Moreover this behaviour is amplified when the packing fraction increases.Finally, the two sets of partial structure factors respectively introduced by Ashcroft-Langreth( S ij ( Q )) and Bhatia-Thornton ( S NN , S NC and S CC ) were analyzed. They are different Fouriertransforms of linear combinations of the P ij ( r ) and can be directly compared with the resultsof diffraction experiments on liquid or disordered materials. In particular the diagonal S ii ( Q )give the partial structure factors of the partial aggregates made from i coordinated spheres andthe variable shapes of these S ii , especially their prepeak, first and second Q oscillations, weredescribed. Furthermore the concentration-concentration partial structure factor of Bhatia gavethe mutual or chemical order between spheres with different coordinations i and j and put forwardhetero-attractions or segregation of spheres within the different aggregates, confirming that severalkind of segregation may indeed exist.From all these results we could conclude that the random irregular polytetrahedral aggregatesstudied here are made from two basic components, namely fully random component (without regu-lar tetrahedra) and regular polytetrahedral component (only built with regular tetrahedra), whoseproportion decreases with increasing packing fraction. This composite nature of the aggregatesproduces prepeaks (at very small Q ) in the aggregate structure factors and the main featuresdifferentiating the RP component from the counterpart FR component are the following: • A high proportion of regular Delaunay tetrahedra (as expected) • an increase in δ -peaks noticeable in: distributions of Delaunay tetrahedra distortion indexes,angular bond distributions, global or partial pair distribution functions • an increase of the second mode in the distributions of Delaunay tetrahedra distortion pa-rameters κ D and L max • a spreading of contact coordination numbers distributions η ij • strong variations in the radial dependency of the average contact coordination number ofspheres • a collapse of the continuum of partial bond angle distributions • a more similar quasi-first-neighbours part of partial pair distribution functions of low andhigh contact coordination number spheres7 • a more similar radial dependency of packing fraction for spheres of any contact coordinationnumber • a reduction of topological discontinuities in pair distribution function at r = 2 √ • a splitting of the structure factor’s second peak ( Q ≈ [1] B. Alder, T. Wainwright, Studies in molecular dynamics, Journal of chemical physics 31 (1959) 459.[2] B. D. Lubachevsky, F. H. Stillinger, Geometric properties of random disk packings, J. Stat. Phys. 60(1990) 561–583.[3] W. S. Jodrey, E. M. Tory, Computer simulation of close random packing of equal spheres, Physicalreview A 32 (1985) 2347–2351.[4] N. Karayiannis, M. Laso, Monte carlo scheme for generation and relaxation of dense and nearly jammedrandom structures of freely jointed hard-sphere chains., Macromolecules 41 (2008) 1537–1551.[5] G. Y. Onoda, E. G. Liniger, Random loose packings of uniform spheres and the dilatancy onset, Phys.Rev. Lett. 64 (1990) 2727–2730. [6] K. J. Dong, R. Y. Yang, R. P. Zou, A. B. Yu, Role of interparticle forces in the formation of randomloose packing, Phys. Rev. Lett. 96 (2006) 145505.[7] A. Anikeenko, N. Medvedev, Polytetrahedral nature of the dense disordered packings of hard spheres,Physical review letters 98 (2007) 235504.[8] R. Jullien, P. A., P. Meakin, Random packings of spheres built with sequential models, J. Phys A:Math. Gen. 25 (1992) 4103–4113.[9] M. Bletry, J. Bletry, Fluctuations, structure factor and polytetrahedra in random packings of stickyhard spheres, Journal of Non-Crystalline Solids 411 (2015) 85–100.[10] Hereafter, unless otherwise mentioned, all distances are expressed in r s as the unit of length ( r s = 1).Similarly, r − s is the unit in reciprocal space.[11] Yet another family of algorithms has been developed, for which the newly added sphere is positionedby controlling the value of the building tetrahedra distortion. However, this strategy leads to the exactsame P ( r ) as the one obtained for RA built with the positioning in function of the hole size and hasthe same maximum packing fraction: these aggregates will not be further studied here.[12] J. Bl´etry, Diffraction des neutrons polaris´es par les phosphures de cobalt amorphes. interpr´etation pardes mod`eles structuraux applicables au m´etaux et aux liquides., Ph.D. thesis, INPG (1979).[13] N. Medvedev, E. Pilyugina, Three-dimensional packing of perfect tetrahedra, in: K. Sugihara, D. Kim(Eds.), Voronos Impact on Modern Science, Vol. 4 of Proc. of the 5th International Symposium onVoronoi Diagrams in Science and Engineering, 2008, pp. 144–156.[14] E. W. Weisstein, Sphere-sphere intersection, A Wolfram Web Resource Edition, MathWorld.URL http://mathworld.wolfram.com/Sphere-SphereIntersection.html [15] L. E. Silbert, Jamming of frictional spheres and random loose packing, Soft Matter 6 (2101) 2918–2924.[16] B. Delaunay, in: Proceeding of the mathematical congress, 1924, pp. 695–700.[17] CGAL, Computational Geometry Algorithms Library.URL [18] C. Jamin, S. Pion, M. Teillaud, 3D triangulations, in: CGAL User and Reference Manual, 4.9 Edition,CGAL Editorial Board, 2016.URL http://doc.cgal.org/4.9/Manual/packages.html [19] A. Anikeenko, N. Medvedev, T. Aste, Structural and entropic insights into the nature of the rendom-close-packing limit, Physical review E 77 (2008) 031101.[20] M. Skoge, A. Donev, F. Stillinger, S. Torquato, Physical Review E 74 (2006) 041127.[21] N. Karayiannis, K. Foteinopoulou, M. Laso, The structure of random packings of freely jointed chainsof tangent hard spheres, The Journal of Chemical Physics 130 (16) (2009) 164908.[22] G. Fournet, Bul. Soc. Fr. Min. 74 (1951) 37.[23] J. Bl´etry, Sphere and distance models for binary disordered systems, Philosophical Magazine B 62(1990) 469.[24] J. G. Kirkwood, F. P. Buff, The statistical mechanical theory of solutions. i, The Journal of Chemical Physics 19 (6) (1951) 774–777.[25] G. Etherington, A. Wright, J. Wenzel, J. Dore, J. Clarke, S. R.N., Journal of Non-Crystalline Solids48 (1982) 265–289.[26] P. Gaskell, A. Saeed, P. Chieux, D. McKenzie, Physical Review Letters 67 (10) (1991) 1286–1289.[27] K. Laaziri, S. Kycia, S. Roorda, M. Chicoine, J. Robertson, J. Wang, S. Moss, Physical Review Letters82 (17) (1999) 3460–3463.[28] A. Narten, M. Danford, H. Levy, Discussions of the Faraday Society 43 (1967) 97.[29] A. H. Narten, C. G. Venkatesh, S. A. Rice, The Journal of Chemical Physics 64 (3) (1976) 1106–1121.[30] D. Polk, Journal of Non-Crystalline Solids 5 (1974) 365.[31] A. Bhatia, D. Thornton, Physical Review B 2 (1970) 3004.[32] J. Bl´etry, Z. Naturforsch. 31 a (1976) 960–966.[33] J. Bl´etry, Zeitschrift Naturforschung Teil A 33 (1978) 327.[34] R. Zallen, The physics of amorphous solids, Wiley, New York, 1983.[35] G. Scott, Packing of spheres: Packing of equal spheres, Nature 188 (1960) 908.[36] J. Bernal, J. Mason, Co-ordination of randomly packed spheres, Nature 188 (1960) 910.[37] J. Bernal, The structure of liquids, Proceedings of the Royal Society of London 280 (1962) 299.[38] R. D. Mountain, T. Ruijgrok, Monte-carlo study of the maier-saupe model on square and trianglelattices, Physica A: Statistical Mechanics and its Applications 89 (3) (1977) 522 – 538.
ACKNOWLEDGMENTS
P. C´en´ed`ese is greatly acknowledged for his help in openMP parallelization of some codes.Anonymous referees are thanked for their observations and suggestions that helped in the processof writing the present paper.Anonymous referees are greatly acknowledged for their critics and suggestions that helped im-proving the quality of this paper.
Appendix A: Supplementary information1. Isotropy of the aggregates
Isotropic disordered aggregates are characterized by an isotropic distribution of ˆ r ij bonds, whereˆ r ij are the unit vectors joining the (i,j) sphere centers. This requirement provides a convenientcheck for the building method and/or the minimum aggregate size beyond which the isotropy isreached. Such a distribution is characterized by a uniform random distribution of azimuthal angles ϕ and a polar angles θ distribution following a probability density P b ( θ ) = sin( θ ) / θ –distribution presenting distinguishable peaks; thus the features of the θ distribution appears tobe a convenient tool to control both the aggregates randomness and isotropy.In order to go beyond a simple qualitative characterization, we can quantify the deviation fromisotropy by looking if a favoured direction emerges among the { ˆ r ij } distribution. For this, we followthe usual method of the framework of the nematic liquid statistics [38] according to which, fromthe diagonalization of the two rank tensor¯ Q = 1 N (cid:88) ij
12 (3ˆ r ij ˆ r ij − ¯ I ) (A1)where ¯ I is the identity tensor, a nematic order parameter, say S is obtained as the largest eigen-value λ max of ¯ Q . In equation (A1), the sum running over all the { i, j } pairs, can be limited to thesum over the bonds ˆ r ij , which can be replaced by a sum over the unit vectors carried by the (cid:126)r i ,with respect to a fixed reference (cid:126)r o , which avoids the surface effects when both i and j are locatedat the aggregate surface.The results for the polar angle distributions relative to the axis ˆ x , ˆ y and ˆ z for the most andleast dense aggregates obtained by MAX-1 algorithm are presented respectively in figures 20.a and20.b. For the aggregates studied in this work, the deviation from isotropy as measured by the valueof λ max is found to decrease when the packing fraction increases. A reliable determination of thedependence of λ max with respect to γ is beyond the scope of this paper, all the more that the rangeof λ max values is quite small. The largest eigenvalue of ¯ Q is smaller than 3 10 − . Such a value isvery small leading us to conclude that these building methods lead indeed to isotropic aggregates.To illustrate quantitatively this point, figure 20.c displays the result of the θ distribution in termsof λ for the following model: P b ( θ ) = sin( θ )2 (cid:104) e ( − θ / (2 σ )) + e ( − ( π − θ ) / (2 σ )) (cid:105) (A2)where the polar axis is a favoured direction according to the variance σ . Clearly, the θ distributionscharacterized by λ values lower than a few 10 − can be considered isotropic.1 (a) (b)(c) FIG. 20. Pair angle distributions with respect to x, y and z for the (a) most and (b) least dense aggregatesproduced by MAX-1 algorithm. c) Polar angle distribution relative to the favorable direction in the uniaxialanisotropy case for different values of the nematic order parameter, λ defined as the largest eigenvalue ofthe tensor ¯ Q , equ.(A1) P ii ( r ) (a) (b)(c) (d)(e) (f) FIG. 21. Sample P ( r ) and P ii ( r ) obtained for aggregates generated by algorithms a) MAX-1 ( γ = 0 . γ = 0 . κ D = 1 .
771 d) RAN1 ¯ κ D = 1 .
771 e) RAN4 f) RMIN-MAX1 . P ij ( r ) (a) (b)(c) (d) FIG. 22. A sample of P ij ( r ) for various packing fractions, algorithm MAX-1 – a) P , ( r ), b) P , ( r ), c) P , ( r ), d) P , ( r )
4. AL structure factors (a) (b)(c)
FIG. 23. Sample Ashcroft Langreth structure factors for MAX-1 aggregates. a) S , b) S , c) S .
5. Bhatia-Thornton formalism for binary mixtures of spheres