Relationship between two-level systems and quasi-localized normal modes in glasses
RRelationship between two-level systems and quasi-localized normal modes in glasses
Dmytro Khomenko, David R. Reichman, and Francesco Zamponi Department of Chemistry, Columbia University, New York, NY 10027, USA Laboratoire de Physique de l’Ecole Normale Sup´erieure, ENS, Universit´e PSL,CNRS, Sorbonne Universit´e, Universit´e de Paris, 75005 Paris, France
Tunnelling Two-Level Systems (TLS) dominate the physics of glasses at low temperatures. YetTLS are extremely rare and it is extremely difficult to directly observe them in silico . It is thuscrucial to develop simple structural predictors that can provide markers for determining if a TLSis present in a given glass region. It has been speculated that Quasi-Localized vibrational Modes(QLM) are closely related to TLS, and that one can extract information about TLS from QLM. Inthis work we address this possibility. In particular, we investigate the degree to which a linear ornon-linear vibrational mode analysis can predict the location of TLS independently found by energylandscape exploration. We find that even though there is a notable spatial correlation between QLMand TLS, in general TLS are strongly non-linear and their global properties cannot be predicted bya simple normal mode analysis.
I. INTRODUCTION
The thermodynamic behavior of low-temperatureamorphous solids has been a topic of considerable in-terest since the seminal experiments of Zeller and Pohldemonstrated marked deviations between the thermalproperties of glasses and those of crystalline solids, nearlyfifty years ago [1]. In particular, in the range of onedegree Kelvin, these experiments and numerous experi-ments that have followed [2–10] have demonstrated thatthe specific heat of a disordered solid is much larger thanthat of a crystal composed of the same material, witha temperature dependence that varies linearly, as op-posed to cubically, and a thermal conductivity that variesquadratically, as opposed to cubically, with temperature.A theoretical framework provided independently by An-derson, Halperin and Varma [11] and by Phillips [12, 13],ascribes the origin of this puzzling behavior to dilute de-fects which tunnel between their lowest lying quantumstates at low temperatures. This two-level systems (TLS)picture has successfully rationalized diverse experimen-tally observed properties, although several outstandingpuzzles lie beyond its reach [13–15]. In particular, under-standing the microscopic nature of the TLS, as well asthe quasi-universal aspects of the thermodynamic data,have remained as outstanding challenges [16–18].Over the last five years, the swap Monte Carlo tech-nique has provided a major advance in computationalglass physics, opening the door to the creation of in sil-ico glasses that have comparable stability properties tothose found in the laboratory [19]. Using this technique,and building on previously developed landscape explo-ration algorithms [20–29], we recently provided a de-tailed computational investigation of the TLS model [30],providing a direct microscopic description of TLS anddemonstrating that their density decreases as the cool-ing rate decreases, similar to what is seen in several re-cent experiments [6–10]. Although this work only con-sidered one model system, and thus the question ofthe quasi-universality [3] of low-temperature thermody- namic anomalies could not be investigated, a detailed de-scription of the nature of TLS was provided. Specifically,tunneling motion in TLS was found to be comprised ofdefect-vacancy-like motion of one or a handful of par-ticles, although occasionally highly collective tunnelingmotion of a large number of particles was observed.A major issue with the landscape exploration methodscurrently used to identify TLS in silico [22, 25, 29, 30]is that they are computationally very expensive. Thisbottleneck is due to (i) the need to accumulate a suf-ficient number of inherent structures (IS) [20, 21] andthen, for at least the most promising pairs of IS (see [30]for details), (ii) the need to identify a relaxation path-way between the two minima in the 3 N -dimensional en-ergy landscape achieved via a computationally expensiveminimization of a path function in the space of possiblepaths [25, 31–34]. Hence, it would be extremely helpful toidentify a priori the glassy configurations (or sub-regionsof them) that are most likely to include TLS with theproper energy splitting, via some sort of simple struc-tural indicator.Over the last two decades a seemingly distinct typeof (partially) localized entity, namely quasi-localized vi-brational modes (QLM) [35], have been intensely scru-tinized. QLM are characterized by a defect-like local-ized core with a power-law decaying elastic background,and are prominently found in the low-frequency wing ofthe density of states of amorphous systems. They havebeen connected to the universal non-Debye behavior ofthe low-frequency density of states [36–43] and to the at-tenuation of sound waves in glassy systems [44–47], tothe dynamical heterogeneity upon approach to the glasstransition from the high temperature side [48, 49], to theplasticity of the glass under strain [50–52], and to thecritical-like behavior in jamming [53, 54].It is natural to assume a connection between QLM andTLS. Indeed, it is known from computer simulations thatas model supercooled liquids are cooled, the concentra-tion of real-space localized cores associated with QLMrapidly decreases [40, 41], as is also the case for TLS.It is thus reasonable to assume that at the glass transi- a r X i v : . [ c ond - m a t . d i s - nn ] D ec local minimum A local minimum Bsaddle point FIG. 1. Schematic illustration of the multidimensional mini-mal energy transition path between two energy minima withina double-well potential. Also illustrated are the displacementvector r AB between the two minima, the displacement vector r AS between the first minimum and the saddle point, and thetangent vector r in the first minimum. tion these QLM cores are “frozen” into the sample andprovide the seeds for low-temperature defects [55]. In-deed, this notion is central to the successful soft potentialmodel [56–60] of low-temperature glasses, which extendsthe models of Anderson-Halperin-Varma and of Phillipsto somewhat higher temperatures by connecting TLS toanharmonic vibrational modes in the glass. If this con-nection is precise, QLM could be used as structural pre-dictors for the location of TLS in glass samples, thus aid-ing the computational search for tunneling states. More-over, establishing this connection more precisely couldhelp validate or invalidate models of low-temperatureglasses based on interacting anharmonic modes [59–63].In this work, leveraging our ability to prepare realis-tically cooled samples and extract detailed informationabout both TLS and QLM, we explore their putativeconnection in detail. In our system we find that, in realspace, TLS and QLM are well correlated, in the sensethat particles that move the most in a TLS are typicallyclose to particles that move the most in a QLM. On theother hand, in phase space, we demonstrate a surpris-ingly weak correlation between TLS and both linear andnon-linear QLM: the 3 N -dimensional vector r AB that en-codes the displacement in phase space of all particles ina TLS is often completely unrelated to the vectors thatdefine QLMs. II. SYSTEM
Our study is based on our previous work [30], in whichwe prepared in silico glasses at cooling rates, from poorlyannealed to ultrastable, and explored their energy land-scape. For completeness, we will give a very brief sum-mary of the methodology, but we refer to Ref. [30] fordetails.Our system is a non-additive polydisperse mixture of N = 1500 particles. In the following, r denotes a 3 N -dimensional phase space vector encoding the position ofall particles, r i denotes the three-dimensional coordinate of particle i = 1 , · · · , N , and r ij = | r i − r j | is the scalardistance between particles i and j . We use the followinginter-particle interaction potential: v ij ( r ij ) = (cid:15) (cid:16) σ ij r ij (cid:17) + (cid:15)F (cid:16) σ ij r ij (cid:17) , r ij < r cut , , r ij > r cut , (1)where r cut = 1 . σ ij , and σ ij is the non-additive inter-action length scale associated to particle pair ij . Thefunction F ( x ) is a fourth-order polynomial which guar-antees the continuity of the potential up to its secondderivative at r cut . We express all dimensional quantitiesin units of energy (cid:15) , length (cid:104) σ (cid:105) = 1 (the average beingover particle pairs), time (cid:112) (cid:15)/ ( m (cid:104) σ (cid:105) ), and the numberdensity is set to ρ = 1 in these units. The mode-couplingtemperature, which sets the onset of strongly glassy dy-namics, is T MCT = 0 .
104 [19].Using the swap Monte Carlo algorithm [19] we pre-pare fully equilibrated configurations at three differ-ent preparation temperatures T f = 0 .
062 (ultra-stableglasses), 0.07 (liquid cooled glasses), 0.092 (poorly an-nealed glasses). Normal molecular dynamics (MD) ini-tialized in these configurations is fully arrested (exceptfor T f = 0 . T f corresponds to Tool’s “fictive temperature” [64], and en-codes the degree of glass stability. For each of the glasses,we explore the glass basin in the energy landscape and de-termine a set of local energy minima, or “inherent struc-tures” (IS), by running MD at a slightly lower temper-ature T MD = 0 .
04 and periodically minimizing the sys-tem’s potential energy [20, 21]. Pairs of energy minimathat are subsequently visited a large enough number oftimes (see [30] for details) are further analysed by theNudged Elastic Band (NEB) method [31, 32] to find theminimum energy pathway connecting them, see Fig. 1,and the associated value of the classical energy barrier.This procedure allows us to obtain a library of distinctDouble Well potentials (DWs) in the high-dimensionalenergy landscape. For each DW, we then solve an ef-fective one-dimensional Schr¨odinger equation to obtainthe quantum tunnel splitting E = E − E between thefirst two energy levels. We find that the relevant DWs,which define active TLS in the quantum regime (namelythose with a tunnel splitting equal to or below the tem-perature T Q which defines the low temperature regimein [30]), have E < . E < . III. NORMAL MODES
Having curated a library of DWs along with associ-ated displacement fields r AB , we can perform a normalmode analysis in each of the two minima A and B (be-cause our search procedure is statistically symmetric, in -3 -2 -1 FIG. 2. Participation ratio of normal modes versus their fre-quency, for one selected energy minimum. Modes shown inred stars are considered to be quasi-localised. In the inset,the mode frequency is shown as a function of the mode index. the following we focus on minimum A without loss ofgenerality), and check if the normal modes overlap withthe DW displacement field, to be defined below. We nowdefine more precisely the linear and non-linear normalmodes and their relationship with the minimum energypath illustrated in Fig. 1.Several displacement fields can be associated with aDW transition: the difference between the two energyminima, r AB , the difference between the first minimumand the saddle point, r AS , and the tangent direction r of the minimal energy path in A (see Fig. 1). The vector r is estimated by a discretization of the minimum en-ergy path, as the difference between the position of thefirst two beads (or images) of the NEB path.Next, we define the tensors M and U as follows: M αβ = ∂ V∂r α ∂r β , U αβγ = ∂ V∂r α ∂r β ∂r γ , (2)where V ( r ) = (cid:80) i We now report the results of the calculations describedabove, for all DWs in the data set obtained in Ref. [30].For T f = 0 . , . , . , , , , 248 are activeusing NiP units, respectively.As a measure of overlap between a displacement field r and a normal mode v , we will focus on two quantities:the simple normalized scalar product of 3 N -dimensionalvectors, c ( v ) = | r · v || r || v | , (7)and the scalar product of N -dimensional vectors obtainedby collecting the absolute values of particle displace-ments, a ( v ) = (cid:80) Ni =1 | r i || v i | (cid:113)(cid:80) Ni =1 | r i | (cid:80) Ni =1 | v i | , (8)where r i and v i are the displacements of particle i in the r and v vectors. The parameter a ignores the polarizationof vectors and compares only the mobility of particles,and hence is an analog of the “softness” field used in [40]for a single mode. A. Tangent vectors are parallel to a soft QLM For a given energy minimum, we can compute all thelinear modes. Their participation ratio as a function offrequency is given in Fig. 2, and an example of a QLM isgiven in Fig. 3.Our first result is that, in the limit of an infinite num-ber of NEB images, ( N b → ∞ ), when the path becomescontinuous, the tangent vector r to the minimal energypath in minimum A coincides with one of the softest lin-ear modes in the same minimum. To prove this, for aselected DW, we show in Fig. 4 the scalar product c α between r and all the vibrational modes in minimum A , here labeled by α = 1 , · · · , N . Note that the firstthree modes are trivial zero modes related to transla-tional invariance, hence the non-trivial modes are labeledby α = 4 , · · · , N . We clearly see that, upon increasing N b from 40 to 600, the overlap with mode 4, which isa QLM, increases while all the other overlaps decrease.Hence, we conclude that in the limit N b → ∞ , r be-comes essentially parallel to a soft QLM. Note that mostof our simulations have been conducted with N b = 40,and that increasing the number of beads to N b = 600makes the NEB calculation very computationally expen-sive, forcing us to restrict this investigation to a smallnumber of DWs.To provide further support for this statement, in Fig. 5we report the probability distribution of the maximumoverlap, c = max α c α , between linear modes and r ,over a subset of DW potentials, as described in the figurecaption. We observe that in all cases, c > . 4, and thatin most cases c is quite close to one, which confirms thatthe tangent vector is indeed parallel to a linear mode.Note that the results of Fig. 5 are for N b = 40, and weexpect c to increase upon increasing N b .The fact that r is parallel to a linear mode, typicallya soft QLM, implies that in real space there is always asoft QLM whose localized core is close to the particlesthat move the most in the DW. Hence, soft QLM aregood predictors of the spatial location of DW potentials. r r AS r AB FIG. 5. Distribution over DW potentials of the largest mode projection, c = max α c α , calculated over all modes for r , r AS and r AB . Here we used N b = 600, and since the NEB calculation is expensive, we have used a smaller subset of all DWs, onlyselecting TLS with energy splitting E < − and participation ratio of the DW transition P R DW < P R forDWs is normalized differently than for normal modes, see [30]) for T f = 0 . -4 -2 FIG. 6. Scatter plot of the tunnel splitting E versus the fre-quency ω in minimum A , for all DWs found at T f = 0 . B. The frequency in a minimum is anti-correlatedwith the tunnel splitting Another interesting observation concerns the relationbetween the curvature ω of the energy profile along thetransition path (which, as discussed in section IV A, co-incides with the frequency of a soft QLM) and the tun-nel splitting E associated to the DW. In Fig. 6 we re-port a scatter plot of these two quantities, which showsa marked anticorrelation.We thus conclude that, although the vector r isstrongly associated with a soft linear mode, its frequencyis not among the softest, and in particular the DWs withlowest splitting are associated to relatively higher fre-quencies. This behavior likely arises due to the fact thatactive TLS typically display a symmetric DW profile,with a relatively high barrier and hence a relatively highfrequency of the two wells. C. Linear modes are poor predictors of theminimal energy path curvature We now discuss whether the linear modes of minimum A are good predictors of the transition path associatedwith a DW, as encoded by the minimum-to-saddle dis-placement r AS and by the minimum-to-minimum dis-placement r AB . In Fig. 5 we compare the statistics ofthe maximum overlap coefficient c = max α c α of linearmodes with r , r AS , and r AB . From these figures weclearly see that the values of c for r AS and r AB are muchsmaller than for r . We thus conclude that the mini-mal energy path between two energy minima, which (asshown in section IV A) tends to start along one of thesoftest modes locally around each minimum, markedlychanges its direction upon approaching the saddle point.We find that linear modes are poor predictors of thischange of direction. D. Non-linear modes are better correlated withthe minimum energy path curvature We next consider whether non-linear modes can be bet-ter predictors of the direction of the minimum energypath around the saddle point S or the arrival minimum B . Because most of the data we will show are quali-tatively similar for r AB and r AS , we will focus on theformer for the rest of this section.In order to find the closest non-linear mode to the min-imum energy path, we use π = r AB / | r AB | as an initialguess for the iterative procedure in Eq. (6), and we it-erate until convergence to the corresponding non-linearmode. We then compute the overlap coefficients c and a between the non-linear mode and r AB , defined respec-tively in Eq. (7) and Eq. (8). In order to provide a directcomparison with linear modes, we find the linear modethat has the maximum overlap ( a or c ) with r AB . Thestatistics of a and c is shown in Fig. 7 for the full set ofavailable DWs at the three T f values. From these plots, FIG. 7. Probability distribution of overlap coefficients c and a for linear and non-linear modes, obtained from a recursiveprocedure, starting from r AB as initial guess. Data are for the three preparation temperatures T f and for the full set of DWs.FIG. 8. Overlap of the softest localized linear modes, and ofthe non-linear modes that are constructed using those linearmodes as an initial condition, with the displacement field r AB .Data are for T f = 0 . 062 and for the full set of DWs. one can see that non-linear modes generally have a muchstronger overlap with r AB than do the linear modes, butnevertheless, for a large fraction of DWs, the overlap re-mains small even for non-linear modes. It is very impor-tant to stress that in this analysis we used the a priori known information encoded in r AB as an initial guessto search for the closest non-linear mode to the reactionpath. Hence, we expect the values reported in Fig. 7 toprovide an upper bound on the possible overlaps.We have repeated the same analysis, without assum-ing any a priori knowledge about the second minimum. We first diagonalize the Hessian matrix in minimum A ,and identify the softest mode, which we use as an initialguess for the non-linear mode search. We then comparethe resulting non-linear mode with r AB by computing theoverlap c , shown in Fig. 8 for the linear mode used as ini-tial guess and for the corresponding non-linear mode. Weobserve in this case that linear and non-linear modes havecomparably poor predictive power. Of course, we cannotexclude that there is another DW starting from mini-mum A and connecting to another minimum B (cid:48) , whichmight be better correlated with these modes, althoughthis scenario seems unlikely.Finally, we investigated whether the energy profilealong the non-linear modes we found using this proce-dure, i.e. v ( s ) = V ( r A + sπ ), displays a DW shape, andwe did not find any DW in approximately 95% of cases.This is consistent with results reported in [65], and illus-trates the complexity of the energy landscape in whichTLS reside. V. CONCLUSIONS In this paper, we have investigated the relationship be-tween QLM and TLS in silico in a model glass, exploit-ing the TLS library constructed in [30]. We find that softQLM are generally associated with the initial direction ofthe minimum energy path connecting two minima, and asa consequence, DW potentials are spatially located closeto a soft QLM. However, the frequency of the QLM isanticorrelated with the tunnel splitting associated to theDW, hence TLS are typically not associated to the softestmodes, which on the contrary are expected to be respon-sible for plasticity [52, 66]. We conclude that QLM withproperly tuned frequency could serve as good predictorsof the spatial location of TLS. However, we also find thatthe minimal energy path is strongly curved within a high-dimensional space, in such a way that the saddle pointand the secondary minimum are uncorrelated with thedirection of the initial tangent vector. We find that lin-ear modes are poor predictors of the minimum-to-saddleor minimum-to-minimum directions.We have also considered non-linear cubic modes [37]and find that one of these modes is often well correlatedwith the minimum-to-minimum direction r AB . However,locating this individual mode is difficult: if the searchis initialized with r AB itself, convergence to the correctmode is facile. If, on the contrary, the search is initial-ized in a soft linear mode, convergence to the correctnon-linear mode does not occur. We conclude that inabsence of some prior information about the direction of r AB , it is difficult to predict the more global displacementfield associated with TLS via either linear or non-linearmodes.The problem of finding good structural predictors forTLS thus remains somewhat open. It is possible that bet-ter search strategies could exploit the information con-tained in linear or non-linear modes more efficiently. Ma- chine learning techniques [67, 68] might be able to exploitthis information (and perhaps additional structural infor-mation) to achieve better performance at contact predic-tion. Exploring this possibility is a clear direction forfuture work. ACKNOWLEDGMENTS We would like to thank Ludovic Berthier, Eran Bouch-binder, Wencheng Ji, Edan Lerner, Felix-Cosmin Mo-canu, Corrado Rainone, Camille Scalliet, PierfrancescoUrbani, Matthieu Wyart for useful discussions and Lu-dovic Berthier and Camille Scalliet for assistance and col-laboration at the early stages of this project.This project has received funding from the EuropeanResearch Council (ERC) under the European Union’sHorizon 2020 research and innovation programme (grantagreement n. 723955 - GlassUniversality), and it wassupported by a grant from the Simons Foundation( [1] R. C. Zeller and R. O. Pohl, Phys. Rev. B , 2029 (1971).[2] M. T. Loponen, R. C. Dynes, V. Narayanamurti, andJ. P. Garno, Physical Review B , 1161 (1982).[3] J. F. Berret and M. Meißner, Zeitschrift f¨ur Physik BCondensed Matter , 65 (1988).[4] A.-M. Boiron, P. Tamarat, B. Lounis, R. Brown, andM. Orrit, Chemical physics , 119 (1999).[5] A. L. Burin, J. M. Leveritt III, G. Ruyters, C. Sch¨otz,M. Bazrafshan, P. Fassl, M. von Schickfus, A. Fleis-chmann, and C. Enss, EPL (Europhysics Letters) ,57006 (2013).[6] D. R. Queen, X. Liu, J. Karel, T. H. Metcalf, and F. Hell-man, Phys. Rev. Lett. , 135901 (2013).[7] T. P´erez-Casta˜neda, R. J. Jim´enez-Riob´oo, and M. A.Ramos, Phys. Rev. Lett. , 165901 (2014).[8] T. P´erez-Casta˜neda, C. Rodr´ıguez-Tinoco, J. Rodr´ıguez-Viejo, and M. A. Ramos, Proceedings of the NationalAcademy of Sciences , 11275 (2014).[9] X. Liu, D. R. Queen, T. H. Metcalf, J. E. Karel, andF. Hellman, Phys. Rev. Lett. , 025503 (2014).[10] D. R. Queen, X. Liu, J. Karel, H. C. Jacks, T. H. Metcalf,and F. Hellman, Journal of Non-Crystalline Solids ,19 (2015).[11] P. W. Anderson, B. I. Halperin, and C. M. Varma, Phi-los. Mag. , 1 (1972).[12] W. A. Phillips, J. Low Temp. Phys. , 351 (1972).[13] W. A. Phillips, Reports on Progress in Physics , 1657(1987).[14] Y. M. Galperin, V. Karpov, and V. Kozub, Advances inPhysics , 669 (1989).[15] A. J. Leggett and D. C. Vural, The Journal of Physical Chemistry B , 12966 (2013).[16] D. Zhou, Journal of Physics: Condensed Matter ,055704 (2019).[17] H. M. Carruzzo and C. C. Yu, Phys. Rev. Lett. ,075902 (2020).[18] C. Artiaco, F. Balducci, and A. Scardicchio,arXiv:2010.03498 (2020).[19] A. Ninarello, L. Berthier, and D. Coslovich, PhysicalReview X , 021039 (2017).[20] F. H. Stillinger and T. A. Weber, Physical Review A ,978 (1982).[21] T. A. Weber and F. H. Stillinger, Physical Review B ,5402 (1985).[22] A. Heuer and R. J. Silbey, Phys. Rev. Lett. , 3911(1993).[23] D. Dab, A. Heuer, and R. Silbey, Journal of Lumines-cence , 95 (1995).[24] A. Heuer and R. J. Silbey, Phys. Rev. B , 609 (1996).[25] F. Demichelis, G. Viliani, and G. Ruocco, PhysChem-Comm , 20 (1999).[26] J. Reinisch and A. Heuer, Phys. Rev. B , 064201(2004).[27] J. Reinisch and A. Heuer, Phys. Rev. Lett. , 155502(2005).[28] A. Heuer, Journal of Physics: Condensed Matter ,373101 (2008).[29] T. Damart and D. Rodney, Phys. Rev. B , 014201(2018).[30] D. Khomenko, C. Scalliet, L. Berthier, D. R. Reichman,and F. Zamponi, Physical Review Letters , 225901(2020). [31] H. J´onsson, G. Mills, and K. W. Jacobsen, in Classi-cal and Quantum Dynamics in Condensed Phase Simu-lations , edited by B. J. Berne, G. Ciccotti, and D. F.Coker (World Scientific, Singapore, 1998).[32] G. Henkelman, B. P. Uberuaga, and H. J´onsson, TheJournal of chemical physics , 9901 (2000).[33] P. G. Bolhuis, D. Chandler, C. Dellago, and P. L.Geissler, Annual review of physical chemistry , 291(2002).[34] E. Vanden-Eijnden, Annual review of physical chemistry , 391 (2010).[35] B. B. Laird and H. R. Schober, Phys. Rev. Lett. , 636(1991).[36] M. Baity-Jesi, V. Mart´ın-Mayor, G. Parisi, and S. Perez-Gaviro, Phys. Rev. Lett. , 267205 (2015).[37] L. Gartner and E. Lerner, SciPost Phys. , 016 (2016).[38] E. Lerner, G. D¨uring, and E. Bouchbinder, Phys. Rev.Lett. , 035501 (2016).[39] H. Mizuno, H. Shiba, and A. Ikeda, Proc. Natl. Acad.Sci. U.S.A. , E9767 (2017).[40] L. Wang, A. Ninarello, P. Guan, L. Berthier, G. Szamel,and E. Flenner, Nature communications , 26 (2019).[41] C. Rainone, E. Bouchbinder, and E. Lerner, Proc. Natl.Acad. Sci. U.S.A. , 5228 (2020).[42] D. Richard, K. Gonz´alez-L´opez, G. Kapteijns, R. Pater,T. Vaknin, E. Bouchbinder, and E. Lerner, Phys. Rev.Lett. , 085502 (2020).[43] P. Das, H. G. E. Hentschel, E. Lerner, and I. Procaccia,Phys. Rev. B , 014202 (2020).[44] H. Mizuno, S. Mossa, and J.-L. Barrat, Phys. Rev. B , 144303 (2016).[45] S. Gelin, H. Tanaka, and A. Lemaˆıtre, Nature materials , 1177 (2016).[46] M. Paoluzzi, L. Angelani, G. Parisi, and G. Ruocco,Phys. Rev. Lett. , 155502 (2019).[47] A. Moriel, G. Kapteijns, C. Rainone, J. Zylberg,E. Lerner, and E. Bouchbinder, J. Chem. Phys. ,104503 (2019).[48] C. Oligschleger and H. R. Schober, Phys. Rev. B , 811(1999).[49] A. Widmer-Cooper, H. Perry, P. Harrowell, and D. R.Reichman, Nature Phys. , 711 (2008).[50] M. L. Manning and A. J. Liu, Phys. Rev. Lett. ,108302 (2011).[51] W. Ji, M. Popovi´c, T. W. J. de Geus, E. Lerner, andM. Wyart, Phys. Rev. E , 023003 (2019).[52] D. Richard and et al., arXiv:2003.11629 (2020).[53] M. Shimada, H. Mizuno, M. Wyart, and A. Ikeda, Phys.Rev. E , 060901 (2018).[54] M. Shimada, H. Mizuno, L. Berthier, and A. Ikeda, Phys.Rev. E , 052906 (2020).[55] V. Lubchenko and P. G. Wolynes, Phys. Rev. Lett. ,195901 (2001).[56] U. Buchenau, Y. M. Galperin, V. L. Gurevich, and H. R.Schober, Phys. Rev. B , 5039 (1991).[57] U. Buchenau, Y. M. Galperin, V. L. Gurevich, D. A.Parshin, M. A. Ramos, and H. R. Schober, Phys. Rev.B , 2798 (1992).[58] L. Gil, M. A. Ramos, A. Bringer, and U. Buchenau,Phys. Rev. Lett. , 182 (1993).[59] V. L. Gurevich, D. A. Parshin, and H. R. Schober, Phys.Rev. B , 094203 (2003).[60] P. Das, H. G. E. Hentschel, E. Lerner, and I. Procaccia,Phys. Rev. B , 014202 (2020). [61] R. K¨uhn and U. Horstmann, Phys. Rev. Lett. , 4067(1997).[62] C. Rainone, E. Bouchbinder, E. Lerner, P. Urbani, andF. Zamponi, arXiv:2010.11180 (2020).[63] S. Albert, G. Biroli, F. Ladieu, R. Tourbot, and P. Ur-bani, arXiv:2010.03294 (2020).[64] A. Q. Tool, Journal of the American Ceramic society ,240 (1946).[65] Edan Lerner, private communication. See the Lec-tures at the workshop Recent progress in glassy sys-tems , Les Houches, February 2020 for related results.Slides available at https://sites.google.com/view/leshouches2020 .[66] S. Karmakar, A. Lemaˆıtre, E. Lerner, and I. Procaccia,Phys. Rev. Lett. , 215502 (2010).[67] D. M. Sussman, S. S. Schoenholz, E. D. Cubuk, and A. J.Liu, Proceedings of the National Academy of Sciences , 10601 (2017).[68] V. Bapst, T. Keck, A. Grabska-Barwi´nska, C. Donner,E. D. Cubuk, S. S. Schoenholz, A. Obika, A. W. Nel-son, T. Back, D. Hassabis, et al. , Nature Physics16