Spatiotemporal intermittency and localized dynamic fluctuations upon approaching the glass transition
J. Ariel Rodriguez Fris, Eric R. Weeks, Francesco Sciortino, Gustavo A. Appignanesi
SSpatiotemporal intermittency and localized dynamic fluctuations upon approachingthe glass transition
J. Ariel Rodriguez Fris , Eric R. Weeks , Francesco Sciortino , , and Gustavo A. Appignanesi INQUISUR, Departamento de Qu´ımica,Universidad Nacional del Sur (UNS)-CONICET,Avenida Alem 1253, 8000 Bah´ıa Blanca, Argentina. Physics Department, Emory University,Atlanta, Georgia 30322, USA. Dipartimento di Fisica, Sapienza Universita’ di Roma,Piazzale A. Moro 5, Roma 00185, Italy. CNR-ISC, c/o Sapienza,Piazzale A. Moro 5, Roma 00185, Italy. (Dated: October 3, 2018)We introduce a new and robust method for characterizing spatially and temporally heterogeneousbehavior within a system based on the evolution of dynamic fuctuations once averaged over differentspace lengths and time scales. We apply it to investigate the dynamics in two canonical systemsas the glass transition is approached: a simulated Lennard-Jones glass-former and a real densecolloidal suspensions. We find that in both cases the onset of glassines is marked by spatiallylocalized dynamic fluctuations originating in regions of correlated mobile particles. By removingthe trivial system size dependence of the fluctuations we show that such regions contain tens tohundreds of particles for time scales corresponding to maximally non-Gaussian dynamics.
I. INTRODUCTION
Glasses are solid materials with disordered liquid-likestructure. These are typically formed by rapidly quench-ing a liquid from a hot to a cold temperature, or com-pressing a liquid from a low to a high pressure [1, 2].How the transition from an equilibrium liquid to anout of equilibrium glass takes places is highly debated,with many different and contrasting interpretations pro-posed; see Refs. [3–8] for reviews. One point is known:the onset of unusual behavior within a sample precedesthe glass transition. “Supercooled liquids”, despite theirmetastable equilibrium, have a markedly higher viscosity η than normal liquids. This step rise of η is associateto the onset of dynamical heterogeneity: diffusive mo-tion takes place in a spatially and temporally heteroge-neous fashion [9–13]. At any given time, some regionswithin the sample are frozen, while other regions arequite mobile. The mobile regions are characterized by“cooperative” motion where localized groups of moleculeshave nearly simultaneous large displacements [11, 14–16].Over time, the locations of faster and slower dynamicschange, such that at any given position the dynamics aretemporally heterogeneous as well [4, 17, 18].In the last two decades, a variety of methods have beenproposed and implemented to characterize such dynam-ical heterogeneities. Early work studied simulations ofsoft particles or Lennard-Jones particles and identifiedsubsets of particles that had large displacements ∆ r [14–16, 19, 20], showing that these formed spatially local-ized clusters. A key result was the identification of thenon-Gaussian time scale ∆ t ∗ as an important time scalerelated to these clusters [14, 20, 21]. This time scale isidentified by examining the behavior of the non-Gaussian parameter α (∆ t ) = 3 (cid:104) ∆ r (cid:105) / (5 (cid:104) ∆ r (cid:105) ) −
1, a quantityderived from the moments of the displacement distribu-tion P (∆ r ) [22]. ∆ t ∗ is the time scale for which α (∆ t )is maximal. α is zero when P (∆ r ) is a Gaussian, whichholds to a good approximation for simple liquids. Fordynamically heterogeneous supercooled liquids, α > α value is expected to provide information on the co-operatively moving clusters [14–16] and their structuraland dynamic properties. Early studies [14–16] used var-ious somewhat arbitrary criteria to define mobile parti-cles. Later work examined spatial correlation functionsaveraged over all particles in various ways attempting toidentify the length and time scales of dynamical hetero-geneity [11, 16, 23–29].In this manuscript, we present a new and robust anal-ysis method to characterize spatial and temporal dynam-ical heterogeneity that does not requires any a priori def-inition of particle mobility. In particular, we here use thesystem mean square displacement as a “null hypothesis”for particle motion, and quantify spatially and tempo-rally localized deviations of particle motion away fromthis null hypothesis. We apply this method to the Kob-Andersen Lennard-Jones glassforming system [30] and tocolloidal supercooled liquid data [31]. Our results showthat dynamical heterogeneity is most obvious for subsys-tems comprised of tens to hundreds of particles, with thesize growing as the glass transition is approached. Ad-ditionally, we examine how dynamical heterogeneity be-comes averaged out at larger length scales. As a byprod-uct we confirm that ∆ t ∗ is the time scale of maximumheterogeneity. While our method of localized fluctuations a r X i v : . [ c ond - m a t . d i s - nn ] J a n is applied to particle displacements, the idea is generaliz-able to other quantities which may have spatiotemporalfluctuations such as structure [32, 33]. An advantage ofour technique is that it is applicable to small data setssuch as the experimental colloidal data we use. II. SIMULATION AND EXPERIMENTALDETAILSA. Simulation
We use LAMMPS to simulate the Kob-Andersen bi-nary Lennard-Jones glassforming system [30]. Briefly,this is an 80:20 mixture of A and B particles. The par-ticles interact via the Lennard-Jones potential [34] U αβ ( r ) = 4 (cid:15) αβ (cid:20)(cid:16) σ αβ r (cid:17) − (cid:16) σ αβ r (cid:17) (cid:21) (1)with α, β ∈ A, B . A and B particles have the samemass. The energy scales are (cid:15) AA = 1 . (cid:15) AB = 1 .
5, and (cid:15)BB = 0 .
5. The size scales are σ AA = 1 . σ AB = 0 . σ BB = 0 .
88, chosen so that A and B particles areencouraged to mix rather than segregate, and thus crys-tallization is frustrated [30]. For most of our results wepresent data with 8000 total particles, and for our analy-sis we consider the N = 6400 A particles. To verify thatour analysis is not biased by finite size, in a few cases wecompare with a N = 8 × data set. Periodic boundaryconditions were used with a cubical box. B. Experiment
Colloids have long been used as model systems tostudy the glass transition [35–43]. We reanalyze pre-viously published data from experiments using confocalmicroscopy to observe dense colloidal samples [31]. Thesamples were sterically stabilized colloidal poly-(methylmethacrylate) for which the key control parameter is thevolume fraction φ [35]. The glass transition for this ex-periment occurred at φ g = 0 .
58. Here we examine datawith φ = 0 . , . , .
56 with N ≈ , , µ m and a polydispersity of 0.045 [44] and wereslightly charged. Confocal microscopy and particle track-ing was used to follow the positions of the particles inthree dimensions [45, 46]. The imaging volume was rect-angular with an aspect ratio roughly 5 : 5 : 1; see Ref. [31]for details. III. RESULTS
We aim at characterizing the growing spatial and tem-poral fluctuations on approaching the glass transition without introducing any arbitrary cut-off quantity. Fol-lowing prior work [28, 29, 47], we start by defining thedistance matrix ∆ ( t (cid:48) , t (cid:48)(cid:48) ), an object that represents theaverage of the squared particle displacements betweentime t (cid:48) and t (cid:48)(cid:48) of a collection of the N particles belong-ing to a predefined set S ( S may be the entire system orsome subsystem):∆ ( t (cid:48) , t (cid:48)(cid:48) ) ≡ N N (cid:88) i =1 | (cid:126)r i ( t (cid:48) ) − (cid:126)r i ( t (cid:48)(cid:48) ) | (2)= (cid:104)| (cid:126)r i ( t (cid:48) ) − (cid:126)r i ( t (cid:48)(cid:48) ) | (cid:105) i ∈ S (3)where the angle brackets indicate an average over the N particles in S . Further averaging ∆ ( t (cid:48) , t (cid:48)(cid:48) ) over all pairs t (cid:48) and t (cid:48)(cid:48) such that t (cid:48)(cid:48) − t (cid:48) = ∆ t produces the well-knownaverage mean square displacement M (∆ t ) of particlesin S . More precisely M (∆ t ) = (cid:104) ∆ ( t (cid:48) , t (cid:48)(cid:48) ) (cid:105) t (cid:48)(cid:48) − t (cid:48) =∆ t (4)where the average is over t (cid:48) , t (cid:48)(cid:48) with fixed time inter-val ∆ t = | t (cid:48)(cid:48) − t (cid:48) | and also over all of the particlesin S . Assuming stationary dynamics (true as long asthe system is not aging), for a sufficiently large ∆ t ,lim ∆ t →∞ ∆ ( t (cid:48) , t (cid:48) + ∆ t ) = M (∆ t )For small systems under glassy relaxation conditions,∆ has temporal fluctuations, as shown in Fig. 1(a) fora sub-system of N = 125 particles. Darker regions indi-cate time intervals ( t (cid:48) , t (cid:48)(cid:48) ) over which this subsystem hasrelatively little particle motion. Clearly, there are spe-cific times for which this subsystem undergoes fairly largechanges, signalled by larger displacements and substan-tially different particle positions. It is expected [48] thanon increasing the sample size well beyond any dynamiccorrelation length, different regions of the system will in-dependently display such burst motion, such that ∆ forlarger systems will appear much smoother, as shown inFig. 1(b) for the fully system of N = 8000 particles. Asfor the long time limit, for a sufficiently large system S ,lim N →∞ ∆ ( t (cid:48) , t (cid:48) + ∆ t ) = M (∆ t )The question we turn to is how the large system limitis reached. In particular, we wish to use the approachto the large system limit to characterize the spatial scaleof dynamical heterogeneities. The obvious features ofFig. 1(a) are the large fluctuations that differentiate itfrom Fig. 1(b). This motivates us to consider the nor-malized difference between ∆ and the expectation for alarge system, defined byΩ S ( t (cid:48) , t (cid:48)(cid:48) ) = [∆ ( t (cid:48) , t (cid:48)(cid:48) ) − M (∆ t )] [ M (∆ t )] (5)with the convention ∆ t = | t (cid:48)(cid:48) − t (cid:48) | . Ω S , a measure ofthe dynamic intermittency, represents the matrix of nor-malized squared deviations from the mean value for theparticles squared displacements and will be equal to zerowhen ∆ is calculated for sufficiently large systems, forwhich time averages and space averages are equivalentFIG. 1: (a) Contour plot of the distance matrix ∆ ( t (cid:48) , t (cid:48)(cid:48) ) for a binary Lennard-Jones system at T = 0 .
50 within acubical subsystem containing 125 particles. The legend indicates the values of the gray scale. (b) Contour plot ofthe same system for the full 8000 particle simulation. Images taken from Ref. [29], permission pending. (c) Ω( t (cid:48) , t (cid:48)(cid:48) ).The data correspond to the same subsystem as in (a). In the image shown, darker points indicate larger values. (d)The values of ∆ as a function of ∆ t = | t (cid:48) − t (cid:48)(cid:48) | for the data shown in (a). For small ∆ t this tends to 0, and forintermediate ∆ t the scatter indicates the temporal heterogeneity in this subsystem. The overall increase mirrors themean square displacement M (∆ t ). (e) Ω (∆ t ) for the data shown in (c). The scatter at small ∆ t reflects a functionof ∆ t , which has much larger fluctuations at short ∆ t , reflecting non-Gaussian statistics of the displacements atthose time scales.and ∆ = M . Otherwise, Ω S > (local in both spaceand time) and the expectation for a large system (thatis, M , a quantity averaged over all space and all time).An example of Ω S ( t (cid:48) , t (cid:48)(cid:48) ) is shown in Fig. 1(c), wherethe darker regions indicate time periods for which themean motion within the sub-volume is anomalouslylarger or smaller compared to the expectation from M .Fig. 1(d) and 1(e) show scatter plots of the values of∆ and Ω S as functions of t (cid:48)(cid:48) − t (cid:48) , taken from the dataof Fig. 1(a) and 1(c) respectively. ∆ starts at the ori-gin and rises, consistent with the idea that on averageit should behave similar to the mean square displace-ment M (∆ t ). At intermediate time scales the scatterin the data of Fig. 1(d) indicates the temporal fluctua-tions of the motion. In contrast, Ω S in Fig. 1(e) has largefluctuations at the shortest time scales, indicating largefluctuations of the motion relative to M (∆ t ) on thosetime scales. At larger time scales, the temporal averagingreduces Ω S toward zero.As defined, Ω S is local in space and time. To focuson the space dependence of the fluctuations we need toevaluate a single nondimensional scalar quantity Ω( N )characterizing the mobility fluctuations for subsystemsof size N (the ratio of the dispersion to the average [49]for the particle squared displacements). To do so we par-tition the system into distinct cubical boxes containing N particles each and evaluate the sum of Ω S ( t (cid:48) , t (cid:48)(cid:48) ) over alltime pairs ( t (cid:48) , t (cid:48)(cid:48) ) [i.e. the sum over all points entering inthe scatter-plot as the one shown in Fig. 1(e)] divided bythe number of such pairs for each of the boxes. We thenaverage the resulting number over all boxes and finallytake the square root of the result. This procedure yieldsthe desired scalar quantity Ω( N ). Note that the specificvalues of Ω( N ) will depend on the total time studied, that is, the maximum of | t (cid:48)(cid:48) − t (cid:48) | that is studied. As is appar-ent from Fig. 1(c) and 1(e), at large | t (cid:48)(cid:48) − t (cid:48) | , Ω S decaysto zero, and the more of this included in the average, thesmaller Ω( N ) will be. However, for a given data set, whatwill matter is the N -dependence which is insensitive tothe total time studied, as long as that time is sufficientto capture the temporal fluctuations seen in Fig. 1(c). Inpractice, we ensure that our data sets have a total dura-tion of ≈ t ∗ where ∆ t ∗ depends on the temperature(for the Lennard-Jones simulations) or the volume frac-tion (for the colloidal experiments). This then will allowfor a sensible comparison between different data sets.Ω( N ) for the Lennard-Jones system is plotted inFig. 2(a), and for the colloidal experiments in Fig. 2(b).In both cases, we see how the dynamical fluctuations av-erage out for larger subsystem sizes. Notably, the systemscloser to the glass transition require larger subsystems be-fore the dynamical fluctuations are averaged out – coldersystems for the LJ data (a), and higher volume fractionsystems for the colloidal data (b).An interesting point emerges when analyzing the func-tional form of the decay of Ω( N ) with N . From thehigh temperature data of Fig. 2(a) we observe the spa-tially localized dynamic fluctuations measured by Ω dis-play the usual N − / size scaling dependence. This re-flects that particle motion is nearly spatially uncorrelatedwithin a subsystem and so the average of ∆ ( N ) con-verges to the large-system limit M as N − / . However,for glassier systems a clear departure from this trivialbehavior is seen, and the decay of Ω( N ) gets slower.The fact that the localized dynamical fluctuations arehigher than expected, persisting at large system sizes,speaks of the existence of regions of correlated mobileparticles, an effect that is more pronounced upon super-cooling [11, 14–16, 25, 31]. We truncate our calculationFIG. 2: Ω(N) and Ω R (N) as a function of subsystem size N for the binary Lennard-Jones system [(a) and (c)] fordifferent temperatures as indicated and for the experimental colloidal suspension [(b) and (d)] at different volumefractions φ as indicated. In the case of Ω(N) the N particles are part of the same subsystem. In the case of Ω R (N)the N particles are selected randomly among all particles in the system.FIG. 3: Subtracting Ω calculated for randomly chosenparticles from Ω calculated for compact subsystems. (a)For the LJ system; (b) for the colloidal suspensions.The black curve in (a) shows the same quantityevaluate for a system of 1M particles at T = 0 . N = N max /
10, as we desire at least 10 subsystems toevaluate a reasonable Ω( N ). We note that both for theLennard-Jones system and the colloidal suspensions, thebehavior of the curves of Fig. 2 at very large system sizesis due to lack of statistics (the subsystems are not smallas compared to the large system and, thus, we can onlyaverage over a few of them).As evident from Fig. 2, for the smallest possible sub-systems ( N = 1), the increasing value of Ω as the glasstransition is approached reflects the well-known increas-ing non-Gaussian nature of the displacement distribu-tion [30]. In fact this points out a limitation of Ω, in thatlarge values of Ω can reflect either spatial fluctuations inthe dynamics or simply a non-Gaussian distribution ofdisplacements. To remove the latter influence (that is,to remove the trivial system size dependence and thus tohighlight the local correlations) we separately computeΩ R ( N ) based not on compact subsystems of size N buton N randomly chosen particles. Here the subscript R indicates an average over many such randomly chosensubsets. For the Lennard-Jones system, this is plottedin Fig. 2(c), showing different behavior from Fig. 2(a)for the colder temperature data. In fact, now the ran- domly distributed dynamical fluctuations quantified byΩ R ( S ) display the typical N − / decay at all temper-atures. Similar behavior is found in Fig. 2(d) for thecolloidal suspensions, as compared to Fig. 2(b). Again,the N − / scaling is recovered at all volume fractions.To understand the differences between the spatiallylocalized and the randomly distributed dynamic fluctu-ations, in Fig. 3(a) we show the result of subtractingΩ R ( N ) from Ω( N ). For N = 1 particles the result is zeroas there is no distinction between the two calculations.Likewise for N → A particles),the two calculations are identical. At intermediate num-bers of particles, nonzero values are found in Fig. 3(a)indicating nontrivial spatially localized values of Ω. Inparticular, for colder temperatures the dynamical fluctu-ations are larger [higher curves in Fig. 3(a)] and the sub-system size with the largest fluctuations grows slightly[peak position shifts rightward in Fig. 3(a)]. This lastobservation is quantified in Fig. 4 which shows the peakposition of Fig. 3(a) as a function of T . The peak oc-curs for larger subsystems at colder temperatures. Thisindicates the size of regions with maximally variable dy-namics contain about 50 particles for the coldest samples.For hotter samples, the maximum shifts to smaller sys-tem sizes; for the hottest data, Ω for compact subsystemsis nearly indistinguishable from Ω calculated for randomparticles, and we cannot identify a maximum. In turn,Fig. 3(b) shows the behavior for the colloidal suspensions,which is similar to that found for the Lennard-Jones sys-tem. Namely, there is a peak in Fig. 3(b) at a specificsubsystem size, and the position and height of this peakis larger for volume fractions closer to the glass transitionvolume fraction ( φ g ≈ .
58 [31]). An intriguing differenceis that the peak position indicates that larger subsystemscontaining a few hundred particles are maximally hetero-geneous, as compared to Fig. 3(a) which peaks at sub-system sizes containing a few tens of particles. Earlierwork with the Lennard-Jones system found mobile clus-ters containing 10-30 particles [50] or 40 particles [28]at the coldest temperatures, in agreement with our re-sult. Likewise, earlier analyses of the same colloidal datafound the largest mobile clusters contained ∼
50 particleson average for φ = 0 .
562 [31]. The new result seen here ishow the sample behaves over larger length scales. For ex-ample, Fig. 3(a) shows that there is still nontrivial spatialheterogeneity for subsystems containing N = 1000 parti-cles, more than an order of magnitude larger than the N corresponding to the peak. This is strong evidence thatthe dynamically heterogeneous regions of size N ∼
50 arenot randomly distributed throughout the sample but arethemselves spatially clustered.FIG. 4: (a) Value of N that maximizes Ω( N ) inFig. 3(a) as a function of T . (b) Comparison of the timescale of the maxima from Fig. 3(a) with the time scale∆ t ∗ .In the preceding analysis, Ω S is averaged over timescales to highlight the N dependence of the dynamicalfluctuations. We now turn to the complementary case,and average Ω S over subsystem sizes N to find the timescale ∆ t of the same fluctuations. Since the value of N covers several order of magnitude, we average Ω S ( t (cid:48) , t (cid:48)(cid:48) )over subsystems with sizes N picked to be evenly dis-tributed in log( N ), with N ranging from 1 to the systemsize. We first average over all subsystems of a given size N , and then average Ω ( N, t (cid:48) , t (cid:48)(cid:48) ) over N and ( t (cid:48) , t (cid:48)(cid:48) ) with fixed time interval ∆ t = t (cid:48)(cid:48) − t (cid:48) to result in Ω(∆ t ). Thisis shown in Fig. 4(b) where the peak of each curve ismarked with a red circle. For comparison, the time scale∆ t ∗ is indicated for each data set with a yellow circle;it appears the peak of Ω(∆ t ) is always close to the non-Gaussian time scale ∆ t ∗ . As expected, the time scaleof maximum dynamical heterogeneity [as measured byΩ(∆ t )] grows as the system approaches the glass transi-tion. IV. CONCLUSIONS
In summary, we have constructed a new measure ofspatial and temporal dynamic heterogeneity. The mea-sure does not require defining subsets of mobile or im-mobile particles, but rather looks for fluctuations awayfrom the large system behavior. Additionally, it allowsus to examine dynamical heterogeneity on a variety oflength scales, showing that the approach to the large sys-tem limit is slower than would be expected for randomlydistributed fluctuations in the dynamics. The methodcan be straightforwardly applied to experimental systemssuch as the dense colloidal solution we examine; it doesnot require finite-size scaling, for example. While we fo-cus on particle motion where the mean square displace-ment is the null hypothesis, the method can be gener-alized to any other spatially and temporally fluctuatingquantities, as long as there is a well-defined null hypoth-esis based on the large system limit.The work of E.R.W was supported by a grant from theNational Science Foundation (DMR-1609763). GAA andJARF acknowledge suport form CONICET, UNS andANPCyT(PICT2015/1893). Inspiration for this workcame from past conversations with Prof. Walter Kob. [1] C. A. Angell, Science , 1924 (1995).[2] C. A. Angell, K. L. Ngai, G. B. McKenna, P. F. McMil-lan, and S. W. Martin, J. App. Phys. , 3113 (2000).[3] J. S. Langer, Rep. Prog. Phys. , 042501 (2014).[4] D. Chandler and J. P. Garrahan, Ann. Rev. Phys. Chem. , 191 (2010).[5] G. Biroli and J. P. Garrahan, J. Chem. Phys. ,12A301 (2013).[6] M. D. Ediger and P. Harrowell, J. Chem. Phys. ,080901 (2012).[7] V. Lubchenko and P. G. Wolynes, Ann. Rev. Phys. Chem. , 235 (2007).[8] A. Cavagna, Phys. Rep. , 51 (2009).[9] H. Sillescu, J. Non-Cryst. Solids , 81 (1999).[10] M. D. Ediger, Ann. Rev. Phys. Chem. , 99 (2000).[11] S. C. Glotzer, Physics of Non-Crystalline Solids 9 , J.Non-Cryst. Solids , 342 (2000).[12] E. Hempel, G. Hempel, A. Hensel, C. Schick, andE. Donth, J. Phys. Chem. B , 2460 (2000).[13] R. Richert, J. Phys.: Condens. Matter , R703 (2002).[14] W. Kob, C. Donati, S. J. Plimpton, P. H. Poole, and S. C. Glotzer, Phys. Rev. Lett. , 2827 (1997).[15] C. Donati, J. F. Douglas, W. Kob, S. J. Plimpton, P. H.Poole, and S. C. Glotzer, Phys. Rev. Lett. , 2338(1998).[16] C. Donati, S. C. Glotzer, and P. H. Poole, Phys. Rev.Lett. , 5064 (1999).[17] J. P. Garrahan and D. Chandler, Proc. Nat. Acad. Sci. , 9710 (2003).[18] A. S. Keys, L. O. Hedges, J. P. Garrahan, S. C. Glotzer,and D. Chandler, Phys. Rev. X , 021013 (2011).[19] M. M. Hurley and P. Harrowell, Phys. Rev. E , 1694(1995).[20] M. M. Hurley and P. Harrowell, J. Chem. Phys. ,10521 (1996).[21] A. H. Marcus, J. Schofield, and S. A. Rice, Phys. Rev.E , 5725 (1999).[22] A. Rahman, Phys. Rev. , A405 (1964).[23] E. Flenner and G. Szamel, J. Phys.: Condens. Matter , 205125 (2007).[24] B. Doliwa and A. Heuer, Phys. Rev. E , 6898 (2000).[25] E. R. Weeks, J. C. Crocker, and D. A. Weitz, J. Phys.: Condens. Matter , 205131 (2007).[26] N. Laˇcevi´c, F. W. Starr, T. B. Schrøder, and S. C.Glotzer, J. Chem. Phys. , 7372 (2003).[27] A. S. Keys, A. R. Abate, S. C. Glotzer, and D. J. Durian,Nature Phys. , 260 (2007).[28] G. A. Appignanesi, J. A. Rodriguez Fris, R. A. Montani,and W. Kob, Phys. Rev. Lett. , 057801 (2006).[29] G. A. Appignanesi and J. A. Rodriguez Fris, J. Phys.:Condens. Matter , 203103 (2009).[30] W. Kob and H. C. Andersen, Phys. Rev. E , 4626(1995).[31] E. R. Weeks, J. C. Crocker, A. C. Levitt, A. Schofield,and D. A. Weitz, Science , 627 (2000).[32] C. P. Royall and S. R. Williams, Phys. Rep. , 1(2015).[33] C. P. Royall, S. R. Williams, T. Ohtsuka, and H. Tanaka,Nat Mater , 556 (2008).[34] J. E. Lennard-Jones, Proc. Roy. Soc. London A: Math-ematical, Physical and Engineering Sciences , 463(1924).[35] G. L. Hunter and E. R. Weeks, Rep. Prog. Phys. ,066501 (2012).[36] L. Marshall and C. F. Zukoski, J. Phys. Chem. , 1164(1990).[37] P. N. Pusey and W. van Megen, Phys. Rev. Lett. , 2083 (1987).[38] P. N. Pusey and W. van Megen, Nature , 340 (1986).[39] W. H¨artl, Curr. Op. Coll. Int. Sci. , 479 (2001).[40] W. van Megen and P. N. Pusey, Phys. Rev. A , 5429(1991).[41] E. Bartsch, V. Frenz, S. Moller, and H. Sillescu, PhysicaA , 363 (1993).[42] A. Meller and J. Stavans, Phys. Rev. Lett. , 3646(1992).[43] H. M. Lindsay and P. M. Chaikin, J. Chem. Phys. ,3774 (1982).[44] R. Kurita, D. B. Ruffner, and E. R. Weeks, NatureComm. , 1127 (2012).[45] A. D. Dinsmore, E. R. Weeks, V. Prasad, A. C. Levitt,and D. A. Weitz, App. Optics , 4152 (2001).[46] J. C. Crocker and D. G. Grier, J. Colloid Interface Sci. , 298 (1996).[47] I. Ohmine and H. Tanaka, Chem. Reviews , 2545(1993).[48] E. La Nave and F. Sciortino, J. Phys. Chem. B ,19663 (2004).[49] D. Chandler, Introduction to Modern Statistical Mechan-ics, by David Chandler, pp. 288. Oxford University Press,Sep 1987. , 288 (1987).[50] C. Donati, S. C. Glotzer, P. H. Poole, W. Kob, and S. J.Plimpton, Phys. Rev. E60