<figure><img src="image_1.png" /></figure>

# Revealing the three-dimensional structure of liquids using four-point correlation functions

Zhen Zhangª and Walter Kobª,ᵇ,ᵇ

ªLaboratoire Charles Coulomb, University of Montpellier, CNRS, F-34095 Montpellier, France; and ᵇInstitut Universitaire de France, 75231 Paris Cedex 05, France

Edited by Michael L. Klein, Temple University, Philadelphia, PA, and approved May 4, 2020 (received for review March 25, 2020)

Disordered systems like liquids, gels, glasses, or granular materials are not only ubiquitous in daily life and in industrial applications, but they are also crucial for the mechanical stability of cells or the transport of chemical and biological agents in living organisms. Despite the importance of these systems, their microscopic structure is understood only on a rudimentary level, thus in stark contrast to the case of gases and crystals. Since scattering experiments and analytical calculations usually give only structural information that is spherically averaged, the three-dimensional (3D) structure of disordered systems is basically unknown. Here, we introduce a simple method that allows probing of the 3D structure of such systems. Using computer simulations, we find that hard sphere-like liquids have on intermediate and large scales a simple structural order given by alternating layers with icosahedral and dodecahedral symmetries, while open network liquids like silica have a structural order with tetrahedral symmetry. These results show that liquids have a highly nontrivial 3D structure and that this structural information is encoded in nonstandard correlation functions.

liquids | three-dimensional structure | Lennard-Jones | silica |
computer simulations

The microscopic structure of many-particle systems is usually determined from scattering experiments that give access to the static structure factor $S(\vec{q})$, $\vec{q}$ is the wave vector, and for crystals, such measurements allow to obtain a complete knowledge of the structure of the material (1–4). This is not the case for disordered materials such as liquids, foams, and granular materials, since on the macroscopic scale, these are isotropic, and hence $S(\vec{q})$ depends only on the norm $q = |\vec{q}|$, i.e., the whole three-dimensional (3D) structural information is projected onto a single function $S(q)$. This projection entails a huge loss of structural information, which subsequently has to be recovered, at least partially, from physical arguments on the possible arrangement of the particles. Such arguments exist for the first few nearest-neighbor shells of the particles (5–17) but not for the arrangements on larger scales.

Microscopy on colloidal systems and computer simulations have shown that for hard sphere-like systems, the local structure can be surprisingly varied, in particular, if the liquid is constituted by more than one type of particle. The geometry of these locally favored structures depends on packing fraction and is rather sensitive to parameters like the composition of the system, the size of the particles, or the interaction energies (7, 11). As a consequence, it has so far not been possible to come up with a universal description of the structure on the local scale, and it is unlikely that such a universal description exists.

In view of this difficulty, it is not surprising that very little effort has been made so far to investigate the structure of disordered systems on length scales beyond the first few nearest neighbors (11, 16). A further reason for this omission is the fact that the characterization of the structure on larger scales seems to be a daunting task, since already the classification of the local structure is highly complex. Other studies have therefore focused on the possible existence of orientational order that extents to larger

distances (10, 18, 19), but no such order was found beyond the scale of a few particle diameters (19). However, whether or not disordered systems have indeed a structural order that extends beyond a few particle diameters is an important question since it is, e.g., related to the formation of the critical nucleus for crystallization or the possible growth of a static length scale that is often invoked for rationalizing the slow dynamics in glass-forming systems (3, 10, 20–24). In the present work, we use an approach to reveal that liquids do have a highly nontrivial 3D structure that is surprisingly simple at length scales beyond the first few neighbors.

In order to show the generality of our results, we will consider two systems that have a very different local structure: a binary mixture of Lennard-Jones particles (BLJM), with 80% A particles and 20% B particles (25), and silica (*Materials and Methods*). The former liquid has a close-packed local structure that is similar to the one of a hard sphere system, while the latter is a paradigm for an open network liquid with local tetrahedral symmetry (3).

We study the equilibrium properties of the BLJM in a temperature range in which the system changes from a very fluid state to a moderately viscous one, i.e., $5.0 \ge T \ge 0.40$ (25). Silica is studied at 3,000 K, a temperature at which the liquid is relatively viscous (26). To probe the 3D structure of the BLJM, we introduce a local coordinate system as follows (Fig. 1A): take any three A particles that touch each other, i.e., they form a triangle with sides that are less than the location of the first minimum in the radial distribution function $g_{AA}(r)$, i.e., $\approx 1.4$ (SI Appendix, Fig. S1A). We define the position of particle 1 as the origin, the direction from particles 1 to 2 as the z axis, and the plane containing the three particles as the z-x plane. (For the case of $SiO_2$, we use a Si atom as the central particle and two nearest

## Significance

The structure of disordered systems like liquids, gels, granular materials, etc. is considered to be isotropic, and hence very few studies exist that have investigated the three-dimensional structure of such systems. Here, we introduce a method that allows characterization of this structure. Considering the examples of a hard sphere-like liquid and of silica, an open network liquid, we show that the local three-dimensional arrangement of the particles is highly anisotropic up to distances of several particle diameters and shows a simple symmetry.

Author contributions: Z.Z. and W.K. designed research; Z.Z. and W.K. performed research; Z.Z. analyzed data; and Z.Z. and W.K. wrote the paper.

The authors declare no competing interest.

This article is a PNAS Direct Submission.

Published under the PNAS license.

Data deposition: The data discussed in the paper have been deposited in Zenodo (https://zenodo.org/record/3783469#.XsfaDxDZS.qN).

This article contains supporting information online at https://www.pnas.org/lookup/suppl/doi:10.1073/pnas.2005638117/-/DCSupplemental.

First published June 8, 2020.

<footnote n="1">¹To whom correspondence may be addressed. Email: walter.kob@umontpellier.fr.</footnote>
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Distribution of particles in three dimensions for the BLJM. (A) The definition of the local coordinate system involves three particles that are nearest neighbors to each other. (B) Radial distribution function $g_{AN}(r)$ for different temperatures (N = A + B). For the sake of clarity, the different curves have been shifted vertically by multiples of 0.3. (C) An icosahedron is the dual polyhedron of a dodecahedron and vice versa. (D–F) Density distribution $ρ(θ, φ, r)$ for different values of $r$, i.e., the distribution of the particles that are in a spherical shell of radius $r$ and thickness 0.4 around the central particle. The temperatures are 2.0, 0.8, and 0.4 for D, E, and F, respectively.</figcaption></figure>

neighbor oxygen atoms to define the three coordinate axes.) This local reference frame allows to introduce a spherical coordinate system $θ, φ, r$ and to measure the probability of finding any other particle at a given point in space, i.e., to measure a four-point correlation function. Note that this coordinate system can be defined for all triplets of neighboring particles, and these density distributions can be averaged to improve the statistics. Since this coordinate system is adapted to the configuration by the three particles, it allows to detect angular correlations that are not visible in $g(r)$ or in previously considered structural observables.

For the BLJM, Fig. 1D shows the 3D normalized distribution $ρ(θ, φ, r)$ of the particles on the sphere of radius $r$ centered at an A particle. The temperature is $T = 2.0$, i.e., above the melting point of the system, which is around $T = 1.0$ (27). We recognize that $ρ(θ, φ, r)$ has a noticeable angular dependence not only at small distances but also at intermediate ones, i.e., $r = 4.5$, which corresponds to the fifth nearest-neighbor shell in $g(r)$ (Fig. 1B). Here, we denote by $g(r)$ the partial radial distribution function $g_{AN}(r)$, where N stands for A + B (see also SI Appendix, Fig. S1). If temperature is decreased to $T = 0.8$ (Fig. 1E), the angular signal can be easily detected up to $r = 5.9$ and for $T = 0.40$ (Fig. 1F), even at $r = 8.0$, i.e., the ninth nearest-neighbor shell.

These snapshots show that this liquid has a nontrivial angular correlation that extends to distances well beyond the first few nearest-neighbor shells (see Movies S1 and S2 for a dynamical presentation of these results).

Furthermore, one notes that $ρ(θ, φ, r)$ has a highly symmetric shape: for $r = 1.65$, corresponding to the distance between the first minimum and the second nearest neighbor peak in $g(r)$, one observes a dodecahedral-like symmetry (Fig. 1C). For $r = 2.2$, the distribution has instead an icosahedral symmetry. This result can be understood by recalling that an icosahedron is the dual of a dodecahedron, and vice versa (Fig. 1C), and hence the local dips formed by particles in the first minimum will be occupied by the particles in the subsequent shell. As shown below, this “duality mechanism” works even at large distances, thus leading to a nontrivial angular correlation in which, as a function of $r$, density distributions with icosahedral symmetry alternates with distributions with dodecahedral symmetry. Fig. 1E and F show that with decreasing temperature the intensity of the signal at intermediate and large distances increases, indicating an enhanced order at low $T$.

Also for the case of the network liquid SiO$_2$, one finds a pronounced anisotropy of the density correlation function (Fig. 2).
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Distribution of particles in three dimensions for silica ($T = 3,000$ K). The thickness of the shell is 1 Å. Depending on $r$, the high-/low-density regions show a tetrahedral symmetry with two different symmetry axes.</figcaption></figure>

In contrast to the hard sphere-like liquid one finds here that the spherical shells with a pronounced orientational order all show a tetrahedral symmetry, which makes sense since the dual of a tetrahedron is again a tetrahedron. We emphasize that for geometrical reasons, at large $r$, a region with a high value of $\rho(\theta, \phi, r)$ is not a single particle but a structure that grows linearly with $r$ and hence is a whole collection of particles, i.e., for fixed $r$, the structure is given by patches with a high density of particles that alternate with patches with low density.

To analyze these findings in a quantitative manner, we use the standard procedure to decompose the signal on the sphere into spherical harmonics $Y_l^m$, $\rho(\theta, \phi, r) = \sum_{l=0}^{\infty} \sum_{m=-l}^{l} \rho_l^m(r) Y_l^m(\theta, \phi)$, where the expansion coefficients $\rho_l^m$ are given in Materials and Methods, and to consider the square root of the angular power spectrum $S_\rho(l, r) = [(2l+1)^{-1} \sum_{m=-l}^{l} |\rho_l^m(r)|^2]^{1/2}$. For the BLJM, the component with $l=6$ has the largest amplitude (SI Appendix, Fig. S2A), independent of $r$, a result that is reasonable in view of the icosahedral and dodecahedral symmetries that we find in the density distribution. For SiO$_2$, it is the component $l=3$ that has the strongest signal (SI Appendix, Fig. S2B) since this mode captures well the tetrahedral symmetry of the density field.

In Fig. 3, we show the $r$-dependence of $S_\rho(l, r)$, and one sees that for both systems, the signal decays quickly with increasing $r$. For the BLJM, Fig. 1F shows that at distance $r = 5.85$, the density distribution has a pronounced structure, while from Fig. 3A, one sees that at this $r$, the absolute value of $S_\rho(l, r)$ is small. This smallness is due to the fact that $S_\rho(l, r)$ is not only sensitive to the angular dependence of the distribution but also to the amplitude of the signal. In order to probe whether or not the density distribution has a pronounced symmetry, it is therefore useful to consider a normalized density distribution $\eta(\theta, \phi, r) = [\rho(\theta, \phi, r) - \rho_{\min}(r)]/[\rho_{\max}(r) - \rho_{\min}(r)]$, where $\rho_{\max}(r)$ and $\rho_{\min}(r)$ are the maximum and minimum of $\rho(\theta, \phi, r)$, respectively (at fixed $r$). The square root of the angular power spectrum of $\eta(\theta, \phi, r)$, $S_\eta(r)$, is included in Fig. 3A as well. We see that for the BLJM, $S_\eta$ oscillates around a constant value, which demonstrates that for this system, the density distribution has a pronounced orientational order even at large distances. For distances larger than a threshold $\xi_\eta(T)$, $S_\eta(r)$ starts to decay before it reaches at large $r$ a value that is determined by the noise of the data. (See SI Appendix for a precise definition of $\xi_\eta$ and its $T$-dependence [SI Appendix, Fig. S5].) For distances smaller than two- to three-particle diameters, there is no direct correlation between $S_\rho(r)$ and $g(r)$ since at these $r$ values, the local packing is determined also by energetic considerations (SI Appendix, Fig. S3).

Most remarkable is the observation that for the case of the BLJM, the height of the local maxima in $S_\eta(r)$ shows a periodic behavior in that a high maximum is followed by a low one. A visual inspection of $\rho(\theta, \phi, r)$ reveals that these high/low maxima correspond to distances at which the distribution has a pronounced icosahedral/dodecahedral symmetry, demonstrating that these two Platonic bodies are present not only at short distances but also at large ones, in agreement with the snapshots

in Fig. 1. One thus concludes that for hard sphere-like systems, the distribution of the particles in three dimensions is given by shells in which particles are arranged in an icosahedral pattern,

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Quantitative characterization of the structural order. (A and B) The angular power spectra and radial distribution function for the BLJM at $T = 0.4$ (A) and for silica at $T = 3,000$ K (B). The power spectra $S_\rho(6, r)$ and $S_\rho(3, r)$ (magenta curves) show an exponential-like decay as a function of the distance $r$. The power spectra for the normalized density distribution, $S_\eta(6, r)$ and $S_\eta(3, r)$ (red curves), stay large even at intermediate $r$. For the BLJM and $r \gtrsim 4.0$, the high/low maxima in $S_\eta(r)$, labeled I and D, coincide with the minima/maxima (labeled M) in $|g(r) - 1|$ (blue line). This up-down behavior is related to the alternating icosahedral/dodecahedral symmetry in the distribution of the particles when $r$ is increased. For SiO$_2$, the arrows indicate the distances at which $g_{isi}(r) = 1$. For better visibility, $S_\eta(r)$ and $|g(r) - 1|$ have been shifted vertically.</figcaption></figure>
followed by a shell in which there is a dodecahedral pattern. For distances larger than $r \approx 4$, one finds that the radial positions of these two geometrical arrangements match perfectly the locations of the minima/maxima in $g(r)$ (Fig. 3A). This observation can be rationalized by the fact that a dodecahedron has 20 vertices [i.e., regions in which $\rho(\theta, \phi, r)$ has high values] and an icosahedron only 12, thus making that the former structure corresponds to the maxima of $g(r)$ and the latter to the minima. Similar results are found for $T = 2.0$ (SI Appendix, Fig. S4), demonstrating that these findings hold also for liquids that are not supercooled.

In contrast to the BLJM, we find that for silica (Fig. 3B), the locations of the maxima in $S_ρ(3, r)$ do not correspond to the ones in $|g_{SiSi}(r) - 1|$ but are instead close to distances at which $g(r) = g_{SiSi}(r) = 1$ (as indicated by the arrows in the graph), i.e., corresponds to a distance at which one expects no structural order. (See SI Appendix, Fig. S1D for the Si–O partial radial correlation function.) This shows that for liquids that have an open network structure, the distances at which one finds the highest orientational symmetry is not associated with a dense packing of particles, in contrast to hard sphere-like systems. Finally, we note that for both systems, the decay of $S_ρ(r)$ matches very well the one of $g(r)$. This indicates that the two functions are closely related to each other, i.e., the loss in the symmetry of the density field in three dimensions leads to the decay in the structure as measured by $g(r)$, a result that is reasonable since the angular integral of $\rho(\theta, \phi, r)$, if not normalized, is proportional to $g(r)$.

Fig. 3B shows that also for silica, $S_η(r)$ is high for small and intermediate distances, but even in this range, it decreases slowly, indicating that for this network liquid, the orientational symmetry is gradually lost with increasing $r$. This result might be related to the fact that the coordination number of the silica network is

lower than the one of the densely packed hard sphere-like liquid; hence, the former structure is more flexible, and therefore it is more difficult to propagate the orientational order in space to large distances.

Since we have found that the distribution of the particles around a central particle is anisotropic, it is of interest to consider also the radial distribution functions in which one probes the correlations in a specific direction with respect to the local coordinate system shown in Fig. 1A. This type of information can be obtained for colloidal systems from confocal microscopy experiments and, more indirectly, from scattering experiments (8). *Insets* in Fig. 4A and C show the directions we considered for the two type of liquids: for the BLJM, the directions corresponding to the vertices of the icosahedra/dodecahedra and the directions given by the midpoints between these two type of vertices; for silica, the directions of the vertices of the tetrahedra, the points given by the midpoints of the faces of the tetrahedra, and the directions given by the midpoints between the two former directions.

In Fig. 4A, we show for the BLJM, the radial distribution functions for these special directions, and one recognizes that the amplitude of the signal depends indeed strongly on the direction considered. For the directions of the icosahedra and of the dodecahedra, $g_I(r)$ and $g_D(r)$, respectively, we find for intermediate and large distances that $g_D(r)$ oscillates in phase with $g(r)$, whereas $g_I(r)$ has oscillations that are in antiphase. These observations are coherent with the aforementioned argument that the number of vertices in the dodecahedra exceed the ones for the icosahedra.

Furthermore, Fig. 4A shows that the amplitudes of the oscillations in $g_I(r)$ and $g_D(r)$ are significantly larger than the ones found in $g(r)$, a result that is reasonable since the latter

<figure><img src="image_5.png" /><figcaption>**Fig. 4.** Anisotropic radial distribution functions. (A) BLJM. (A, Inset) I and D are the directions defined by the vertices of the icosahedra/dodecahedra. N is the midpoint between these two directions. (A) Radial distribution function as measured in the directions I, D, and N. (B) $|g(r) - 1|$ on logarithmic scale. For the sake of clarity, only the maxima in the curves are shown. (C) Silica. (C, Inset) T1 and T2 are the directions defined by the two interlocked tetrahedra. N is the midpoint between these two directions. (C) Radial distribution function as measured in the direction T1, T2, and N. (D) $|g(r) - 1|$ on logarithmic scale. For the sake of clarity, only the maxima in the curves are shown.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Three-dimensional representation of the density field. The shown layers correspond to distances at which $S_ρ(r)$ has a local maximum. Only regions with high density (covering 35% area of the sphere) are depicted. (A) BLJM at $T = 0.4$. The bluish/reddish colors correspond to the locations of the high/low maxima in $S_η(r)$ and thus to shells with icosahedral/dodecahedral symmetry. (B) Silica at $T = 3,000$ K. The high-density regions form interlocked zones with a tetrahedral symmetry.</figcaption></figure>

function is a weighted average of the two former ones and hence will be affected by cancellation effects. The distribution function in the direction that corresponds to the midpoint of the line connecting two neighboring vertices of an icosahedron and a dodecahedron, $g_N(r)$, shows significantly smaller oscillations than $g(r)$, a result that is expected since one probes the structure in a direction that does not pass close to the locations that correspond to the vertices of the icosahedra/dodecahedra. Fig. 4B shows these distribution functions on a logarithmic scale. (For the sake of clarity, only the maxima and minima of the functions are shown.) One notices that the slope of the curves for $g(r)$, $g_I(r)$, and $g_D(r)$, i.e., the length scale over which the correlation decays, is basically independent of the function considered, demonstrating that they are indeed closely related to each other. In contrast to this, the data from $g_N(r)$ decay faster, showing that in this direction, the correlation length is smaller.

For the case of silica, the connection between the extrema in $g(r)$ and the ones obtained from the radial distribution functions in the special directions T1 and T2 (Fig. 4C, Inset) is not straightforward. One finds that the peaks in $g_{T1}(r)$ ($g_{T2}$) are where $g(r)$ rises (decreases) quickly (Fig. 4C). In fact, the extrema of $g_{T1}(r)$ are very close to the distances at which $g(r)$ becomes 1.0, i.e., an $r$ at which the Si density corresponds to the one expected for an ideal gas. The reason for this is presently not known, and thus it will be interesting to determine whether this is a general feature for liquids that have an open network structure.

The radial distribution function for the “neutral” direction $N$ has a signal that is in phase with $g(r)$, and its amplitude is smaller than the one of $g(r)$. The latter result is expected since one measures the density field in a direction in which the fluctuations between the interlocked tetrahedra basically cancel each other. Fig. 4D shows on logarithmic scale the maxima of $|g(r) - 1|$ that correspond to the $g(r)$ curves in Fig. 4C. One recognizes that all of them decay in the same exponential manner with a slope that is independent of the direction.

To give a comprehensive view of the particle arrangement in three dimensions, we present in Fig. 5 the density distribution of the two systems. The colored regions correspond to the zones in which the particle density is high, and, by construction, they cover 35% of the sphere. For the BLJM at intermediate and large distances, one recognizes clearly the presence of high density zones with icosahedral symmetry (bluish color) interlocked with zones with dodecahedral symmetry (reddish color). The directions in which the blue and red regions touch each other correspond to the neutral direction $N$ defined above and in which the particle correlation is weak. For silica, one finds instead interlocked tetrahedra at all distances (Fig. 5B). Again, the neutral direction corresponds to the one in which the blue and yellow regions touch.

In conclusion, we have demonstrated that liquids have nontrivial structural symmetries at beyond short range that have gone unnoticed so far. This result has been obtained by using a method that takes into account the 3D angular dependence of the structure and which can be readily applied to any system for which the particle coordinates are accessible, such as colloidal and granular systems, or materials in which some of the particles have been marked by fluorescence techniques (9, 28–31). Since we find that the nature of the orientational order in three dimensions depends on the system considered, the method allows to make a more precise classification of the structure of liquids, an aspect that should trigger the improvement of experimental techniques that probe this structural order.

## Materials and Methods

**System and Simulations.** The BLJM we study is a 80:20 mixture of Lennard-Jones particles (types A and B) with interactions given by $V_{α,β}(r) = 4ε_{α,β}[(σ_{α,β}/r)^{12} - (σ_{α,β}/r)^6]$, where $α, β \in \{A, B\}$, $σ_{AA} = 1.0$, $ε_{AA} = 1.0$, $σ_{AB} = 0.8$, $ε_{AB} = 1.5$, $σ_{BB} = 0.88$, and $ε_{BB} = 0.5$ (25). Here, we use $σ_{AA}$ and $ε_{AA}$ as the units of length and energy, respectively. We set the mass of all particles equal to $m = 1.0$ and the Boltzmann constant is $k_B = 1.0$. We simulate a total of $10^5$ particles at constant volume (box size, 43.68) and temperature. At the lowest temperature, $T = 0.40$, the run was $1.4 \cdot 10^8$ time steps (step size is 0.005) for equilibration and the same length for production, time spans that are sufficiently large to completely equilibrate the system, i.e., the mean squared displacement (MSD) of the particles was larger than 1.0. For higher $T$ values, the MSD was significantly larger. For the analysis of the data, we used 8 and 20 configurations for $S_ρ$ and $g(r)$, respectively. Previous studies have shown that this system starts to show a tendency to crystallization at temperatures around 0.4 if the system is simulated for several $α$-relaxation times (32). However, our simulations lasted only a few $α$-relaxation times (the MSD has reached 1 to 2), which is long enough for equilibrating the liquid but not long enough for allowing the system to crystallize (as indicated by the static structure factor; SI Appendix, Fig. S1). Also, all of our results show a completely smooth dependence on temperature, and hence it is unlikely that they are affected by the presence of crystalline order.

For the simulation of silica, we use a recently optimized interaction potential proposed by Sundararaman et al. and which has been shown to be able to describe reliably the properties of real silica (33). Although this potential is based only on pair interactions, it has been found to be able to describe better the structural and mechanical properties than other potentials for silica, including potentials with three-body interactions (33, 34). A cubic simulation box containing 120,000 atoms was used, which corresponds at room temperature and zero pressure to a box size of about 120 Å. The simulation was carried out in the isothermal-isobaric ensemble at 3,000 K for $10^6$ time steps (step size is 1.6 fs). This time span is sufficiently long to equilibrate the liquid. (The analysis of the MSD indicates that the liquid has become diffusive and the mean displacement of Si is larger than 18 Å). After equilibration, we collected 8 configurations spaced by $10^5$ time steps for the subsequent structural analysis. Note that the melting temperature of silica is around 2,000 K; thus, our simulations are far above the melting point, and the system is in its stable liquid phase. For both systems, we used the Large-scale Atomic/Molecular Massively Parallel Simulator software (35) to carry out the simulations.

**Angular Power Spectrum.** The coefficient $ρ_l^m$ for the expansion of the density distribution into spherical harmonics is given by

$$ρ_l^m = \int_0^{2\pi} d\phi \int_0^{\pi} d\theta \sin\theta \rho(\theta, \phi, r) Y_l^{m*}(\theta, \phi), \quad [1]$$

where $Y_l^{m*}$ is the complex conjugate of the spherical harmonic function of degree $l$ and order $m$. In practice, this integration was done for the BLJM by sampling the integrand over up to $2 \cdot 10^9$ points for each shell of width 0.4. The corresponding numbers for SiO$_2$ are $10^8$ points and 1.0 Å.

**Data Availability.** The data discussed in the paper have been deposited in Zenodo (https://zenodo.org/record/3783469#.XsfaDxDZS.qN).

**ACKNOWLEDGMENTS.** We thank D. Coslovich, G. Monaco, M. Ozawa, and K. Schweizer for discussions. Part of this work was supported by China Scholarship Council Grant 201606050112 and Agence Nationale de la Recherche (ANR) Grant ANR-15-CE30-0003-02.
1. N. W. Ashcroft, N. D. Mermin, *Solid State Physics* (Holt, Rinehart and Winston, New York, NY, 1976).

2. J. P. Hansen, I. R. McDonald, *Theory of Simple Liquids* (Elsevier, Amsterdam, The Netherlands, 1990).

3. K. Binder, W. Kob, *Glassy Materials and Disordered Solids: An Introduction to Their Statistical Mechanics* (World Scientific, Singapore, 2011).

4. P. S. Salmon, Decay of the pair correlations and small-angle scattering for binary liquids and glasses. *J. Phys. Condens. Matter* **18**, 11443–11469 (2006).

5. H. Jónsson, H. C. Andersen, Icosahedral ordering in the Lennard-Jones liquid and glass. *Phys. Rev. Lett.* **60**, 2295–2298 (1988).

6. D. B. Miracle, A structural model for metallic glasses. *Nat. Mater.* **3**, 697–702 (2004).

7. Y. Cheng, E. Ma, Atomic-level structure and structure–property relationship in metallic glasses. *Prog. Mater. Sci.* **56**, 379–473 (2011).

8. P. Wochner et al., X-ray cross correlation analysis uncovers hidden local symmetries in disordered matter. *Proc. Natl. Acad. Sci. U.S.A.* **106**, 11511–11514 (2009).

9. C. Xia et al., Origin of noncubic scaling law in disordered granular packing. *Phys. Rev. Lett.* **118**, 238002 (2017).

10. D. R. Nelson, F. Spaepen, "Polytetrahedral order in condensed matter" in *Solid State Physics* (Elsevier, 1989), Vol. 42, pp. 1–90.

11. C. P. Royall, S. R. Williams, The role of local structure in dynamical arrest. *Phys. Rep.* **560**, 1–75 (2015).

12. A. Malins, J. Eggers, C. P. Royall, S. R. Williams, H. Tanaka, Identification of long-lived clusters and their link to slow dynamics in a model glass former. *J. Chem. Phys.* **138**, 12A535 (2013).

13. A. J. Dunleavy, K. Wiesner, R. Yamamoto, C. P. Royall, Mutual information reveals multiple structural relaxation mechanisms in a model glass former. *Nat. Commun.* **6**, 6089 (2015).

14. D. Coslovich, G. Pastore, Understanding fragility in supercooled Lennard-Jones mixtures. I. Locally preferred structures. *J. Chem. Phys.* **127**, 124504 (2007).

15. M. Leocmach, H. Tanaka, Roles of icosahedral and crystal-like order in the hard spheres glass transition. *Nat. Commun.* **3**, 974 (2012).

16. X. Fang, C. Wang, Y. Yao, Z. Ding, K. Ho, Atomistic cluster alignment method for local order mining in liquids and glasses. *Phys. Rev. B* **82**, 184204 (2010).

17. X. Fang et al., Spatially resolved distribution function and the medium-range order in metallic liquid and glass. *Sci. Rep.* **1**, 194 (2011).

18. P. J. Steinhardt, D. R. Nelson, M. Ronchetti, Bond-orientational order in liquids and glasses. *Phys. Rev. B* **28**, 784–805 (1983).

19. R. M. Ernst, S. R. Nagel, G. S. Grest, Search for a correlation length in a simulation of the glass transition. *Phys. Rev. B* **43**, 8070–8080 (1991).

20. G. Adam, J. H. Gibbs, On the temperature dependence of cooperative relaxation properties in glass-forming liquids. *J. Chem. Phys.* **43**, 139–146 (1965).

21. K. Kelton, A. L. Greer, *Nucleation in Condensed Matter: Applications in Materials and Biology* (Elsevier, Amsterdam, The Netherlands, 2010).

22. C. P. Royall, W. Kob, Locally favoured structures and dynamic length scales in a simple glass-former. *J. Stat. Mech. Theor. Exp.* **2017**, 024001 (2017).

23. D. Chandler, J. P. Garrahan, Dynamics on the way to forming glass: Bubbles in space-time. *Annu. Rev. Phys. Chem.* **61**, 191–217 (2010).

24. X. Xia, P. G. Wolynes, Fragilities of liquids predicted from the random first order transition theory of glasses. *Proc. Natl. Acad. Sci. U.S.A.* **97**, 2990–2994 (2000).

25. W. Kob, H. C. Andersen, Testing mode-coupling theory for a supercooled binary Lennard-Jones mixture: The van Hove correlation function. *Phys. Rev. E* **51**, 4626–4641 (1995).

26. J. Horbach, W. Kob, Static and dynamic properties of a viscous silica melt. *Phys. Rev. B* **60**, 3169–3181 (1999).

27. U. R. Pedersen, T. B. Schrøder, J. C. Dyre, Phase diagram of kob-andersen-type binary Lennard-Jones mixtures. *Phys. Rev. Lett.* **120**, 165501 (2018).

28. W. K. Kegel, A. van Blaaderen, Direct observation of dynamical heterogeneities in colloidal hard-sphere suspensions. *Science* **287**, 290–293 (2000).

29. E. R. Weeks, J. C. Crocker, A. C. Levitt, A. Schofield, D. A. Weitz, Three-dimensional direct imaging of structural relaxation near the colloidal glass transition. *Science* **287**, 627–631 (2000).

30. J. F. Sherson et al., Single-atom-resolved fluorescence imaging of an atomic mott insulator. *Nature* **467**, 68–72 (2010).

31. B. Kou et al., Granular materials flow like complex fluids. *Nature* **551**, 360–363 (2017).

32. T. S. Ingebrigtsen, J. C. Dyre, T. B. Schrøder, C. P. Royall, Crystallization instability in glass-forming mixtures. *Phys. Rev. X* **9**, 031016 (2019).

33. S. Sundararaman, L. Huang, S. Ispas, W. Kob, New optimization scheme to obtain interaction potentials for oxide glasses. *J. Chem. Phys.* **148**, 194504 (2018).

34. Z. Zhang, S. Ispas, W. Kob, The critical role of the interaction potential and simulation protocol for the structural and mechanical properties of sodosilicate glasses. *J. Non-Cryst. Solids* **532**, 119895 (2020).

35. S. Plimpton, Fast parallel algorithms for short-range molecular dynamics. *J. Comput. Phys.* **117**, 1–19 (1995).