Structure of sticky-hard-spheres random-aggregates: the viewpoint of contact coordination and tetrahedra
M. Blétry V. Russier E. Barbé J. Blétry





[Corresponding author ]bletry@icmpe.cnrs.fr



ICMPE-CNRS, 2-8 rue Henri Dunant, 94320 Thiais, France

Honorary Professor, Bahia Blanca University, Argentina


We study more than 10^4 random aggregates of 10^6 monodisperse sticky hard spheres each, generated by various static algorithms. Their packing fraction varies from 0.370 up to 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 on pair distribution functions. 
Distortion of structural (Delaunay) tetrahedra is studied by two parameters, which show some similarities and some differences in terms of overall tendencies. Isotropy of aggregates is characterized by the nematic order parameter.
The overall structure is then studied by distinguishing spheres in function of their contact coordination number (CCN). 
Distributions of average CCN around spheres of a given CCN value show trends that depend on packing fraction and building algorithms. 
The radial dependency of the average CCN turns out to be dependent upon the CCN of the central sphere and shows discontinuities that resemble those of the pair distribution function.
Moreover, it is shown that structural details appear when the CCN is used as pseudo chemical parameter, such as various angular distribution of bond angles, partial pair distribution functions, Ashcroft-Langreth and Bhatia-Thornton partial structure factors. These allow distinguishing aggregates with the same values of packing fraction or average tetrahedral distortion or even similar global pair distribution function, indicative of the great interest of paying attention to contact coordination numbers to study more precisely the structure of random aggregates.



Structure of sticky-hard-spheres random-aggregates: the viewpoint of contact coordination and tetrahedra
    Blétry, J
    
========================================================================================================






§ INTRODUCTION

Random aggregates of monodisperse spheres are of great interest to simulate various physical systems, such as amorphous solids, liquids, powders, etc. In particular, a fine characterization of their local or longer range structural properties is of interest to classify the various types of random aggregates of spheres that can exist and study their properties. 

Many different approaches allow building random aggregates of spheres. Broadly, they may be divided into two categories. The first one consists of algorithms for which all spheres are introduced at once and then the system relaxes towards some more or less disordered structural state (eg molecular dynamics <cit.>, Lubachevsky-Stillinger algorithm <cit.>, Jodrey-Tory algorithm <cit.> or Monte-Carlo relaxation of chains of hard-sphere <cit.>). This family of algorithms produces aggregates whose properties mimic what was found in many disordered systems, notably the dependency of the contact coordination number with packing fraction and they are also able to produce aggregates with packing fraction equal or even superior to the random close packing (RCP) value (γ_RCP≈ 0.64) or equal to or lower than the random loose packing value (RLP, γ_RLP≈ 0.555 see <cit.>). However, γ_RLP presents a strong dependency on interparticles interactions <cit.> which suggests that no geometrical property must be associated to the RLP transition, whereas the RCP state might have structural or geometrical constraints that set an upper bound to γ_RCP, such as the proportion of spheres involved in quasi perfect tetrahedra, as observed by Anikeenko and Medvedev <cit.>.

The second large category of algorithms introduces spheres in the aggregate one by one and sets their position definitely at once (see, e.g., <cit.>). Aggregates built by such methods have two important differences with those of the first category: first, their average contact coordination number remains close to 6, second, to the best of our knowledge, it is impossible to produce aggregates with packing fraction higher than 0.6. Hence, the latter family of algorithms gives access to random aggregates with somewhat different structural properties from the former. Moreover, using sticky hard spheres, contact neighbours are rigorously defined by Dirac δ functions and this allows studying various properties related to the coordination number of spheres in a simple and non ambiguous manner.

In <cit.>, a first study of several hundreds sequentially built aggregates was conducted. It was found that for the lowest packing fractions, δ peaks appear on the pair distribution function, which corresponds to the formation of a disordered polytetrahedral structure in the random packing aggregate. 

The present study extends the results obtained in <cit.> by looking at structural tetrahedra, isotropy, contact coordination numbers, partial distribution functions as well as partial structure factors, 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 building tetrahedra, which appears essential to understand the composite nature of other sequentially built random aggregates.



§ STUDIED AGGREGATES

More than 10000 aggregates, each one containing 10^6 spheres with radius r_s=1 (and, hence, diameter d=2) were built and studied[Hereafter, unless otherwise mentioned, all distances are expressed in r_s as the unit of length (r_s = 1). Similarly, r_s^-1 is the unit in reciprocal space.]. They fall into two broad families: random irregular polytetrahedral aggregates (RIPA) and random regular polytetrahedral aggregates (RRPA). All of them were built by adding spheres one by one at their final position, tangentially to three existing ones.



 §.§ Random irregular polytetrahedral aggregates

Most aggregates were of RIPA type. They were built using algorithms that have been presented in details in <cit.> and are only summarized here. 
A seed of three spheres forming an equilateral triangle is used. Each new sphere P is positioned tangentially to at least three already present spheres (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 spheres at once (according to the index NINS, 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 packing fraction of the aggregates. The larger the value of a, the higher the packing fraction. It is varied between about 3.4 and 8 as, for a<3.4, no aggregate can be generated and for a>8, the maximum packing fraction is reached and no evolution of the generated aggregates is noted for higher values of a (see figure 2 in <cit.>).
More than 300 aggregates of 10^6 spheres were built, by varying a for each family of these algorithms. 

An additional modification with regard to the aggregates studied in <cit.> was to choose the origin O as close as possible from the center (0,0,0) of the growing aggregate (RMIN algorithms), instead of a purely random fashion. This change has the effect of increasing the maximum packing fraction, reached for (MAX-1, a>3.5), of about 1 %, from 0.586 to 0.593. All other aggregates generated by RMIN-MAX algorithms have a slightly higher packing fraction than their MAX counterparts (i.e. same values of a and NINS) with no significant changes concerning the structural results presented in <cit.>. This modification also allows for a more homogeneous growth of the aggregate. Only the results obtained for RMIN-MAX-1 will be used hereafter[Yet another family of algorithms has been developed, for which the newly added sphere is positioned by controlling the value of the building tetrahedra distortion. However, this strategy leads to the exact same P(r) as the one obtained for RA built with the positioning in function of the hole size and has the same maximum packing fraction: these aggregates will not be further studied here.].

Finally, as it turned out that the previous families of aggregates had some inhomogeneities of their packing fraction close to (0,0,0) (see hereafter), some aggregates were generated by using as 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 from an aggregate with the same value of other parameters (a and NINS), far from the origin (0,0,0), which has the effect to remove the central area with higher packing fraction. Typically, the number of spheres in the seed is between 30 and 600. However, this modification entails only slight changes of pair distribution functions or structure factors.



 §.§ Random regular polytetrahedral aggregates

A last algorithm has been used, that produces aggregates with only regular building tetrahedra <cit.>. In this case, the newly inserted sphere P forms a regular tetrahedron with the three already contacting spheres O, A and B, i.e. PO = PA = PB = OA = OB = AB = 2. For these algorithms, a has virtually no impact on packing fraction, but the number of inserted spheres around a given local origin (NINS) does. 

Once again, the local origin O can be chosen at random, or as the closest one from the center of the aggregate. When it is chosen randomly (RRPA), the maximum packing fraction (0.418) is reached for NINS = 3 and a minimum of 0.408 is observed for NINS = 1 and beyond 4. 
When O is chosen as the closest possible origin from the aggregate center (RMIN-RRPA), the maximum packing fraction is 0.456 and is reached for NINS>4. 
Similar aggregates, so-called saturated polytetrahedra, have been studied by Medvedev and Pilyugina <cit.>. They found a packing fraction of 0.435 for aggregates consisting of about 576 000 spheres. This value falls in between the maximum ones obtained for aggregates choosing O randomly (lower bound) and those taking O as the closest sphere from the origin (higher bound).




It should be noted that RIPA and RRPA distinguishes aggregates based on their building algorithms, not their structure, which will be discussed in details below.



The isotropy and the randomness of all aggregates has been systematically checked through the distribution 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.



§ PACKING FRACTION



 §.§ Basic relations

The radius R_m of a large spherical aggregate centered in (0,0,0) and made of N spheres centered in R⃗_i is given to a good approximation by <cit.>:

    = √(^2)

where  is the average quadratic radius of all spheres in the aggregate:

    ^2 = 5/31/N∑_i=1^N 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  of the sphere cut into the aggregate bulk 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 the origin (0,0,0) with another sphere , of radius , centered in (0,0,0) is <cit.>:

    V_s(r,)    = 4/3π r_s^3                                                    r≤-r_s
    
    V_s(r,)    = π(+r_s-r)^2(r^2+2rr_s-3r_s^2+2r+6 r_s-3^2)/12r   -r_s < r < +r_s
    
    V_s(r,)    = 0                                                                 r≥+r_s

These relationships can also be used for any sphere  centered in x, y and z by a mere change of reference frame. 

Hence, the packing fraction of the sphere  can be directly determined for any radius , as:

    γ () = ∑_i=1^N V_s(r_i,)/4/3π^3

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 the portion of spheres outside of the shell, according to the volume complementary of relation <ref>.



 §.§ Packing fraction of spheres inscribed in the aggregate

The evolution of the packing fraction of spheres inscribed in the aggregate as a function of their radius γ = f(/) (figure <ref>) shows that, for all aggregates, whatever their building algorithm, the packing fraction decreases slightly when  increases.

A seed effect appears as for MAX-1 and RMIN-MAX-1 aggregates (3 spheres forming an equilateral triangle as seed), a first regime of fast decrease is observed for < 0.2, i.e. for a number of spheres below approximately 8000, then the packing fraction tends to plateau whereas for N-RMIN-MAX-1 aggregates (seed consisting of spheres taken far from the origin in a previously generated aggregate), this initial decrease is much faster. Nevertheless, the exact range of effect of the 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 are extremely rare in high packing fraction aggregates. As a matter of fact, this seed dependency disappears 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 same behaviour for larger values of .



For every aggregate, a second regime is observed when →, i.e. when  reaches the limit of the aggregate: logically, the packing fraction decreases faster. For a perfectly spherical aggregate, when >+r_s, then the "packing fraction" of  should decrease as ^-3. When this decrease is at first more progressive, it shows that the aggregate has an imperfect shape and either 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 than the corresponding MAX-1 aggregates, as well as a sharper decrease of γ when →, which shows that RMIN-MAX-1 aggregates have a more regular surface than MAX-1 aggregates. For the latter, the thickness of "imperfect aggregate" is about 1.5d for a=3.5 (γ = 0.586) and about 5d for a=1.79 (γ = 0.370). For RMIN-MAX-1 aggregates, the thickness of imperfect aggregate is roughly d/4 for a=3.5 (γ = 0.593) and about 1.5d for a=1.78 (γ = 0.378).



 §.§ Packing fraction of shells



Figures <ref>.a and b represent the variation of packing fraction in shells with thickness w=0.1 for the densest aggregates produced by algorithms MAX-1, RMIN-MAX-1 and N-RMIN-MAX-1. Globally, the packing fraction of such shells oscillates with r. Its average value decreases from a higher value near the seed, to a smoother behaviour when r increases. N-RMIN-MAX-1 aggregate presents virtually no effect of the seed: the packing fraction of shells reaches the average behaviour for very small values of r, which seems logical as for these aggregates the seed consists of a set of spheres with the average structure. RMIN-MAX-1 and N-RMIN-MAX-1 converge for r∈ [0.1,0.15]: the effect of the initial equilateral-triangle seed of RMIN-MAX-1 aggregates seems then to act on about 10^3 spheres in the whole aggregate, consisting of 10^6 spheres, i.e. significantly less than suggested above by the comparison, in figure <ref>.b, of γ = f(/) for the two same aggregates. Figure <ref> 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 former aggregate has a brighter atypical central area in the region of the seed, denoting unusual structural properties with respect to the rest of the aggregate, whereas the latter displays a seed area much more similar to the rest of the aggregate.

Moreover, for aggregates with lower packing fraction (i.e. for aggregates built using lower values of a and in which regular polytetrahedra appear), the range of aggregate affected by the seed 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 seed ceases to be special in comparison with the rest of the structure. 
At large values of r (figure <ref>.b), oscillations can still be detected in the packing fraction of shells, however with a much smaller amplitude. A slight decrease of the average value is noticeable: the farther a shell is from the center of the aggregate, the lower its packing fraction, on average. The origin of this phenomenon is not obvious to us. For r>, the packing fraction of shells falls rapidly to 0.





§ TETRAHEDRAL STRUCTURE

The structure of random packings of spheres is commonly assessed via the tetrahedra connecting sphere centers, forming the so-called Delaunay tesselation <cit.> (these tetrahedra are noted by the subscript _D in what follows). For the present study, Delaunay tessellations were built using the cgal library <cit.>. 

In <cit.>, 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 new sphere P tangentially to the three other ones, O, A and B. 
It should be noted that such tetrahedra may or may not belong to Delaunay's tessellation.
The distortion of building tetrahedra was shown to be a very significant structural parameter, allowing the correlation of various structural traits of the aggregates. In this section, we focus on two distortion parameters of Delaunay tetrahedra.



 §.§ Distortion parameters

The first tetrahedral distortion parameter has been defined for the characterization of building tetrahedra <cit.> by relation:

    κ_BT = 3d^2 + OA^2 + OB^2 + AB^2/6d^2

where O, A and B are the three sphere centers used to add the new sphere P and the term 3d^2 corresponds to the three necessary sphere contacts PO, PA, PB imposed to building tetrahedra by the algorithm. The maximum value of κ_BT is 2 and is obtained for a centered equilateral triangle with side d√(3) of three spheres, with the additional sphere 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:

    = ∑_i∑_j>id_ij^2/6d^2

where i and j are the vertices of the tetrahedron and d_ij the vertices length. The smallest distance possible between sphere centers is d_ij=d=2, hence the smallest possible value is obtained for a regular tetrahedron and is = 1.

The last distortion parameter to be studied hereafter is the longest edge length () of the considered tetrahedron. The smallest possible value of , is d, which is found in the case of a regular tetrahedron. The behaviour of  has already been studied, along with others, notably by Anikeenko et al <cit.> on aggregates built using Jodrey-Tory (JT) algorithm <cit.> and Lubachevsky-Stillinger (LS) algorithm <cit.>. Anikeenko et al <cit.> have found that , in spite of its simplicity, shows a great consistency when compared with two other parameters, namely the edge differences and the procrustean distance.




 §.§ Distributions of distortion parameters




Globally, the distribution of distortion parameters suggest that two limiting aggregates exist, the densest one, produced by (RMIN)-MAX-1 algorithms (here for γ = 0.586) and the RRPA. The first limiting aggregate will be called fully random (FR) component, and the second one, regular polytetrahedra (RP) component. The notation FR and RP components refer to structural traits of the aggregates studied here. As a matter of fact, RRPA aggregates are fully RP, whereas RIPAs may share features of these two basic components in a variable proportion, depending on their packing fraction and building algorithm.
Figure <ref> shows local structures of (a) the RIPA aggregate with the smallest proportion of RP component and (b) a RRPA aggregate. The RRPA aggregate presents larger holes in the structure, however the polytetrahedral nature of the latter is quite difficult to visualize. On the other hand, it is clear that the RIPA with high packing fraction presents a much more homogeneous structure than the RRPA.





Figure <ref>.a presents various distributions of  obtained for aggregates built using MAX-1 algorithm (RIPA) and one RRPA, and figure <ref>.b represents the distribution of  for the same aggregates. Concerning MAX-1, these distributions present very similar behaviours: their maximum decreases with packing fraction while their full width at half maximum (FWHM) increases when packing fraction decreases. A bimodal component appears for the lowest packing fraction aggregates on both distributions, respectively centered around ≈ 2 and ≈ 3.5. 
These two values may be related in this way: assuming that tetrahedra with ≈ 3.5 have their 5 other edge lengths regularly distributed in the range [2 ; 3.5], leads to a  value of approximately 1.96, i.e. close to 2, suggesting that these two modes are indeed associated.
δ peaks appear for the lowest packing fraction, which is consistent with the existence of well defined distances observed on pair distribution function for the same aggregates, noting the existence of regular polytetrahedra and recurrent configurations of spheres in the structure (see <cit.>). 
The RRPA, on the other hand, appears as a limiting case. Indeed, the first peak observed for  and  distributions in the case of 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  present a discontinuity at ≈ 1.33, and the distributions of  have one at =2.827, which is likely associated. 
The distributions of the lowest packing-fraction aggregates have a change of slope for =4.
Additional work is needed to analyse precisely the configurations corresponding to these discontinuities.



 §.§ Average values of distortion parameters







Although average distortion parameters are basic structural parameters, they do not determine the packing fraction of the aggregate because they only involve the short range order of spheres (through the tetrahedral description) but do not take into account longer range order. It is therefore worthwhile studying the variation of packing fraction with distortion parameters for all aggregates.
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 high RP component in their structure and as RP means a regular tetrahedral basis.

Figure <ref>.a presents the dependency of  with packing fraction, which decreases when γ increases and turns out to be linear for each family of algorithm.
More specifically, the various algorithms are roughly distinguished as, for MAX-i aggregates, the average distortion increases with i for a given packing fraction. RMIN-MAX-1 stands on its own, with a slightly different slope. The same is observed for RAN-i aggregates, which span a larger interval of  for a given packing fraction (approximately 3 times as wide as that of MAX aggregates) and a certain overlap is observed as MAX-1 to 4 are between RAN-4 and RAN-6 for γ<0.5, which is the maximum packing fraction that RAN aggregates can reach. 

The behaviour of =f(γ) is globally the same (figure <ref>.b) as it also decreases when packing fraction increases. However, its behaviour deviates more from linearity: in the case of MAX-1 and RMIN-MAX-1 aggregates a change of slope is observed between 0.47 and 0.5 and they behave very differently from the rest of the other aggregates, with a higher value of  than any other aggregate for a given packing fraction. Furthermore, there is no overlap between RAN and MAX aggregates. RAN aggregates appear more dispersed than MAX aggregates and, for the latter, they more or less converge on the same curve (that of MAX-9) with the noticeable exception of MAX-1 aggregates.

Hence, surprisingly, the average distortions measured by both indicators do not agree as, for example, in the case of , RMIN-MAX-1 aggregates appear as the least distorted and as the most distorted according to .



 §.§ Proportion of regular and quasi regular Delaunay tetrahedra

Using either distortion parameter, it is possible to evaluate the volume fraction of regular Delaunay tetrahedra (Φ_V, with =1 or = 2), which globally decreases when packing fraction 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.073 for γ = 0.415. These values remain rather low and prompt the question: what is the geometrical upper bound of the volume fraction of regular tetrahedra in random aggregates? 

The results for random packings of RIPA type are presented in figure <ref>.a, which represents the variation of Φ_V with packing fraction. The highest proportion is obtained for RAN-6 aggregates (Φ_V = 0.141, γ = 0.422). RAN-1 turns out to be the family of aggregates with the smallest fraction of perfect tetrahedra for γ∈ [0.45;0.48] and outside of this interval, RMIN-MAX-1 outside of this interval. The proportion of perfect tetrahedra goes to 0 at the highest packing fractions.


Anikeenko and Medvedev <cit.> have studied the volume fraction of quasi regular Delaunay tetrahedra (PQRT, i.e. with <2.3) within aggregates generated with JT and LS algorithms (see figure <ref>.b of the present article). They have found that this volume proportion increases with packing fraction in their studied interval of packing fraction, i.e. roughly between 0.53 and 0.71. Figure <ref>.b superimposes their results with the ones found in the present study. 
Interestingly, the curves for MAX-4 and MAX-3 match over a rather narrow packing fraction interval, around γ=0.55, and all MAX-i and RMIN-MAX aggregates show a similar increase of the volume fraction of tetrahedra with <2.3 for γ∈ [0.56;0.59], however with lower satisfying quantitative agreement. 
The highest value of the volume fraction of quasi regular tetrahedra found here is 0.127 and is obtained for RMIN-RRPA aggregates.

This comparison shows that regular and quasi regular tetrahedra are two distinct populations: Φ_V decreases when packing fraction increases and goes to 0 beyond a threshold that depends on the algorithm but is roughly γ = 0.57, whereas PQRT goes to a minimum for a given packing fraction that depends on the algorithm (between γ = 0.55 and γ = 0.57) and then increases with packing fraction beyond this value of γ.






§ EFFECT OF SPHERE COORDINATION


The contact coordination number (CCN) of each sphere is determined when the aggregate is built: 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.



 §.§ Contact coordination number



  §.§.§ Partial distributions of contact coordination numbers
 
Let η_ij be the number of contacts between spheres with CCN i and j respectively. 
The (normalized) distributions of η_ij can be easily determined from the sphere positions for all values of i and j. 
Figure <ref> introduces distributions of η_ij for three aggregates built by RMIN-MAX-1 algorithms, from the highest to the lowest packing fraction (fig. <ref>.a to c) and a RRPA aggregate (with the highest packing fraction among RRPA, fig. <ref>.d). 



These distributions show a progressive shift of the maxima from high to low packing fraction (fig. <ref>.a to c). At high packing fraction, spheres with low CCN tend to be surrounded by spheres with higher CCN thus reducing local fluctuations of the CCN. Then, as the packing fraction decreases, all η_ij curves more or less collapse, meaning that in this regime, all spheres–whatever their contact coordination number η_i–have the same η_ij distribution, centered on the average contact coordination number, hence a very similar environment in terms of contact neighbours. Finally, for even lower packing fraction, an inversion is observed and high CCN spheres are preferentially surrounded by high CCN spheres, which corresponds to a contact segregation effect. 
At the same time, the FWHM of the distributions widen as the packing fraction decreases and they become less symmetrical.
For the lowest packing fraction (fig. <ref>.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. <ref>.d) as low packing fraction RIPAs, going through a minimum for η_i5. However, the order of the various η_ij distributions appears inverted in the case of the RRPA as, on the high j end, the spheres with the highest proportion of contact with high CCN spheres are spheres with lower CCN (i), i.e. η_3,11>η_4,11>...>η_12,11, whereas for RIPA with low packing fraction, the opposite situation is observed, i.e. η_12,11>η_11,11>...>η_3,11 (fig. <ref>.c).




  §.§.§ Evolution of  around i sphere

, the average value of η_ij around spheres with CCN i, can easily be determined from the distributions presented above. Figure <ref> presents the evolution of  for all values of i as a function of packing fraction for MAX-1 and RMIN-MAX-1 algorithms. The inversion suggested by the shift of the maximum of distributions in the previous section appears clearly for RMIN-MAX-1 aggregates (figure <ref>.b) as they form a crossover for γ≈ 0.52 which separates a low and a high packing fraction regimes. In the low packing fraction regime, high coordination spheres tend to be surrounded by spheres with higher CCN than the spheres surrounding low CCN sphere i.e. ⟨η_12,j⟩  >  ⟨η_11,j⟩  > ... >  ⟨η_4,j⟩, with the exception of the limit case ⟨η_3,j⟩, corresponding to the 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 with higher CCN, thus reducing density fluctuations.

In the case of aggregates generated by MAX-1 algorithm, this crossover is not captured but might take place at higher values of packing fraction – unaccessible by this algorithm – as all curves begin to collapse for γ>0.55. In this case, only the low packing fraction regime is observed.





  §.§.§ Radial dependency of the average contact coordination number


The radial dependency of the average contact coordination number ⟨ CCN ⟩ of spheres within [r;r+dr] from an i coordinated sphere has been determined for all aggregates and all values of i. Figure <ref> introduces two examples obtained for the RMIN-MAX-1 aggregates with the highest and lowest packing fractions (figure <ref>.a and b).




At high packing fraction (see figure <ref>.a), the average value of the CCN of spheres at a distance r<3.5 from a low coordinated sphere is higher than that of spheres surrounding a sphere with higher CCN. Then, as the packing fraction decreases, the relative positions of the various curves for r<3.5 are progressively inverted: they superimpose for γ≈ 0.54, which incidentally matches the packing fraction of the crossover of the various  curves in figure <ref>.b for the same aggregates. At even lower packing fraction, the inversion is complete as it is exemplified in figure <ref>.b: low coordinated sphere are, on average, surrounded by quasi-first-neighbour spheres with low CCN and vice versa. This regime corresponds to a segregation effect in the range of quasi first neighbours instead of contact neighbours. At even larger values of r (r>4), this segregation effect gets inverted.

All curves end up superimposing one another beyond some value r=r_e, however r_e increases when packing fraction decreases (from r_e≈ 4 for γ=0.593 up to r_e≈ 6.5 for γ=0.378), suggesting that structural inhomogeneities extend over larger and larger scales as γ decreases.

Discontinuities are observed at r=d√(3) and r=2d, which matches discontinuities of the various PDF of the same aggregates (cf. section <ref> below and <cit.>). The amplitude of these discontinuities increases with the CCN of the central sphere.

At low packing fraction (fig <ref>.b), strong local increases of ⟨ CCN ⟩ are noticeable for values of r which matches δ-peaks on the PDF, hence corresponding to distances characteristics of the presence of regular polytetrahedra in the aggregate (cf. section <ref> below). 



 §.§ Bond angle distributions


Two spheres i and j form a bond when they are contact neighbours. The bond angles α around a sphere i are defined as the angles formed between all possible vectors R⃗_ij between contacting neighbours. Hence, for a sphere i with n contact neighbours there are n(n-1)/2 bond angles.



  §.§.§ Global bond angle distributions

Bond angle distributions α have been calculated. The smallest possible bond angle is between three contacting spheres, i.e. α = π/3, and the largest is, of course, π.

Figure <ref>.a presents bond angle distributions for the most and least dense aggregates obtained by 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 to the configuration of three spheres in contact with each other. The value here is so high that it is outside the scale of the figure. The second discontinuity occurs for α = 2.093± 0.002, close to 2π/3, which corresponds to the situation where four spheres form two coplanar equilateral triangles sharing a common side.
In <cit.> Karayiannis et al obtained a very similar bond distribution in their structures 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 <cit.>, which we interpret as a finite size broadening since the maximum number of spheres considered in <cit.> is 54000 instead of 10^6 used here.





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 the presence of regular polytetrahedra in the disordered structure (see <cit.>). The addition of another sphere on top of three contacting ones does not introduce a new bond-angle. The addition of a fifth sphere forms then a trigonal bi-pyramid, with α≈ 1.9106 rad. This angle should correlate with the δ-peak at r=d√(8/3) observed on the pair distribution function of low density aggregates (see <cit.>). It shows the same smooth dependency with the irregularity index of building tetrahedra (see figure <ref>.b) as P(r=d√(8/3)).



  §.§.§ Partial bond angle distributions


Bond angle distributions for specific contact coordination numbers are shown in figure <ref>. 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 depends strongly on the CCN. Low CCN spheres have a stronger asymmetry between α = π/3 and α = 2π/3 than spheres with high CCN. For high packing fraction aggregates, distribution of low CCN spheres show a depletion of lower bond angles and an excess of high angles and the distributions get more even when the CCN increases. 
When the packing fraction decreases, low and high CCN spheres tend to have more similar angular environment.

The appearance of regular polytetrahedra, associated with δ peaks, like in the case of the global bond angle distribution, is logically correlated with a decrease of the continuum component of the distribution. This continuum disappears completely for the RRPA (figure <ref>.d).



 §.§ Partial pair distribution functions

The structure of random sphere packings is usually characterized by the probability per unit volume 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 structural 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 sphere with contact coordination number i, normalized to 1 at large r.



  §.§.§ 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 a spherical aggregate with radius R, P_ij(r) write:

    P_ij(r) = V^2/N_iN_jΔ N_ij/(2-δ_ij)S(r)Δ r

where N_i and N_j are the number of spheres with CCN i and j, respectively, in the volume V=4π R^3/3 of the aggregate; 
Δ N_ij is the number of sphere pairs of CCN i and j respectively, lying in the interval [r : r+Δ r];
δ_ij is the usual Kronecker symbol;
S(r) is the spherical shape factor of the aggregate <cit.>:

    S(r) = π^2/6r^2(2R-r)^2(4R+r)


P_ii(r) PPDF describe the arrangement of i coordinated spheres, while P_ij(r) PPDF with i≠ 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 <cit.>:

    P_ij(d) = η̅_ijV/4π N_j d^2δ(r-d) i.e. numerically η_ij = 3N_j/R^2d^2 2σ P(d)

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 range of the P_ij(r) figures presented below.

One can then define the probability P_i(r), to find a sphere with any coordination number at a distance r from a sphere with coordination number i. It writes:

    P_i(r) = ∑_j=3^12 C_j P_ij(r)

where C_j = N_j/N is the concentration of j coordinated spheres (N=∑_j=3^12N_j).

These P_i(r) characterize the global arrangement of spheres around a sphere with coordination i.

Finally, the global PDF writes:

    P(r) = ∑_i=3^12C_i P_i(r)





  
Random regular polytetrahedral aggregates
Figure <ref> presents P_ii(r) curves (with i=4, 6 and 8) of an RRPA with γ = 0.418, as well as its global PDF (see <cit.> for the calculation procedure of the latter).




Concerning the global PDF, a striking difference from what was observed for RIPA (see <cit.>), is the disappearance of the topological discontinuities at r = √(3)d and r=2d. The continuous structure of the PDF observed in RIPAs almost disappears and is replaced by a set of polyhedral δ peaks (<cit.>), which are due to a large (virtually infinite) regular polytetrahedron. The positions of these peaks, which are due to precise configurations of spheres with well defined distances, are identical to those studied by Medvedev et al <cit.>. The continuous structure grows almost linearly as a function of r and goes to 1 for r≈ 3.5d = 7.

The partial PDFs (P_ii(r)) show a systematic tendency, as their continuous regime in the region of the quasi first neighbours (QFN) starts from 0 for spheres with coordination number i equal or superior to 8 and increases like the global P(r). The lower coordination numbers (i<8) lead to P_ii(r) starting from higher values in the QFN area. This can be qualitatively understood: the more contacting first neighbours a sphere has, the less quasi first neighbours it can accept.


  
Random irregular polytetrahedral aggregates

The global PDF of RIPA was studied in <cit.>, with the exception of RMIN aggregates. However, the latter bring no qualitative differences to these results. 



    P_ii(r)
 The PPDF P_ii(r) together with the corresponding PDF are presented in figures <ref>.a and b for MAX-1 aggregates with the two most extreme packing fraction, γ = 0.586 and γ = 0.370.
First and foremost, for a given packing fraction, the number of quasi first neighbours (corresponding to values of r close to d) decreases when the coordination number (i) increases, like in the case of RRPA and for the same reasons: the most contacting neighbours a sphere has, the less quasi first neighbours it can accept. 




On the other hand, the comparison between figure <ref>.a and b, shows that the number of quasi first neighbours decreases with packing fraction. At low packing fractions, they form plateaus whose level increases when the coordination number decreases, whereas they have a very distinct behaviour for high packing fractions: low coordination number spheres possess a lot of quasi first neighbours while high coordination number ones have a very limited number of QFN.

Both topological discontinuities (at r = √(3)d and r=2d) increase when the coordination number increases and when packing fraction increases.

Finally, (see supplementary informations for figures) "iso-packing-fraction" and "iso-" aggregates obtained by different building algorithms show significant differences in their global and partial PDF, confirming that parameter γ and  are greatly insufficient for a full structural description of random aggregates. Aggregates sharing similar values of these two parameters (γ and ) still present significant differences in structural properties. Hence, even when used together, γ and  are not satisfying predictors of structural properties of disordered systems. Besides, it has been impossible to find aggregates with similar  and γ to compare them in a similar fashion leaving the question open for a possible combination of these two parameters as good predictor of random aggregates global structure. Conversely, it turns out that PPDFs allow a more sensitive distinction between aggregates built by different algorithms than the corresponding global PDF and are thus interesting structural descriptors.




    P_ij(r) with i≠ j
 
A sampling of P_ij(r) curves is shown in supplementary information. They show 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,6 and 9) obtained for aggregates built by MAX-1 algorithm with packing fraction 0.370 and 0.586 are displayed in figure <ref>.
These PPDF appear to depend strongly on 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√(3) and 2d are more intense for high coordination numbers;





On the other hand, the main differences between low and high packing fraction aggregates are that:


  * 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 quasi first 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 the lowest packing fraction aggregate.




 §.§ Radial evolution of local packing fraction from pair distribution functions



  §.§.§ Principle

Knowing pair distribution functions, it becomes possible to study the variation of the local packing fraction around an average sphere as a function of r.
Using relation <ref>, the packing fraction of a sphere  with arbitrary diameter >d situated within the aggregate writes:

    γ () = 4π r_s^3/3 + ρ∫_0^+r_s P(r)4π r^2 V_s(r,) dr/4π^3/3

with ρ=N/V, the number of spheres per unit volume in the aggregate;
the 4π r_s^3/3 term corresponds to the sphere in r=0; V_s(r)  is defined by equation <ref>; P(r) can be a total or partial pair distribution function. It also comes γ (r≤ r_s)=1.

It is also possible to remove the contribution of the central sphere and its contact first neighbours to packing fraction, which writes: 

    γ_CN() = 4π r_s^3/3 + iV_s(r=d,)/4π^3/3

where i is the number of contacting neighbours of the central sphere. 
It then comes:

    γ_WCN() = γ () - γ_CN()

γ_WCN is the contribution to packing fraction of quasi-first and further neighbours around an average sphere.



  §.§.§ Local packing fraction around average sphere


Using the global PDF, the formalism just introduced gives access to the variation of the local packing fraction as a function of r around an average sphere. Figure <ref>.a represents γ() and figure <ref>.b represents γ()/γ, which goes to 1 as  goes to ∞. All aggregates behave differently depending on their packing fraction. For the densest ones, the packing fraction falls below the average value, before converging more rapidly at large r than the least dense ones. On the other hand, in the case of the least dense ones, packing fraction remains above the average value and converges more slowly. 
Oscillations of local packing fractions are damped at about r=5d for the lowest packing fraction and around r=3d for the highest packing fraction.






  §.§.§ 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 an average sphere with coordination number i, γ_i.
Figures <ref>.c and d present the sets of corresponding curves 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. Their amplitude is damped for the aggregate with the lowest packing fraction. Figures <ref>.e and f show the evolution of packing fraction without the contribution of contact neighbours (γ_WCN) and the central 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 numbers leads to a significantly lower packing fraction due to QFN than in the case of spheres with low contact coordination number. 

In the case of the aggregate with the highest packing fraction (fig. <ref>.e), the various curves converge for ≈ 3.5d and the highest difference is observed between i=10 and i=3 for ≈ 2.53, γ_WCN,3-γ_WCN,10≈0.28 = 0.48γ_ag.  Between i=10 and i=6, the maximum is found at ≈ 2.57, with γ_WCN,6-γ_WCN,10≈0.16 = 0.27γ_ag.

There is a much lower difference in contribution of the QFN between spheres with different contact coordination number in the case of the aggregate with the lowest packing fraction (in terms of absolute as well as relative value). However, the maximum of the difference between γ_WCN,3 and γ_WCN,10 is obtained for approximately the same value ≈ 2.55 and γ_WCN,3-γ_WCN,10≈0.06 = 0.16γ_ag. Between i=10 and i=6, the maximum is found at ≈ 2.51, with γ_WCN,6-γ_WCN,10≈0.06 = 0.16γ_ag, which is remarkably similar to the former. In this case, figure <ref>.e shows that γ_WCN() for coordination numbers below 7 behave very similarly, suggesting that these types of spheres have more or less the same environment in terms of number of QFN, irrespective of their contact coordination number. 



 §.§ Global and partial structure factors

The previous structural description  achieved by the PPDF analysis of random packings leads naturally to the study of the corresponding global and partial structure factors which can be compared with the results of diffraction experiments on disordered materials.



  §.§.§ Global and Ashcroft Langreth partial structure factors

Let us first introduce the global structure factor:

    S(Q) = 1 + N/V∫_0^∞(P(r)-1) sin(Qr)/Qr4π r^2 dr

and the partial Ashcroft Langreth (AL) structure factors, defined by:

    S_ij(Q) = δ_ij + √(N_iN_j)/V∫_0^∞(P_ij(r)-1)sin (Qr)/Qr4π r^2 d r

where δ_ij is the Kronecker symbol. 

The diagonal terms, S_ii, represent the structure factors of the partial aggregates formed by spheres with contact coordination i and will be studied first. The non diagonal terms S_ij with i≠ j describe the mutual or "chemical" arrangement between i and j coordinated spheres and will be studied later by the more suited Bhatia-Thornton formalism. 

S(Q) and S_ij(Q) values for Q ≪ 1/R cannot be calculated from relations <ref> and <ref> owing to the finite radius R of the aggregates <cit.>. However S(Q=0) values are determined by the statistical density fluctuations according to relation <cit.>:

    S_ij (0) = N_iN_j - N_i N_j/√(N_i N_j)

where upper bars represent averages.

Fluctuations were directly derived from the positions of the sphere centers. To that purpose, a large cube with edge 50d centered on the aggregates origin has been subdivided into 1000 sub-cubes with edge 5d, each of them containing N_i and N_j sphere centers with respective coordination i and j. Average values of N_i, N_j and N_iN_j were taken over the 1000 sub-cubes to get the statistical fluctuations involved in relation <ref>.



  
Random regular polytetrahedral packing – global structure factor
Figure <ref> introduces the global structure factor of two RRPA: <ref>.a (small Q) and b (large Q or asymptotic regime).
It turns out that the structure factors of RRPA exhibit a prepeak, whose position depends on the used algorithm or packing fraction. The position and intensity of this pre-peak corresponds to the folding of the polytetrahedra which controls the space correlation of density fluctuations. It accounts 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) corresponding to the distances between vertices in the infinite RRPA and is only weakly affected by small changes in packing fraction (depending on the used algorithm).







    Comparison with the structure of tetravalent elements

The polytetrahedral structure factor can be usefully compared with the experimentally measured structure factors of pure tetravalent elements (C <cit.>, Ge <cit.>, Si <cit.>) or with the tetrahedral structure of oxygen atoms in glassy water <cit.>. 
As a matter of fact the structure of these elements is based on body centered regular tetrahedra (BCT) connected by their vertices <cit.>, while RRPA's are made from simple regular tetrahedra connected by their faces. 
As a consequence, a "pseudo superstructure" peak for such hard-sphere based centered-tetrahedra should correspond to Q_s/Q_1 = d_1/d_s = √(6)/4 = 0.612, where d_1 is the distance between contacting atoms and d_s is the edge length of the regular centered tetrahedron. 
Experimentally, Etherington finds Q_s/Q_1=0.572 <cit.> for the structure factor of amorphous Ge, Laaziri et al. find 0.583 for the same ratio in a-Si <cit.> and Gaskell et al. find 0.540 for a-C <cit.>. 
The discrepancy from 0.612 for these three values could be due to the error introduced by hard sphere approximation of tetrahedral bond angle between first neighbours of softer 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 another small angle peak at significantly lower value is found around Q_s/Q_1 ≈ 0.25, varying with packing fraction. The latter can be attributed to the pseudo periodic correlation of density fluctuations present in the RRPA model. These correlations happen over larger distances than those of the BCT, hence the smaller Q values of the corresponding pre-peak. It corresponds roughly to the end of significant oscillations of pair distribution functions.

The reason why this longer correlation-distances pre-peak does not appear for tetravalent elements has probably to do with the fact that tetrahedra forming their structures are connected by vertices, allowing for many more degrees of freedom and preventing longer distances correlations than in the case of RRPA aggregates where tetrahedra are connected by their faces, which is geometrically more constraining. 
On the other hand, the second peak of the structure factor of the RRPA splits into two subpeaks in a completely analogous manner with what is found in amorphous tetravalent elements. 

Finally, for RRPA aggregates, characteristic distances are strictly defined, giving rise to δ peaks on the PDF resulting in aperiodic large Q behaviour, also absent from the structure factor of tetravalent elements. In contrast, the lack of such well defined distances, because of higher degrees of freedom in the orientations of tetrahedra relatively to one another in the case of tetravalent element, results in a single δ peak on their PDF, which accounts for their periodic behaviour in the large Q regime <cit.>.





  
Random irregular polytetrahedral aggregates

The AL partial structure factors of random irregular polytetrahedral aggregates are presented in figure <ref>.a to c, along with the partial structure factors of the RRPA (<ref>.d). 
These structure factors all have in common that the global amplitude of their oscillations increases from i=3 up to i=6 and then decreases as i increases beyond 6.
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 like sin (Qd)/Qd (not shown here).

Some peculiarities are also observed. A shoulder appears on the left of the first peak of S_44 (around Q=3.8). 
For S_88 and S_44 an intense pre-peak around Q=0.8 is noted, whose intensity increases with the packing fraction for S_44 and decreases when the packing fraction increases for S_88. This pre-peak is due to correlated density fluctuations between low or high density regions separated by an average distance of 2π/0.8 r_s ≈ 8r_s lying beyond the main P_ii oscillations; they are observed whatever the packing fraction, even if there is no noticeable pre-peak on the global structure factor.
In the case of S_66, S_66(0) decreases when γ increases and barely any pre-peak may be observed at all. The first peak remains symmetrical; its position depends little on packing fraction. The first peak of S_88 becomes strongly asymmetrical for the highest packing fraction. 

Finally, a splitting of the second peak is observed on the AL structure factors with high coordination numbers and low packing fraction. This splitting has been noticed above for the RRPA and it is consistent with the existence of a RP structure embedded in a "continuum" random structure (respectively RP and FR structural components), as was already concluded in <cit.> for these aggregates.






  §.§.§ Bhatia Thornton partial structure factor

The partial structure factors introduced by Bhatia and Thornton <cit.> are associated with density and concentration correlations in binary alloys. They have been extended to alloys consisting of more than two components <cit.> and can then be used by considering random aggregates as alloys of spheres with various coordination numbers.



  
Formalism
The partial structure factor corresponding to:


  * density-density correlation function, S_NN, writes:

    S_NN(Q) = 1+ρ∫_0^∞(P(r) - 1) sin(Qr)/Qr 4π r^2 dr


  * concentration-concentration correlation function, S_C_iC_j, writes:

    S_C_iC_j(Q) = 1+ρ∫_0^∞(-P(r) + P_i(r) + P_j(r) - P_ij (r)) sin (Qr)/Qr4π r^2 d r


  * density-concentration, S_NC_i, writes <cit.>:

    S_NC_i(Q) = ρ∫_0^∞(P_i(r) - P(r)) sin (Qr)/Qr4π r^2 d r


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 supplementary informations. 



  
Results
Between 15 (for high packing fraction) and 45 sets of Bhatia-Thornton structure factors corresponding to different coordination pairs i-j were calculated for 500 aggregates. In the 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 and 8, hence the lower number of BT structure factors calculated for them.

A sample of Bhatia Thornton partial structure factors is provided in supplementary informations for three MAX-1 aggregates, with maximum, minimum and intermediate packing fraction along with those of one RRPA aggregate.


Their main characteristics are the following:



    S_N_iN_j
 The S_N_iN_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≥ 4)



    S_NC_i
 structure factors exhibit significant oscillations whose intensity increases with the coordination difference i-j. Furthermore, for low coordination differences these oscillations decrease when the coordination numbers increase. All these oscillations come from the fact that the average environment of sphere i and j characterized by P_i and P_j differ from the average global environments P(r) (see relation <ref>). This behaviour differs completely from the case of binary "substitution" alloys (with equal atomic diameters) for which S_NC does not oscillate at all because P_i(r)≈ P(r), independently of the chemical order of the two alloy components (which only affects S_CC <cit.>).




    S_C_iC_j
 Most interesting is the study of S_C_iC_j, which characterizes the "chemical" order between i and j coordinated spheres through its dependence on P_ii-P_ij (see equation <ref>). For small coordination differences (at most 2) and all packing fractions, S_C_iC_j oscillates weakly around 1 and there 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 ≈ 0.5 Q_1. 
This pre-deep indicates a segregation effect <cit.>, between spheres with large CCN differences. It varies non uniformly with packing fraction (see figure <ref>.a) and its maximum amplitude lies in the packing fraction interval γ∈[0.5;0.55].
This segregation effect occurs over distances larger than contact neighbours, as it was suggested by the radial dependency of  seen in section <ref>. As a matter of fact, the maximum segregation measured by this pre-deep occurs in a range of packing fraction where there is virtually no contact segregation in the case of RMIN-MAX-1 aggregates (by comparing fig <ref>.b and fig <ref>.a), suggesting that there can exist contact and longer-range segregation effects between spheres of various coordination numbers and that these two types of segregation are not necessarily concomitant.




Furthermore, for these high coordination differences a new S_CC prepeak appears at very small Q values, around 0.2 Q_1, which seems to vanish in the aggregates with the highest packing fraction (fig <ref>.b), with the notable exception of MAX-4. Accordingly, this S_CC prepeak is also observed in the RRPA.
This peak should be associated with the existence of segregation effects extending to larger r range, i.e. to the formation of isocoordinated aggregates in an "average matrix". This suggests the existence of regions of low or high coordination in the aggregate matrix. 



§ 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 <cit.> with contact coordination number 2 <cit.>, 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 <cit.>.

In this paper we focused on the detailed structural characterization of such packings, using a wide variety of static aggregates-building algorithms. 

The randomness or isotropy of the 10^4 aggregates (containing 10^6 spheres each one) was checked first by studying  the angular distribution of pairs of spheres separated by a given distance. The eigenvalues of the corresponding nematic order tensor were shown to be always less than 2.10^-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 the fifth neighbour range, i.e. are limited to the first 8 000 spheres, if the seed consists of sphere arrangements that are rare in the aggregate. This seed effect can be totally removed by using as seed 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 the aggregates depends on the building algorithm and packing fraction, and is limited to a thickness of approximately 6r_s in the "worst case". Finally, the accuracy on the (average) packing fraction increases with the sphere number of the aggregate which must reach  10^6 in order to get a 10^-3 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 tetrahedra connecting sphere centers. The distributions of these tetrahedra have been characterized by two distortion parameters,  and . Their distributions show a bimodal character that varies with packing fraction and appear related, they both include discontinuities whose origin remains unclear. Their average values decrease with increasing packing fraction. 
However, they provide complementary characterization of the aggregates produced by the various studied algorithms, suggesting that beyond their similarities, these distributions also present subtle differences that will be the object of a future study. 

Special attention was paid to the populations of regular tetrahedra (formed by 4 mutually contacting spheres, =2) and quasi regular Delaunay tetrahedra (< 2.3) which were shown to behave quite differently. As a matter of fact, the volume fraction of regular tetrahedra has been shown to decrease with increasing packing fraction and reaches a maximum value of 0.165 for RPPA aggregates (only built with regular construction tetrahedra). 
This raises a fundamental question: 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 regular Delaunay tetrahedra goes to a minimum around γ = 0.56 and then increases with increasing packing fraction.

New structural characterization methods could then be introduced by taking advantage of the unequivocal definition of sphere contacts and, hence, contact coordination numbers.

First, partial characterization was carried for short distances, i.e. contacting neighbours. The distributions of pairs of spheres with respective CCN i and j (η_ij) were first studied. 
Their FWHM increases when packing fraction decreases and their average values () show distinct behaviours with respect to γ for low and high values of i respectively: a non uniform decrease with γ is obtained only in the case of low i values.
The evolution of these distributions shows that contact segregation can exist between spheres of various CCN, whereas the evolution of  for spheres of various coordination shows that another segregation may exist, on the basis of CCN, over a larger 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 ones present a depletion of low angle bonds and have a higher proportion of high angle bonds. These differences are particularly distinct in high packing fraction aggregates.

This structural description was then extended to all distances by introducing the partial pair distribution 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 and this detailed analysis allowed a clear cut distinction between different random packing structures which cannot be distinguished by their packing fraction and/or  parameters. The question remains of whether these aggregates could be characterized by characterizing them through the combination of packing fraction γ and distortion parameter . 

Distinguishing (with Bernal <cit.>) contacting from quasi-contacting spheres, it was shown that high CCN spheres have few quasi contacting neighbours (QCN) while low CCN spheres have a higher 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 Fourier transforms of linear combinations of the P_ij(r) and can be directly compared with the results of 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 and the variable shapes of these S_ii, especially their prepeak, first and second Q oscillations, were described. Furthermore the  concentration-concentration partial structure factor of Bhatia gave the mutual or “chemical” order between spheres with different coordinations i and j and put forward hetero-attractions or segregation of spheres within the different aggregates, confirming that several kind of segregation may indeed exist.

From all these results we could conclude that the random irregular polytetrahedral aggregates studied here are made from two basic components, namely fully random component (without regular tetrahedra) and regular polytetrahedral component (only built with regular tetrahedra), whose proportion decreases with increasing packing fraction. This composite nature of the aggregates produces “prepeaks” (at very small Q) in the aggregate structure factors and the main features differentiating 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 parameters  and 

  * a spreading of contact coordination numbers distributions η_ij

  * strong variations in the radial dependency of the average contact coordination number of spheres

  * a collapse of the continuum of partial bond angle distributions

  * a more similar quasi-first-neighbours part of partial pair distribution functions of low and high contact coordination number spheres

  * a more similar radial dependency of packing fraction for spheres of any contact coordination number

  * a reduction of topological discontinuities in pair distribution function at r=2√(3) and 4

  * a splitting of the structure factor's second peak (Q≈ 7), similar to what is experimentally found in amorphous tetrahedral materials. 


The present results on statically built sticky hard sphere aggregates cover the packing fraction interval 0.370-0.593. It could be interesting to extend this study:

Firstly to the low packing fraction range 0.15 (condensation or percolation limit)-0.370, which involves contact coordination number lower than 3 and cannot be reached by the building algorithms used here.

Secondly, to the high packing fraction range 0.593-0.64 which needs a reduction of the fluctuations of the local packing fraction at all length scales and cannot be reached by our algorithms which only minimize the fluctuations of the local packing fraction up to the second neighbours distances.
In particular, it would be of great interest to study dynamically built aggregates, such as the ones produced by Lubachevsky-Stillinger algorithm, Jodrey-Tory algorithm or molecular dynamics.
The first problem to be resolved would then be to define the distance corresponding to contacting spheres, as, for these aggregates, the first neighbour peaks are spread over a finite length interval. Nonetheless, valuable structural informations might be at hand through the techniques studied in the present work, that might shed new light on various long standing questions concerning random aggregates as well as allowing new comparisons between statically and dynamically built random aggregates.




elsarticle-num




§ ACKNOWLEDGMENTS

P. Cénédèse 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 process of writing the present paper.

Anonymous referees are greatly acknowledged for their critics and suggestions that helped improving the quality of this paper.





§ SUPPLEMENTARY INFORMATION




 §.§ 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 convenient check for the building method and/or the minimum aggregate size beyond which the isotropy is reached. 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(θ)/2. 

It is worth mentioning that a lattice like distribution of spheres is characterized by a θ–distribution presenting distinguishable peaks; thus the features of the θ distribution appears to be 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 from isotropy by looking if a favoured direction emerges among the {r̂_ij} distribution. For this, we follow the usual method of the framework of the nematic liquid statistics <cit.> according to which, from the diagonalization of the two rank tensor

    Q̅ = 1/N∑_ij1/2(3r̂_ijr̂_ij - I̅)

where I̅ is the identity tensor, a nematic order parameter, say S_1 is obtained as the largest eigenvalue λ_max of Q̅.
In equation (<ref>), the sum running over all the {i,j} pairs, can be limited to the sum over the bonds r̂_ij, which can be replaced by a sum over the unit vectors carried by the r⃗_i, with respect to a fixed reference r⃗_o, which avoids the surface effects when both i and j are located at the aggregate surface.

The results for the polar angle distributions relative to the axis x̂, ŷ and ẑ for the most and least dense aggregates obtained by MAX-1 algorithm are presented respectively in figures <ref>.a and <ref>.b.
For the aggregates studied in this work, the deviation from isotropy as measured by the value of λ_max is found to decrease when the packing fraction increases. A reliable determination of the dependence of λ_max with respect to γ is beyond the scope of this paper, all the more that the range of λ_max values is quite small.
The largest eigenvalue of Q̅ is smaller than 3 10^-2. Such a value is very small leading us to conclude that these building methods lead indeed to isotropic aggregates. 
To illustrate quantitatively this point, figure <ref>.c displays the result of the θ distribution in terms of λ for the following model:

    P_b(θ) = sin(θ)/2[ e^(-θ^2/(2σ^2)) +
    e^(-(π-θ^2)/(2σ^2))]

where the polar axis is a favoured direction according to the variance σ. Clearly, the θ distributions characterized by λ values lower than a few 10^-2 can be considered isotropic.







 §.§ P_ii(r)






 §.§ P_ij(r)





 §.§ AL structure factors






 §.§ Bhatia-Thornton formalism for binary mixtures of spheres

General BT relations then reduce to:

    S_N_iN_j(Q) = 1 + ρ' ∫(C'_i^2 P_ii + C'_j^2 P_jj + 2C'_iC'_j P_ij -1) sin (Qr)/Qr4π r^2 d r

where ρ' = (N_i+N_j)/V and C'_i = N_i/(N_i+N_j).

    S_C'_iC'_j(Q) = 1 + ρ' C'_iC'_j ∫[(P_ii(r)-P_ij(r))+(P_jj(r)-P_ij(r))]sin (Qr)/Qr4π r^2 d r


    S_NC'_i(Q) = ρ' C'_j ∫(C'_iP_ii(r) - C'_jP_jj(r) + (C'_j-C'_i)P_ij(r))sin (Qr)/Qr4π r^2 d r









