Understanding the glassy dynamics from melting temperatures in binary glass-forming liquids

Abstract

It is natural to expect that small particles in binary mixtures move faster than large ones. However, in binary glass-forming liquids with soft-core particle interactions, we observe the counterintuitive dynamic reversal between large and small particles along with the increase of pressure by performing molecular dynamics simulations. The structural relaxation (dynamic heterogeneity) of small particles is faster (weaker) than large ones at low pressures, but becomes slower (stronger) above a crossover pressure. In contrast, this dynamic reversal never happens in glass-forming liquids with hard-core interactions. We find that the difference of the effective melting temperatures felt by large and small particles can be used to understand the dynamic reversal. In binary mixtures, we derive effective melting temperatures of large and small particles simply from the conversion of units and find that particles with a higher effective melting temperature usually undergo a slower and more heterogeneous relaxation. The presence (absence) of the dynamic reversal in soft-core (hard-core) systems is simply due to the non-monotonic (monotonic) behavior of the melting temperature as a function of pressure. Interestingly, by manipulating the relative softness between large and small particles, we obtain a special case of soft-core systems, in which large particles always have higher effective melting temperatures than small ones. As a result, the dynamic reversal is totally eliminated. Our work provides another piece of evidence of the underlying connections between the properties of non-equilibrium glass-formers and equilibrium crystal-formers.

---

1 Introduction

The utilization of glassy materials in human civilization has spanned several millennia, yet a comprehensive understanding of the nature of the glass transition and glasses still poses great challenges in physics and materials science.^1–8 On approaching the glass transition, supercooled glass-forming liquids exhibit exotic but universal dynamic properties that distinguish them apparently from simple liquids. These universal properties include the super-Arrhenius behavior in the temperature evolution of the structural relaxation, the dynamic heterogeneity,^8–21 and the breakdown of the Stokes–Einstein relation,^22–32 which have undergone extensive investigations over the past decades and are still ongoing subjects of study.

Actually, there have been numerous common dynamic features reported to be shared by a wide swath of glass-formers; however, when going to a quantitative level, one could find that some dynamic properties could exhibit great diversity over glass-formers under different conditions.^33–43 For instance, some dynamic properties of glass-formers with different interaction potentials could exhibit vastly different pressure evolutions: it was reported that the glass transition temperature and fragility both increase monotonically when pressure increases in hard-core potential systems,^33,34 whereas they exhibit a non-monotonic pressure dependence in soft-core potential systems.^35–42,44

What is more complicated is that the dynamics of different particle species in the same glass-forming liquid could differ remarkably, in addition to the temporal–spatial difference in dynamics, i.e., the dynamic heterogeneity. One is that, in binary glass-forming liquids composed of small and large particles, smaller particles usually exhibit faster relaxation rates, which is commonly observed and hence naturally accepted. However, to the best of our knowledge, the reasons behind the faster relaxation of smaller particles remain largely unknown to date. As a result, we cannot yet exclude the possibility that larger particles may have a faster relaxation under certain conditions.

In this work, we systematically compare the dynamics of large and small particles in binary glass-forming liquids. Two widely studied particle interactions, soft-core harmonic and hard-core repulsive Lennard–Jones (RLJ) potentials, are studied. We find that small particles exhibit faster structural relaxations at all pressures in RLJ systems. In harmonic systems, it is the fact, only at low pressures, above a crossover pressure large particles relax unexpectedly faster instead. Interestingly, by constructing effective melting temperatures for large and small particles, respectively, we find that particles with a higher effective melting temperature can always relax more slowly. The correlation between the relaxation and the effective melting temperature of particle species rationalizes the distinct pressure evolution of particle dynamics between harmonic and RLJ systems. Using this correlation, we successfully eliminate the dynamic reversal in harmonic systems by optimizing the relative softness between large and small particles and hence making the effective melting temperature of large particles higher at all pressures.

2 Simulation details

We perform molecular dynamics simulations at constant particle number N, pressure p, and temperature T. A binary mixture of N_A A (large) and N_B B (small) particles with the same mass m is placed in a three-dimensional cubic cell with side-length L. Periodic boundary conditions are applied in all directions. We use two widely studied particle interactions, namely, harmonic (soft-core)and RLJ (hard-core)where r_ij represents the separation between particles i and j, σ_ij = (σ_i + σ_j)/2 with σ_i being the diameter of particle i, ε_ij is the characteristic energy scale of the interaction, and Θ(x) denotes the Heaviside step function. A comparison between the harmonic and RLJ potentials is illustrated in Fig. 1.

Fig. 1 Schematic illustration of the scaled harmonic and RLJ potentials U(r_ij)/ε_ij.

To fight against crystallization, we use a diameter ratio of γ = σ_A/σ_B = 1.4 and a concentration ratio quantified by the composition concentration of B particles c_B = N_B/(N_A + N_B) = 0.5. We define ε_AA = (1 + Δ)ε_AB and ε_BB = (1 − Δ)ε_AB with Δ ∈ [−1,1] in order to adjust the relative softness of particle species.^45,46 We set the units of mass, length, and energy to be m, σ_B, and ε_AB, respectively. The time and temperature are thus in units of σ_Bm^1/2ε_AB^−1/2 and ε_ABk_B⁻¹ with k_B being the Boltzmann constant.

In molecular dynamics simulations, we solve numerically equations of motion of this form:⁴⁷

where v⃑_i and r⃑_i are the velocity and position of particle i, respectively, F⃑_ij = −∇U_ij is the force of particle i exerted by particle j, and λ and ξ are the Lagrange multipliers to maintain constant p and T.⁴⁷

To quantify the structural relaxation, we calculate the self-part of the intermediate scattering function:⁴⁸

where α is A or B specifying particle species, the sum is over all α particles, r⃑_j(t) is the position of particle j at time t, and q⃑ is chosen in the x-direction with q = |q⃑| satisfying the periodic boundary conditions and being approximately the value at the first peak of the static structure factor. The structural relaxation time τ_α is determined as the time when F_s,α(t) decays to e⁻¹.

To characterize the dynamic heterogeneity, we calculate the four-point dynamic susceptibility:^49,50

which evaluates the fluctuation of an overlap function, where a = 0.3 and W_α,a(r) = 1 if r ≤ a and 0 otherwise. For glass-forming liquids, χ_4,α(t) usually exhibits a non-monotonic time dependence with a maximum χ_4,max,α. Here, we use χ_4,max,α to quantify dynamic heterogeneity at various temperatures or pressures.

The local bond-orientational order of particle i is defined as where with N_b(i) being the number of nearest neighbors of particle i determined by the Voronoi tessellation and Y_6l(r⃑_ij) the spherical harmonics. We employ 〈Q₆〉 averaged over particles and configurations to quantify the overall structural order.

The data were collected after the systems had been equilibrated for several to 100 times of τ_α (depending on the temperature and pressure). We will primarily present the results of systems with N = 1000 particles; however, we have verified partially that there are no apparent finite-size effects on our findings.

3 Results

In binary mixtures of large and small particles, the mobility of small particles is intuitively larger than that of larger ones since small particles usually possess larger free volumes. In this section, we will show that this is always true in hard-core RLJ binary glass-forming liquids. In sharp contrast, it is violated in soft-core harmonic systems at high pressures.

Fig. 2(a) and (b) compares the intermediate scattering functions of A and B particles, F_s,A(t) and F_s,B(t), for both harmonic and RLJ systems. We vary the pressure at a fixed temperature above the glass transition. Here we first show results for Δ = 0, so that ε_AA = ε_BB = ε_AB. This setup of the characteristic energy scales of the interaction has been widely used in previous modeling of glass-formers.^40,51–54 As can be seen from Fig. 2(a) and (b), F_s,B(t) (dashed lines) decays faster than F_s,A(t) (solid lines) at low pressures for both harmonic and RLJ systems, so small particles undergo a faster structural relaxation than large ones, as expected. Fig. 2(b) indicates that F_s,A(t) > F_s,B(t) holds at all pressures for RLJ systems. However, Fig. 2(a) shows that, for harmonic systems, F_s,A(t) < F_s,B(t) when the pressure is above p_d ≈ 0.18, so large particles relax counterintuitively faster than small ones. The crossover pressure p_d signals the dynamic reversal between large and small particles, at which F_s,A(t) ≈ F_s,B(t), as illustrated by the comparison of the red curves in Fig. 2(a).

Fig. 2 Intermediate scattering functions of A and B particles, F_s,A(t) (solid lines) and F_s,B(t) (dashed lines), for (a) harmonic systems at T = 0.003 and p = 0.04 (dark blue), 0.18 (red), and 0.3 (dark green), and for (b) RLJ systems at T = 0.01 and p = 0.04 (dark blue), 0.18 (red), and 0.24 (dark green). Here Δ = 0. (c) and (d) Comparison of the structural relaxation times of A and B particles, τ_A(p) (solid symbols) and τ_B(p) (empty symbols), for harmonic and RLJ systems, respectively, with the solid lines being a guide to the eye. In (c), T = 0.05 (dark blue), 0.0045 (dark green), 0.004 (red), and 0.003 (orange). In (d), T = 0.15 (dark blue), 0.015 (dark green), 0.012 (red), and 0.01 (orange). The vertical dot-dashed line in (c) marks the crossover pressure p_d at which τ_A ≈ τ_B at low temperatures.

In order to verify that the results in Fig. 2(a) and (b) are not limited to some specifically-chosen temperatures, in Fig. 2(c) and (d), we compare τ_A(p) and τ_B(p) at different temperatures. At temperatures far beyond the glass transition, the difference between τ_A(p) and τ_B(p) is tiny. With the decrease of temperature, an observable gap between τ_A(p) and τ_B(p) emerges and grows. For RLJ systems, Fig. 2(d) shows that τ_A(p) > τ_B(p) at all temperatures. However, for harmonic systems, there is always a dynamic reversal across a crossover pressure p_d, which is roughly insensitive to temperature, as shown in Fig. 2(c).

Fig. 2 raises a couple of questions: (1) what makes harmonic and RLJ systems exhibit so different pressure evolutions, and (2) what determines the crossover pressure p_d for harmonic systems? Recent studies suggest that the equilibrium melting temperature plays an unexpected but crucial role in understanding the properties of glasses, such as the glass-forming ability⁴⁵ and the effective temperature of aging glasses.⁵⁵ It is interesting to know if the melting temperature can also shed light on the distinct particle dynamics observed here. Next we will show that it actually does and in a robust way.

For monodisperse systems, the equilibrium melting temperature is well-defined and can be determined operationally as the temperature at which there is a discontinuous change in the density ρ(T), the enthalpy per particle H(T), or the average bond-orientational order 〈Q₆〉 under a constant pressure, as illustrated in Fig. 3. The three parameters are consistent in the determination of the equilibrium melting temperature.⁵⁶

Fig. 3 The temperature T dependence of the density ρ, enthalpy per particle H, and average local bond-orientational order 〈Q₆〉 for monodisperse harmonic systems at p = 0.01. The lines are a guide to the eye. The vertical dot-dashed lines mark how we determine the melting temperature, T_m,0.

For binary glass-forming liquids, it seems strange to define melting temperatures and hard to believe such temperature matter. The basic idea is that in binary mixtures large and small particles are under the same pressure p in the common units of ε_ABσ_B⁻³; however, in their own units (ε_AAσ_A⁻³ and ε_BBσ_B⁻³, respectively), A and B particles ‘feel’ different pressure values, p_A and p_B, which can be simply derived from the conversion of units:

leading toLet us denote T_m,0 as the equilibrium melting temperature of monodisperse systems. Assume there is a monodisperse system of A particles at pressure p_A (in units of ε_AAσ_A⁻³). The melting temperature is thus T_m,0(p_A) (in units of ε_AAk_B⁻¹). Applying the same approach, we have T_m,0(p_B) (in units of ε_BBk_B⁻¹) for B particles. By converting T_m,0(p_A) and T_m,0(p_B) to the common units of ε_ABk_B⁻¹, we can define effective melting temperatures ‘felt’ by particles A and B in binary mixtures at the common pressure p, T_m,A(p) and T_m,B(p):

T_m,A(p) × ε_ABk_B⁻¹ = T_m,0(p_A) × ε_AAk_B⁻¹, (12)

T_m,B(p) × ε_ABk_B⁻¹ = T_m,0(p_B) × ε_BBk_B⁻¹. (13)

The combination of eqn (10)–(13) results inApparently, T_m,A(p) and T_m,B(p) can be simply obtained from T_m,0(p). Note that the procedures used to get eqn (14) and (15) as well as eqn (17) have been discussed in our previous work.⁴⁵

For the case with Δ = 0 shown in Fig. 2, T_m,A(p) = T_m,0(γ³p) and T_m,B(p) = T_m,0(p). Therefore, T_m,A(p) can be simply obtained from T_m,0(p) by multiplying p by γ⁻³. Fig. 4(a) and (b) compares T_m,A(p) and T_m,B(p) for harmonic and RLJ systems, respectively, with Δ = 0. The remarkable difference between harmonic and RLJ systems is that T_m,0(p) is monotonic for RLJ, but non-monotonic for harmonic, due to their soft- and hard-core natures, respectively.⁴⁵ Actually, the non-monotonic evolution of the melting temperature with pressure in soft-core potential systems has been studied intensively.^42,56–63 Because of this difference, the relative standing between T_m,A(p) and T_m,B(p) shows rather different pressure evolutions between the two systems. For RLJ systems, T_m,A is always larger than T_m,B at all pressures, simply because T_m,0(p) is monotonic. In contrast, because of the non-monotonicity, the T_m,A(p) and T_m,B(p) of harmonic systems intersect at a crossover pressure p_m. When p < p_m, T_m,A > T_m,B, while T_m,A < T_m,B otherwise in the pressure regime studied here.

Fig. 4 Effective melting temperatures of A and B particles, T_m,A(p) (squares) and T_m,B(p) (circles), for (a) harmonic and (b) RLJ systems with Δ = 0. The vertical dot-dashed line in (a) marks the crossover pressure p_m at which T_m,A(p) = T_m,B(p).

It is interesting that the effective melting temperature difference strongly correlates with the relaxation time difference between the two particle species. For RLJ systems, T_m,A(p) > T_m,B(p) at all pressures, and we always have τ_A(p) > τ_B(p). For harmonic systems, T_m,A(p) > T_m,B(p) and τ_A(p) > τ_B(p) at low pressures, while T_m,A(p) < T_m,B(p) and τ_A(p) < τ_B(p) at high pressures. As shown in Fig. 4(a), T_m,A = T_m,B at p_m ≈ 0.18, in excellent agreement with p_d ≈ 0.18 at which τ_A = τ_B. All these suggest that particle species with a higher effective melting temperature undergo a slower structural relaxation.

Fig. 2 and 4 show the results for Δ = 0. As can be seen from eqn (14) and (15), T_m,A(p) and T_m,B(p) vary with Δ, so does the crossover pressure p_m. In order to see whether there is always a dynamic reversal in harmonic systems and whether p_d always agrees with p_m, we present results for various Δ values in Fig. 5. For two arbitrarily chosen values of Δ, − 0.35 and 0.2, Fig. 5(a), (b), (d) and (e) indicates that T_m,A(p) and T_m,B(p) intersect at a Δ-dependent p_m; meanwhile, the dynamic reversal always occurs at p_d(Δ) ≈ p_m(Δ). The robust agreement between p_d(Δ) and p_m(Δ) convincingly supports our argument about the correlation between effective melting temperature and particle dynamics.

Fig. 5 Manipulation of the dynamic reversal in harmonic systems by tuning the relative softness of A and B particles Δ. The top row compares the melting temperatures T_m,A(p) (squares) and T_m,B(p) (circles) in systems with (a) Δ = −0.35, (d) Δ = 0.2, and (g) Δ = Δ* = 0.466, respectively. The vertical dot-dashed lines in (a) and (c) mark the pressure p_m at which T_m,A(p) = T_m,B(p). The middle row compares the structural relaxation times τ_A(p) (squares) and τ_B(p) (circles) at T = 4 × 10⁻³ in systems with (b) Δ = −0.35, (e) Δ = 0.2, and (h) Δ = Δ* = 0.466, respectively. The vertical dot-dashed lines in (b) and (d) mark the pressure p_d at which τ_A(p) = τ_B(p). The bottom row compares the maxima of the four-point dynamic susceptibility χ_4,max,A(p) (squares) and χ_4,max,B(p) (circles) at T = 4 × 10⁻³ in systems with (c) Δ = −0.35, (f) Δ = 0.2, and (i) Δ = Δ* = 0.466, respectively. The vertical dot-dashed lines in (c) and (f) mark the pressure p_d.

The non-monotonicity causes the intersection between T_m,A(p) and T_m,B(p). However, it is not guaranteed that T_m,A(p) and T_m,B(p) have to intersect. There is a special value of Δ, denoted as Δ*, which can be easily derived from the condition of p_A = p_B, i.e., by equating eqn (10) and (11):

By substituting Δ* into eqn (14) and (15), we have

T_m,A(p) = γ³T_m,B(p). (17)

Apparently, when Δ = Δ*, T_m,A(p) is always larger than T_m,B(p) without any intersection, as shown in Fig. 5(g). Interestingly, Fig. 5(h) indicates that τ_A > τ_B at all pressures and the dynamic reversal is completely eliminated, even though the melting temperature is non-monotonic in pressure. This result is another strong evidence to support the correlation between effective melting temperature and particle dynamics.

We have compared the dynamics of large and small particles by calculating their structural relaxation times and revealed a connection between the difference in the dynamics and the difference in the effective melting temperatures. In addition, the dynamics usually exhibits spatio-temporal fluctuations in glass-forming liquids, which is referred to as dynamic heterogeneity. It has been well documented^40,49 that a slower dynamics usually corresponds to a more heterogeneous dynamics. Since there could be a reversal in the structural relaxation times between large and small particles under certain circumstances, it is interesting to see whether there could be a reversal in the dynamic heterogeneity as well. Furthermore, if the reversal in dynamic heterogeneity exists, does it start at P_m?

We use the maximum of the four-point dynamic susceptibility χ_4,max,A and χ_4,max,B to quantify the dynamic heterogeneity of A and B particles, respectively. The comparison of χ_4,max,A and χ_4,max,B is shown in Fig. 5(c) with Δ = −0.35, Fig. 5(f) with Δ = 0.2, and Fig. 5(i) with Δ = Δ*. Interestingly, the evolution of the difference between χ_4,max,A(p) and χ_4,max,B(p) with pressure follows nearly the same trend as that of the difference between T_m,A(p) and T_m,B(p), in particular χ_4,max,A(p) = χ_4,max,B(p) at p_m where T_m,A(p) = T_m,B(p). Therefore, one could conclude that there is also a good correlation between the effective melting temperature and the dynamic heterogeneity.

It is worth noting that we have excluded that the dynamic reversal observed in harmonic systems is a direct consequence of the particle demixing; see Fig. 6 and related discussions in the Appendix. Moreover, in Fig. 7 in the Appendix, we also find that our major conclusions do not depend on the diameter ratio or the concentration ratio of A and B particles within an appropriate range. However, in extreme cases where the diameter ratio or the concentration ratio is large and hence the good glass-forming ability⁶⁴ could not be maintained anymore, we find that our conclusions could not work. Therefore, we believe more work is needed to study systems in these extreme cases.

Fig. 6 Correlation between the dynamic reversal and particle demixing. (a) Comparison of the structural relaxation times τ_A(p) (squares) and τ_B(p) (circles) at T = 4 × 10⁻³ in systems with Δ = 0. The vertical dashed line in (a) marks the pressure p_d at which τ_A(p) = τ_B(p). (b) Comparison of the particles’ aggregation degree Ω_A (squares) and Ω_B (circles) for the same systems in (a). (c) The evolution of snapshots with pressure for the same systems in (a). The blue and yellow spheres are A (large) and B (small) particles.

Fig. 7 Effect of the particle concentration c_B and diameter ratio γ on the dynamic reversal in harmonic systems with Δ = 0. The top row compares the melting temperatures T_m,A(p) (squares) and T_m,B(p) (circles) in systems with (a) γ = 1.4, (c) γ = 1.6, and (e) γ = 1.8, respectively. The vertical dot-dashed lines in (a), (c) and (e) mark the pressure p_m at which T_m,A(p) = T_m,B(p). (b) Comparison of the structural relaxation times τ_A(p) (squares) and τ_B(p) (circles) at T = 4.3 × 10⁻³ in systems with the same γ = 1.4 but different c_B = 0.4 (red), 0.5 (blue) and 0.6 (orange), respectively. (d) and (f) Comparison of the structural relaxation times τ_A(p) (squares) and τ_B(p) (circles) at T = 4 × 10⁻³ in systems with γ = 1.6 and γ = 1.8, respectively; c_B = 0.5 for (d) and (f). The vertical dot-dashed lines in (b), (d) and (f) mark the pressure p_d at which τ_A(p) = τ_B(p).

It is also worthwhile to note that this work demonstrates a strong correlation between the glassy dynamics and the effective melting temperatures. In our previous work,⁴⁵ the glass-forming ability has been linked to the gap in the effective melting temperatures. To our knowledge, the glassy dynamics and the glass-forming ability are two facets of the glass transition and hence usually studied separately. However, the glass-forming ability could relate closely to the glassy dynamics if combining our results in this work and those in ref. 45, which warrants further studies to investigate.

4 Conclusions

In summary, we compare the dynamics of large and small particles in binary glass-forming liquids with two distinct particle interactions. For soft-core harmonic systems, we observe the unexpected and counterintuitive dynamic reversal at high pressures: small particles relax faster and exhibit weaker dynamic heterogeneity than large ones at low pressures, while large particles relax faster and exhibit weaker dynamic heterogeneity at high pressures. In contrast, the dynamic reversal is absent in hard-core RLJ systems, in which small particles normally undergo a faster relaxation and less heterogeneous dynamics. We find that the difference in effective melting temperatures between large and small particles can rationalize the anomalous dynamic reversal in harmonic systems and the distinct dynamics between the two systems. Our results convincingly suggest that particle species with a higher effective melting temperature will undergo a slower structural relaxation in binary glass-forming liquids. Therefore, the non-monotonic pressure dependence of the melting temperature is the necessary condition for the dynamic reversal to occur. These findings provide us the freedom to manipulate the occurrence of the dynamic reversal and even eliminate it in soft-core systems with a non-monotonic pressure dependence of the melting temperature.

It is surprising that the equilibrium melting temperature of monodisperse crystal-formers can take effect to control the dynamics of binary glass-formers. The correlation between the effective melting temperature and the particle dynamics revealed by this work is another piece of evidence of the underlying connections between non-equilibrium systems and their equilibrium counterparts.^45,55 Here we are concerned about binary mixtures. Even for polydisperse harmonic systems, if we divide all particles into larger and smaller ones, separated by the average diameter, we expect to see the dynamic reversal as well and interpret it in terms of melting temperatures.

Author contributions

N. X. designed the project. Y. N. performed the simulations. Y. N., L. W., P. G., and N. X. analyzed the data and wrote the paper.

Conflicts of interest

There are no conflicts to declare.

Appendix

In order to check whether the reported dynamic reversal observed in harmonic systems is due to the phase separation, we quantify the degree of mixing for A particles and B particles using⁴⁵where the sum is over all α particles, δ_nsi,ni is the Kronecker delta, n_i is the number of nearest neighbors of particle i, and n_si is the number of nearest neighbors which are the same type as particle i. Ω_A and Ω_B are close to 0 when A and B components mix up very well, whereas larger values of Ω_A and Ω_B mean stronger demixing.

Fig. 6 compares the structure relaxation times of A and B, τ_A and τ_B, in (a) and the mixing degree of A and B particles, Ω_A and Ω_B, in (b) for harmonic systems for Δ = 0. In addition, we also show some snapshots at pressures below, around, and above the pressure at the dynamic reversal. We find it is hard to conclude there is a conclusively direct connection between the dynamic difference and the difference in the mixing degree. However, the demixing (or partial phase separation) seems to become stronger with the increase of pressure at very high pressures,⁵⁵ and thus future work should check how this affects the particle dynamics.

In the main text, we have focused on γ = 1.4 and c_B = 0.5. Fig. 7 shows our further results for c_B = 0.4 and c_B = 0.6 at fixed γ = 1.4, and for γ = 1.6 and γ = 1.8 at fixed c_B = 0.5. For each combination of (γ, c_B), the agreement between particle dynamics and effective melting temperature could be reproduced. However, we find that the agreement will disappear when γ or c_B becomes large, which deserves further study.

Acknowledgements

We are grateful to L. Berthier for bringing us attention to the dynamic reversal in harmonic systems. Y. N. and P. G. acknowledge the support from the National Natural Science Foundation of China (Grant no. T2325004, 522161160330 and U2230402). N. X. acknowledges the support from the National Natural Science Foundation of China (Grant no. 12334009 and 12074355). Y. N. acknowledges the support from the China Postdoctoral Science Foundation (Grant no. 2021M690328). L. W. acknowledges the support from the National Natural Science Foundation of China (Grant no. 12374202 and 12004001), Anhui Projects (Grant no. 2022AH020009, S020218016, and Z010118169), and Hefei City (Grant no. Z020132009). We also acknowledge Beijing Computational Science Research Center, the Supercomputing Center of the University of Science and Technology of China, Hefei Advanced Computing Center, Beijing Super Cloud Computing Center, and the High-Performance Computing Platform of Anhui University for providing computing resources.

References

1. C. A. Angell, Science, 1995, 267, 1924.
2. V. Lubchenko and P. G. Wolynes, Annu. Rev. Phys. Chem., 2007, 58, 235.
3. P. G. Debenedetti and F. H. Stillinger, Nature, 2001, 410, 259.
4. J. C. Dyre, Rev. Mod. Phys., 2006, 78, 953.
5. L. Berthier and G. Biroli, Rev. Mod. Phys., 2011, 83, 587.
6. S. Karmakar, C. Dasgupta and S. Sastry, Annu. Rev. Condens. Matter Phys., 2014, 5, 255.
7. D. Chandler and J. P. Garrahan, Annu. Rev. Phys. Chem., 2010, 61, 191.
8. M. D. Ediger, Annu. Rev. Phys. Chem., 2000, 51, 99.
9. E. Flenner, H. Staley and G. Szamel, Phys. Rev. Lett., 2014, 112, 097801.
10. M. Ozawa and G. Biroli, Phys. Rev. Lett., 2023, 130, 138201.
11. M. Adhikari, S. Karmakar and S. Sastry, J. Phys. Chem. B, 2021, 125, 10232.
12. G. Biroli, K. Miyazaki and D. R. Reichman, arXiv, 2022, preprint, arXiv:2209.02825,  DOI:10.48550/arXiv.2209.02825.
13. A. Tahaei, G. Biroli, M. Ozawa, M. Popović and M. Wyart, Phys. Rev. X, 2023, 11, 031034.
14. R. N. Chacko, F. P. Landes, G. Biroli, O. Dauchot, A. J. Liu and D. R. Reichman, Phys. Rev. Lett., 2021, 127, 048002.
15. Y. Li, Y. Yao and M. P. Ciamarra, Phys. Rev. Lett., 2022, 128, 258001.
16. H. Sillescu, J. Non-Cryst. Solids, 1999, 243, 81.
17. L. Berthier, G. Biroli, H.-P. Bouchaud, L. Cipelletti and W. van Saarloos, Dynamical Heterogeneities in Glasses, Colloids, and Granular Media, OUP Oxford, 2011, vol. 150.
18. W. Kob, C. Donati, S. J. Plimpton, P. H. Poole and S. C. Glotzer, Phys. Rev. Lett., 1997, 79, 2827.
19. D. N. Perera and P. Harrowell, Phys. Rev. Lett., 1998, 81, 120.
20. R. Yamamoto and A. Onuki, Phys. Rev. Lett., 1998, 81, 4915.
21. C. Donati, S. Franz, S. C. Glotzer and G. Parisi, J. Non-Cryst. Solids, 2002, 307, 215.
22. E. Rössler, Phys. Rev. Lett., 1990, 65, 1595.
23. F. H. Stillinger and J. A. Hodgdon, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 1994, 50, 2064.
24. M. T. Cicerone and M. D. Ediger, J. Chem. Phys., 1996, 104, 7210.
25. L. Berthier, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2004, 69, 020201.
26. Y. Jung, J. P. Garrahan and D. Chandler, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2004, 69, 061205.
27. S. K. Kumar, G. Szamel and J. F. Douglas, J. Chem. Phys., 2006, 124, 214501.
28. S. R. Becker, P. H. Poole and F. W. Starr, Phys. Rev. Lett., 2006, 97, 055901.
29. S.-H. Chong and W. Kob, Phys. Rev. Lett., 2009, 102, 025702.
30. S. Sengupta, S. Karmakar, C. Dasgupta and S. Sastry, J. Chem. Phys., 2013, 138, 12A548.
31. P. Charbonneau, Y. Jin, G. Parisi and F. Zamponi, Proc. Natl. Acad. Sci. U. S. A., 2014, 111, 15025.
32. A. D. Parmar, S. Sengupta and S. Sastry, Phys. Rev. Lett., 2017, 119, 056001.
33. S. Sastry, Nature, 2001, 409, 164.
34. L. Wang and N. Xu, Phys. Rev. Lett., 2014, 112, 055701.
35. G. Foffi, F. Sciortino, P. Tartaglia, E. Zaccarelli, F. L. Verso, L. Reatto and C. N. Likos, Phys. Rev. Lett., 2003, 90, 238301.
36. E. Zaccarelli, C. Mayer, A. Asteriadi, C. N. Likos, F. Sciortino, J. Roovers, H. Iatrou, N. Hadjichristidis, P. Tartaglia, H. Lowen and D. Vlassopoulos, Phys. Rev. Lett., 2005, 95, 268301.
37. C. Mayer, E. Zaccarelli, E. Stiakakis, C. N. Likos, F. Sciortino, A. Munam, M. Gauthier, N. Hadjichristidis, H. Iatrou, P. Tartaglia, H. Lowen and D. Vlassopoulos, Nat. Mater., 2008, 7, 780.
38. C. M. Likos, Phys. Rep., 2001, 348, 267.
39. L. Berthier, A. J. Moreno and G. Szamel, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2010, 82, 060501(R).
40. L. Wang, Y. Duan and N. Xu, Soft Matter, 2012, 8, 11831.
41. M. Schmiedeberg, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2013, 87, 052310.
42. R. Miyazaki, T. Kawasaki and K. Miyazaki, Phys. Rev. Lett., 2016, 117, 165701.
43. X. Sun, H. Zhang, L. Wang, Z. Zhang and Y. Ma, Chin. Phys. B, 2020, 129, 126201.
44. R. Miyazaki, T. Kawasaki and K. Miyazaki, J. Chem. Phys., 2019, 150, 074503.
45. Y. Nie, J. Liu, J. Guo and N. Xu, Nat. Commun., 2020, 11, 3198.
46. J. Guo, Y. Nie and N. Xu, Soft Matter, 2021, 17, 3397.
47. M. P. Allen and D. J. Tildesley, Computer Simulation of Liquids, Oxford University Press, Oxford, 1987.
48. W. Kob and H. C. Andersen, Phys. Rev. Lett., 1994, 73, 1376.
49. N. Lačević, F. W. Starr, T. B. Schrøder and S. C. Glotzer, J. Chem. Phys., 2003, 119, 7372.
50. D. Coslovich and C. M. Roland, J. Chem. Phys., 2009, 131, 151103.
51. C. Zhao, K. Tian and N. Xu, Phys. Rev. Lett., 2011, 106, 125503.
52. N. Xu, T. K. Haxton, A. J. Liu and S. R. Nagel, Phys. Rev. Lett., 2009, 103, 245701.
53. L. Wang, P. Guan and W. H. Wang, J. Chem. Phys., 2016, 145, 034505.
54. L. Wang, N. Xu and P. Guan, Phys. Rev. Lett., 2018, 120, 125502.
55. J. Zhang, W. Zheng, S. Zhang, D. Xu, Z. Jiang and N. Xu, Sci. Adv., 2021, 7, eabg6766.
56. M. Zu, J. Liu, H. Tong and N. Xu, Phys. Rev. Lett., 2016, 117, 085702.
57. S. Prestipino, F. Saija and P. V. Giaquinta, Phys. Rev. Lett., 2011, 106, 235701.
58. Y. D. Fomin, N. V. Gribova, V. N. Ryzhov, S. M. Stishov and D. Frenkel, J. Chem. Phys., 2008, 129, 064512.
59. J. C. Pamies, A. Cacciuto and D. Frenkel, J. Chem. Phys., 2009, 131, 044514.
60. D. E. Dudalov, Y. D. Fomin, E. N. Tsiok and V. N. Ryzhov, Soft Matter, 2014, 10, 4966.
61. N. Xu, Chin. J. Polym. Sci., 2019, 37, 1065.
62. M. Zu, T. Peng and N. Xu, Nat. Commun., 2017, 8, 2089.
63. W. L. Miller and A. Cacciuto, Soft Matter, 2011, 7, 7552.
64. J. Russo, F. Romano and H. Tanaka, Phys. Rev. X, 2018, 8, 021040.

---