# Configuration correlation governs slow dynamics of supercooled metallic liquids

^{a}Institute of Physics, Chinese Academy of Sciences, 100190 Beijing, China;^{b}University of Chinese Academy of Sciences, 100049 Beijing, China;^{c}Department of Mechanical and Biomedical Engineering, City University of Hong Kong, Kowloon, Hong Kong, China;^{d}Division of Physics and Applied Physics, School of Physical and Mathematical Sciences, Nanyang Technological University, 637371 Singapore, Singapore;^{e}Beijing Computational Science Research Center, 100193 Beijing, China

See allHide authors and affiliations

Edited by Pablo G. Debenedetti, Princeton University, Princeton, NJ, and approved May 9, 2018 (received for review February 7, 2018)

## Significance

The search for a structural origin governing the dynamical slowing down of a supercooled liquid toward glass transition is an active area of the community of amorphous materials. In the past decade, the locally preferred geometrical orderings, that is, those local polyhedral packing clusters extracted from instantaneous atomic configurations, such as icosahedron, have been suggested as the structural origin of slow dynamics in metallic glass-forming liquids. Here, we demonstrate that it is the intrinsic correlation between configurations that captures the structural origin governing slow dynamics. A correlation length extracted from these configurations' correlation plays a more important role than various dynamic correlation lengths in determining the drastic dynamical slowdown of supercooled metallic liquids.

## Abstract

The origin of dramatic slowing down of dynamics in metallic glass-forming liquids toward their glass transition temperatures is a fundamental but unresolved issue. Through extensive molecular dynamics simulations, here we show that, contrary to the previous beliefs, it is not local geometrical orderings extracted from instantaneous configurations but the intrinsic correlation between configurations that captures the structural origin governing slow dynamics. More significantly, it is demonstrated by scaling analyses that it is the correlation length extracted from configuration correlation rather than dynamic correlation lengths that is the key to determine the drastic slowdown of supercooled metallic liquids. The key role of the configuration correlation established here sheds important light on the structural origin of the mysterious glass transition and provides an essential piece of the puzzle for the development of a universal theoretical understanding of glass transition in glasses.

The underlying mechanism of how viscosities (or structural relaxation times) of glass-forming liquids surges by many orders of magnitude on approaching glass-transition points remains one of the most controversial issues in fundamental sciences (1). A large variety of theoretical frameworks involving growing spatial correlations, either dynamic or static, were proposed to explain the origin of the spectacular slowdown (ref. 2 and therein). In their pioneering theory, Adam and Gibbs (3) proposed that the rearrangements of atoms are collective in localized domains during cooling. The size and lifetime of the cooperatively rearranging regions increase with reducing temperature, thereby leading to vanishing configurational entropy and drastic slowing down of dynamics. This picture of heterogeneous dynamics was later supported by advanced experiments and large-scale computer simulations (4, 5), which shows a wide distribution of local dynamics in supercooled liquids with some domains moving significantly faster or slower than the average. The dynamic correlation length quantifying the extent of spatially correlated motions increases remarkably toward the glass transition. Recently, experiments and modeling related to the spatial heterogeneities in metallic glasses (MGs) also support this conception and imply the potential correlation between dynamics and structure (6⇓⇓–9). However, the puzzling fact is that the remarkable dynamic slowdown is not accompanied by any obvious structure changes based on two-point correlators in traditional diffraction and scattering experiments. Thus, the structural origin of the glass transition is still mysterious.

Inspired by Frank’s proposal (10) in 1952 that densely packed icosahedron showing incompatible symmetry with crystallographic symmetries would stabilize a liquid with constituents of equal size, many research works have been pursuing the role of these locally preferred structures in determining glassy properties (11⇓⇓⇓⇓–16). The dynamics of the centers of icosahedral clusters was found to be slow, because of strong constraints imposed by their neighbors. Thus, geometrical frustration induced by icosahedral clusters has long been suggested as the structural origin of slow dynamics in metallic glass-forming liquids (ref. 12 and therein). However, the concentration of icosahedral clusters in MGs is known to be quite low and highly composition dependent (12). What exact role of geometrical orderings determined by local atomic packings plays in MG formation remains arguable. It is urgent to unravel whether there is an intrinsic order-agnostic parameter corresponding to configurations which governs the slow dynamics of supercooled metallic liquids. Furthermore, all previous studies on slowdown of metallic liquids include the influences of changes in temperature, which itself will induce denser packings during cooling and thus cannot provide a direct link between possible structural orderings and dynamical slowdown. Therefore, the mechanism of vitrification in MGs has not yet been well understood. Moreover, although the prominent phenomenon of dynamic heterogeneity has aroused intensive interests and become a central topic in glass physics (17), it is still elusive whether dynamic heterogeneity is the consequence or the primary origin of the dynamical slowdown (18) in metallic glass-forming liquids.

In this work, we studied the slowing down of dynamics in realistic model MGs by molecular dynamics simulations based on embedded-atom method (EAM) potentials (*Materials and Methods*). By including confinement effects in equilibrated supercooled metallic liquids, we observed strong decoupling between local geometrical orderings, such as icosahedron, and dynamical slowdown, which challenges the common belief that slow dynamics is governed by icosahedral clusters in MGs. More importantly, we found that the dramatic slowdown is determined by the correlation length evaluated from configuration correlation rather than by various dynamic correlation lengths stemming from multipoint correlators. The order-agnostic correlation length provides valuable clues to unravel the mystery of glass transition in MGs.

## Results and Discussion

To reduce the interatomic potential and composition influences, two prototypical compositions Cu_{50}Zr_{50} and Cu_{46}Zr_{46}Al_{8} described by different potentials were chosen in this study (19, 20). We show the results of the former in the main text and the latter in *SI Appendix*. The dynamics of these supercooled liquids in the bulk state is shown in *SI Appendix*, Fig. S1. It is obvious that the dynamics slows down dramatically with decreasing temperature with an extended plateau in the self-part of intermediate scattering function

### Dynamics in Confinement.

Unlike all previous studies focused on the temperature drop triggered slowdown, we employed atomic pinning methods (22, 23) to eliminate the effects of changes in temperature and extract the direct link between structure and dynamics (see *Materials and Methods*). Two pinning strategies, “sandwich” pinning and random pinning, were utilized to overcome possible effects of pinning geometries. Similar results were obtained from both geometries so we only show results of the former in the main text.

Fig. 1*A* shows a sandwich-pinning configuration. Within this geometry, two parallel symmetric rough walls were constructed by freezing the atoms in an equilibrated liquid configuration at a designed temperature *T* (23, 24) and allowing the unpinned atoms to evolve at *T*. After imposing the walls, no phase separation and layering happen, which is confirmed by the invariant number density and local composition at different distances to the walls compared with those of the bulk (*SI Appendix*, Fig. S2). We also did not find any influence of the confinement on partial pair-correlation functions. In view of the above results, we were ready to study the confinement effects on the geometrical structures and relaxation dynamics of the supercooled metallic liquids (24) for which such novel simulations have never been carried out before.

To characterize the local dynamics as a function of *z*, the distance of atoms to the walls, we generalized the conventional form of *j* to the closer wall at time *t* = 0 (see *SI Appendix*, Fig. S3 for more details). Here we analyzed the dynamics without distinguishing motions at different directions because the influence of dynamic anisotropy induced by wall construction is negligible and similar results were obtained by only considering the motion parallel to the walls (see *SI Appendix* for details). The *z*-dependent relaxation times *B* depicts _{50}Zr_{50} at 800 K. Obviously, the *z* descends are quite similar to the case during cooling. The classical two-step relaxation appears with a prolonged plateau at smaller *z*. It means that the dynamics near the walls is strongly reduced while it recovers to the bulk dynamics when far away from them. Meanwhile, the slowing down of dynamics near the walls is accompanied by increasing dynamic heterogeneity reflected by the non-Gaussian parameter *z* (Fig. 1*C*). *α*-relaxation part of *SI Appendix*, Fig. S3). These observations strongly demonstrate that dynamical slowdown is intrinsically accompanied by dynamic heterogeneity, even though there is no temperature change. In other words, glassy dynamics can be realized by introducing confinements to exclude the effects of changes in temperature in supercooled metallic liquids. In what follows, we will provide a detailed analysis to understand whether dynamic heterogeneity is the origin of slow dynamics or not. To quantify the slowdown of dynamics with confinement, we extracted the corresponding *D* shows *z* at various temperatures in comparison with *T* higher than the glass-transition temperature of the bulk. Since this dynamical slowdown at a constant temperature is different from the traditional hyperquenching, we can study the pure structure evolution of supercooled metallic liquids correlated with dynamical slowdown in the absence of changes in temperature.

### Decoupling Between Dynamics and Local Structural Orderings.

We used Voronoi tessellation to characterize the local structure (*Materials and Methods*). Given that full icosahedron *f*_{ico}) as function of *z* as an example in Fig. 2*A*. The numbers of the Voronoi index *f*_{ico} increases during cooling, it is unambiguous that its fraction is low and it is almost *z*-independent albeit with some fluctuations. Similar invariant trends are also observed in other atomic clusters (*SI Appendix*, Fig. S4). Similar results were also obtained in random pinning and in Cu_{46}Zr_{46}Al_{8} (*SI Appendix*). Locally preferred clusters have commonly been treated as the key factor controlling dynamical slowdown in MGs (12). If that were the truth, any dynamics change in the glass-forming liquid would be in accord with the change in these structural orderings, including in the confined systems. Fig. 2*B* shows *z* and *f*_{ico}. Obviously, the relaxation times can change hugely without changing *f*_{ico}. The sharp contrast between *D* and 2*A* (and *SI Appendix*, Fig. S4) clearly demonstrates decoupling of slow dynamics and local geometrical orderings in the absence of temperature changes. This delivers a strong message that local structural orderings, like icosahedra, shall not play a determining role in the slowing down of dynamics, at least in the pinning-induced slow dynamics. Our findings contrast the common belief that local geometrical orderings extracted from instantaneous configurations be a universal structural origin of slow dynamics in metallic glass-forming liquids.

### Coupling Between Dynamics and Configuration Correlation.

Now we come to the long-debated question of what order parameter originating from configurations could physically determine the slow dynamics of supercooled metallic liquids in confinement. Instead of just considering one instantaneous configuration, we considered the correlation between two configurations. In the sandwich-pinning geometry (23), we defined an overlap function *l* = 0.68 Å so that there would be no more than one atom in a single box. A binary digit *i*th box is occupied by an atom at *t*, and =0 otherwise. The overlap profile quantifying the similarity between two configurations separated by *t* in the *z* direction is measured by *z* from the walls. In this calculation, we considered two instantaneous configurations at a temperature without finding their inherent structures. Fig. 3*A* shows the decaying behavior of *z*, the curves collapse with the bulk sample indicating negligible effects of the confinement, consistent with Fig. 1 *B* and *D*. At long time limit, these profiles decay to a plateau *z*, *t* increases due to the confinement. To quantify this effect, we fitted the final decay of *A*, τ, β, and *z*. We can then obtain the change tendency of *z* and *T* (*SI Appendix*, Fig. S5). To explore the role of the configuration correlation in dynamics, it is necessary to establish the relation between *SI Appendix*, Fig. S5), indicating that the available configuration states are strongly constrained by the pinning walls. Surprisingly, when we plot *B*, curves of all temperatures collapse into a master curve:*B*). The

Furthermore, as shown in Fig. 3*C*, the relation between *z* can be well fitted by an exponential form *SI Appendix*, Fig. S6). The correlation length *B* and *C* can lead to the universal scaling relation of *D*. The collapse is better at lower temperatures in the FEL-dominated region, consistent with previous studies (27), and was also found in Cu_{46}Zr_{46}Al_{8} (*SI Appendix*). It indicates that the decoupling relation between *z* (Fig. 3*D*, *Inset*) is controlled by *Materials and Methods* and *SI Appendix*, Fig. S7). In contrast, no scaling collapse can be achieved between *SI Appendix*, Fig. S8). These findings deliver a strong message that it is the configuration correlation rather than dynamic heterogeneities measured by the dynamic correlation lengths that controls the slow dynamics in the confined metallic liquids.

### Critical Role of the Correlation Length in Governing Slow Dynamics.

Next, we move further to unravel the role of *SI Appendix*). We calculated _{50}Zr_{50} at a series of system sizes *N* (*Materials and Methods*). As shown in Fig. 4*A* (*Inset*), there is a clear *N* dependence of *Materials and Methods*) is quite small compared with *SI Appendix*, Fig. S9). This further indicates that not only is its fraction low as shown before, the correlation between

Since dynamic heterogeneity has long been thought to play an important role in explaining slow dynamics, we also calculated four-point dynamic correlation length *Materials and Methods*) in addition to *B*. Evidently, dynamic length scales grow faster than *SI Appendix*, Fig. S10) (18). These findings demonstrate that the slow dynamics in MGs is not governed by dynamic heterogeneities.

The finite-size scaling analysis of structural relaxations versus structural and dynamic correlation lengths provides the compelling evidence that dynamical slowdown in the metallic glass-forming liquids is determined by the structural correlation *SI Appendix*) (33, 34).

## Conclusion

In stark contrast to the previous knowledge about MGs, our work suggests that local geometrical orderings such as icosahedra are not the structural origin of slow dynamics, at least in pinned MGs. Instead, we found that the key factor governing slow dynamics of metallic glass-forming liquids, whether in pinned or in cooled, is the length scale defined by the intrinsic configuration correlation. We also demonstrate that the growing dynamic heterogeneity, although intrinsically accompanying dramatic slowing down of the dynamics during glass transition, is not the primary origin but the consequence of slow structural relaxations in metallic glass-forming liquids. The determining role of the configuration correlation gives direct evidence that glass transition is not a purely dynamical process but has a certain thermodynamic origin. The hidden length scale based on the configuration correlation sheds light on the structural origin of the mysterious glass transition. Our findings also offer one essential piece of the puzzle for the RFOT theory and suggest that the order-agnostic length scale may be universal among a wide range of glass-formers. One important step in the future is to transform the physical definition of the length scale to experimentally measurable observations (18, 35⇓⇓–38). In our view, this leads to another unresolved issue, i.e., what geometric structural feature in a supercooled liquid, if there is any, could be responsible for the increase of the correlation length as temperature decreases. This is definitely a good problem to be attacked in the follow-up studies. Finally, it is worth noting that our findings here are based on numerical simulations of supercooled metallic liquids and thus call for further experimental verifications.

## Materials and Methods

### Molecular Dynamics Simulations.

In this study, we performed all of the molecular dynamics simulations using the open source code LAMMPS (39) based on EAM potentials. These empirical potentials describe the relatively common atomic interactions beyond the pure pair interactions, such as Lennard-Jones, Kob–Andersen, and the inverse power-law potentials, by introducing many-body interactions. The many-body nature of the EAM potentials is a result of the embedding energy term. In our simulations, periodic boundary conditions were applied in three directions. For each system, the initial configuration containing 16,000 atoms was first equilibrated at 2,000 K for at least 1.5 ns and then it was equilibrated at desired temperatures in the NPT ensemble (constant number, constant pressure, and constant temperature), during which the cell size was adjusted to give zero pressure. After equilibration at each temperature of interest, the ensemble was switched to the NVT ensemble (constant number, constant volume and constant temperature) for further equilibration and production runs. The time step was set to 0.002 ps and the temperature was controlled using the Nosé–Hoover thermostat.

### Pinning Atoms in Sandwich Geometry.

We studied the influence of amorphous walls on the structure and dynamics of metallic glass-forming liquids using samples containing 16,000 atoms (∼65 Å × 65 Å × 65 Å). Two symmetric rough walls with thickness of 5 Å were created in an equilibrium liquid configuration by freezing the atoms at the boundaries along the *z* axis. Then we let the free atoms relax at the initial temperature and the periodic boundary conditions along the *z* axis were removed. Therefore, the atoms near one wall will not interact with those near the other wall. To prevent atoms from penetrating into the frozen walls, we also introduced an infinitely hard wall at the two interfaces between walls and unpinned atoms. In this geometry, we divided the confined sample into slices with a thickness of 1.5 Å. The results were averaged over the two slices that have the same distance to one of the walls, and over different realizations of the walls. Since the dynamics near the center always recovers to the bulk dynamics in our studies, the system size is sufficiently large and finite-size effects are absent. We considered eight temperatures in the supercooled liquids in Cu_{50}Zr_{50}: *T* = 780; 800; 850; 900; 950; 1,000; 1,100; and 1,200 K. For the Cu_{46}Zr_{46}Al_{8}, the considered temperatures were 850; 870; 880; 900; 950; 1,000; and 1,200 K. The Vogel–Fulcher–Tammann temperatures of the two systems are 612 and 726 K, respectively (*SI Appendix*, Fig. S1*B*).

### Voronoi Tessellation.

A well-accepted method to characterize the structure of MGs is Voronoi tessellation. It divides space into close-packed polyhedra around atoms by constructing bisecting planes along the lines joining the central atom and all its neighbors. The Voronoi index <*n*_{3}, *n*_{4}, *n*_{5}, *n*_{6}> is used to characterize the geometry feature of atomic clusters, where *n*_{i} (*i* = 3, 4, 5, 6) denotes the number of *i*-edged faces of a Voronoi polyhedron. Voronoi tessellation was performed using the voro++ package (math.lbl.gov/voro++/) (40). For the sandwich-pinning and random-pinning geometries, the configurations for Voronoi tessellation were collected after simulations at least 10 times of the corresponding structural relaxation times. In random pinning, both pinned and unpinned atoms were taken into consideration when conducting Voronoi analysis, to avoid artificial holes created by pinned atoms. However, when calculating the fractions, only unpinned atoms were considered.

### Dynamic Correlation Lengths of Nonlinear Reponses.

According to the profiles of *z*, showing the strong effects of the walls on atomic rearrangements. Such influence is more remarkable at lower temperatures, indicating growing dynamic correlations. To extract the dynamic correlation lengths *SI Appendix*, Fig. S7). By fitting the scaled timescales at large *z* to the exponential forms,

### Four-Point Correlation Lengths.

To accurately measure the four-point dynamic correlation length *n* = 462,150 to calculate *Molecular Dynamics Simulations*. To save computer time, the timestep was chosen as 0.005 ps. To characterize dynamic heterogeneity, we first defined an overlap function *d*_{Zr}; *d*_{Zr} is the atomic diameter of Zr), and =0 otherwise. The fluctuation of *SI Appendix*, Fig. S11). This demonstrates the accuracy of the estimation of *SI Appendix*, Fig. S9).

### Finite-Size Scaling.

We used the model system Cu_{50}Zr_{50} for finite-size scaling analysis. The system size changes from *n* = 100 to *n* = 1,600. The simulation method is similar to the description in *Molecular Dynamics Simulations*. We also included the large system with *n* = 462,150 for analysis. Ten independent simulations were carried out for ensemble average. To obtain *N* dependence of *n* = 462,150 as

## Acknowledgments

We are grateful to Dr. Lijin Wang for helpful discussion. We also acknowledge the computational support from the Beijing Computational Science Research Center. This work is supported by National Natural Science Foundation (NSF) of China Grants 11790291 and 51271195, Ministry of Science and Technology of the People’s Republic of China 973 Program Grant 2015CB856800, Key Research Program of Frontier Sciences Grant QYZDY-SSW-JSC017, and Strategic Priority Research of the Chinese Academy of Sciences Grant XDPB06. P.-F.G. is also supported by NSF of China Grants 51571011 and U1530401. Y.Y. acknowledges the financial support provided by the Research Grant Council–NSF of China joint fund with the Grant of N_CityU116/14.

## Footnotes

- ↵
^{1}To whom correspondence may be addressed. Email: ychu0213{at}gmail.com or pguan{at}csrc.ac.cn.

Author contributions: Y.-C.H., P.-F.G., and W.-H.W. designed research; Y.-C.H. performed research; Y.-C.H., Y.-W.L., Y.Y., P.-F.G., H.-Y.B., and W.-H.W. analyzed data; and Y.-C.H., Y.Y., P.-F.G., and W.-H.W. wrote the paper.

The authors declare no conflict of interest.

This article is a PNAS Direct Submission.

Data deposition: Data relevant to this paper are available at https://github.com/yuanchaohu/PNAS2018.

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

Published under the PNAS license.

## References

- ↵
- ↵
- ↵
- ↵
- ↵
- Berthier L, et al.

- ↵
- ↵
- ↵
- ↵
- ↵
- Frank F

- ↵
- ↵
- ↵
- Hirata A, et al.

- ↵
- ↵
- ↵
- Hu YC,
- Li FX,
- Li MZ,
- Bai HY,
- Wang WH

- ↵
- Berthier L,
- Biroli G,
- Bouchaud J-P,
- Cipelletti L,
- van Saarloos W

- ↵
- Karmakar S,
- Dasgupta C,
- Sastry S

- ↵
- ↵
- ↵
- ↵
- ↵
- Kob W,
- Roldan-Vargas S,
- Berthier L

- ↵
- Scheidler P,
- Kob W,
- Binder K

- ↵
- ↵
- Kirkpatrick TR,
- Thirumalai D

- ↵
- ↵
- Karmakar S,
- Procaccia I

- ↵
- Malins A,
- Eggers J,
- Royall CP,
- Williams SR,
- Tanaka H

- ↵
- ↵
- Hu YC, et al.

- ↵
- Flenner E,
- Szamel G

- ↵
- Kirkpatrick TR,
- Thirumalai D,
- Wolynes PG

- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- Plimpton S

- ↵

## Citation Manager Formats

## Article Classifications

- Physical Sciences
- Physics