A novel approach to numerical measurements of the configurational entropy in supercooled liquids
The configurational entropy is among the key observables to characterize experimentally the formation of a glass. Physically, it quantifies the multiplicity of metastable states in which an amorphous material can be found at a given temperature, and its temperature dependence provides a major thermodynamic signature of the glass transition, which is experimentally accessible. Measurements of the configurational entropy require, however, some approximations which have often led to ambiguities and contradictory results. Here we implement a novel numerical scheme to measure the configurational entropy in supercooled liquids, using a direct determination of the free energy cost to localize the system within a single metastable state at temperature . For two prototypical glass-forming liquids, we find that disappears discontinuously above a temperature , which is slightly lower than the usual estimate of the onset temperature for glassy dynamics. This observation is in good agreement with theoretical expectations, but contrasts sharply with alternative numerical methods. While the temperature dependence of correlates with the glass fragility, we show that the validity of the Adam-Gibbs relation (relating configurational entropy to structural relaxation time) established in earlier numerical studies is smaller than previously thought, potentially resolving an important conflict between experiments and simulations.
pacs:05.10.-a, 05.20.Jj, 64.70.Q-
The configurational entropy (or complexity) plays an important role in descriptions of the glass transition because it quantifies the temperature evolution of the free energy landscape accompanying changes in thermodynamic and dynamic properties of supercooled liquids. It represents both a major experimental signature of the glass transition ediger_supercooled_1996 () and a fundamental quantity within a number of theoretical approaches berthier_theoretical_2011 ().
The configurational entropy is traditionally measured by subtracting a ‘vibrational’ contribution to the total entropy of the system: . While is well-defined, the vibrational contribution requires some approximation. Experiments ediger_supercooled_1996 (); kauzmann_nature_1948 (); richert_dynamics_1998 (); martinez_thermodynamic_2001 (); angell_specific_2002 () use for instance the entropy of the crystalline or glass states to estimate . In simulations, the above decomposition relies on the assumption that the system vibrates around a given ‘state’, further assumed to be equivalent to a local energy minimum, or inherent structure stillinger_hidden_1982 (). A thermodynamic formalism was developed to determine numerically the configurational entropy, and applied to a large number of models sciortino_inherent_1999 (); Scala_nature_2000 (); Sastry_2001 (); mossa_dynamics_2002 (); angelani_configurational_2007 (). These studies have additionally revealed that the Adam-Gibbs relation adam_temperature_1965 ()
between and the structural relaxation time is obeyed over a broad temperature window. In Eq. (1), is an energy scale and a microscopic timescale. Equation (1) is an important relation for supercooled liquids, as its validity would directly establish that the viscosity increase near the glass transition is caused by the temperature evolution of a complex free energy landscape.
Available numerical methods are however not fully satisfactory from both theoretical and experimental viewpoints. Firstly, the identification of metastable states with energy minima within the inherent structure formalism has been questioned berthier_theoretical_2011 (); biroli_inherent_2000 (). Because energy minima exist at all , the inherent structure exists at arbitrarily high temperatures, where the free energy landscape is in fact featureless. In theoretical approaches berthier_theoretical_2011 (); Wolynes_1997 (), is the entropic contribution stemming from the multiplicity of metastable states proliferating at low . While this definition is also plagued by ambiguities, as discussed below, specific calculations show that appears discontinuously below a temperature corresponding (within mean-field approximations) to the mode-coupling transition temperature berthier_theoretical_2011 (); cavagna_supercooled_2009 (). Secondly, the Adam-Gibbs relation in Eq. (1) was numerically found to be valid over the entire supercooled regime Scala_nature_2000 (); Sastry_2001 (); mossa_dynamics_2002 (); sengupta_dependence_2011 (). Experiments report instead that it only holds at low temperatures below the mode-coupling temperature richert_dynamics_1998 (), in a regime not accessible in simulations. These experimental findings are physically sensible because it is only at such low temperatures that the free energy landscape can possibly control the dynamics, but they directly contradict simulations.
We propose and implement a novel numerical method to measure the configurational entropy, which fully resolves these issues. The proposed methodology does not require precise definitions of a free energy landscape and metastable states. Our results show that the configurational entropy appears discontinuously at a characteristic low temperature, and that the Adam-Gibbs relation is not valid above the mode-coupling temperature. Therefore, this alternative approach provides a numerical estimate of the configurational entropy that is conceptually closer to theory, and yields quantitative results which agree better with experiments.
The proposed numerical method is directly inspired by statistical mechanics approaches, where the configurational entropy can be computed from the thermodynamic properties of constrained cloned systems berthier_theoretical_2011 (); cavagna_supercooled_2009 (); monasson_structural_1995 (); mezard_thermodynamics_1999 (); mezard__???? (). The physical idea is that constraining a system to reside ‘close’ to a single state has a free energy cost , because it represents the entropic loss due to an incomplete exploration of the configurational space. To bypass the difficulty of defining metastable states rigorously, we obtain a numerical estimate of by measuring a free energy difference between two thermodynamic phases that can be well-defined. In practice, we estimate from the thermodynamic properties of a system comprising two copies, 1 and 2, of the considered liquid thermalized at temperature . As described in more detail in the Supplementary Information (SI), we conduct equilibrium simulations of these two coupled copies and carefully measure the probability distribution function of their mutual overlap, , where brackets indicate an equilibrium average. We define the overlap as , where is the Heaviside function, denotes the position of particle within configuration 1, is a length comparable to the particle diameter (we take ), and the particle number in each copy. Note that this ‘collective’ overlap is insensitive to particle exchanges. We define the ‘effective potential’ , which is by definition the constrained equilibrium free energy of the total system when the average value of the overlap is franz_phase_1997 (); cardenas_constrained_1999 ().
While was introduced long ago in theoretical calculations franz_phase_1997 (), it was only recently realized that it can be accurately determined in computer simulations by applying tools first devised to study equilibrium phase transitions cammarota_phase-separation_2010 (); berthier_overlap_2013 (); parisi_liquid-glass_2013 (). For a particular model liquid, we have shown berthier_overlap_2013 () that is convex above a critical temperature below which it develops a linear part, corresponding to a strongly non-Gaussian . This observation implies that a thermodynamic field conjugated to the overlap induces, for , an equilibrium first-order transition between a low- and a high- phase franz_phase_1997 (). This first-order transition line ends at a second-order critical point at , as explicitly demonstrated in berthier_overlap_2013 (); parisi_liquid-glass_2013 (). The existence of two phases below suggests to estimate the configurational entropy as:
where denotes the position of the global minimum of , and is determined from the position of the peak in at coexistence, see Fig. 1. Equation (2) states that represents the free energy difference between the low- phase where the two copies independently explore the free energy landscape and the high- phase where they remain close to one another. This free energy difference originates from the fact that one of the copies cannot freely explore the configuration space, and this precisely costs . (Details pertaining to quenched and annealed complexities are discussed below.) While the complexity in Eq. (2) emerges naturally in mean-field calculations franz_phase_1997 (), our work is the first to implement this approach to estimate the configurational entropy in finite dimensional liquids.
The definition (2) shows that the measurement of does not rely on an explicit definition of a free energy landscape and of metastable states, and does not stem in the present approach from an enumeration of states. Instead, by measuring the thermodynamic properties of the high- localized phase, we let the system itself define the extent of a ‘state’. This provides a direct determination of which requires neither an approximate estimate of a vibrational contribution, nor a detailed investigation of the potential energy landscape. This approach, which relies on the direct measurement of a free energy difference, is conceptually much closer to theoretical calculations.
Another consequence of Eq. (2) is that is only defined when two distinct phases can be distinguished, i.e. for . For , is featureless and does not exist, see Fig. 1. In this regime, the entropy cannot be decomposed in configurational and vibrational parts. This is qualitatively consistent with specific theoretical calculations Wolynes_1997 (); monasson_structural_1995 (); mezard_thermodynamics_1999 (); franz_phase_1997 (). Physically, it means that the free energy landscape of the high temperature liquid has a simple topography for which the concept of configurational entropy is not relevant. A discontinuous emergence of the configurational entropy at low is naturally obtained within the present calculations, whereas it is missed by previous methods sciortino_inherent_1999 (); Sastry_2001 (); angelani_configurational_2007 ().
We studied two models of glass-formers using Monte-Carlo simulations berthier_monte_2007 (). The first model is a 50:50 binary mixture of harmonic spheres of diameter ratio 1.4 berthier_compressing_2009 (); berthier_glass_2009 (). Within reduced units kob_non-monotonic_2012 (), this quasi-hard sphere system has an onset temperature flenner_dynamic_2013 (), a mode-coupling temperature kob_non-monotonic_2012 (), and a Vogel-Fulcher temperature flenner_dynamic_2013 () (obtained with low reliability as the system is weakly fragile at this density berthier_compressing_2009 (); berthier_glass_2009 ()). We used , 108 and 256, finding that finite size effects for are small (see SI). We show data for . The second model is a 80:20 binary mixture of Lennard-Jones particles kob_testing_1995 (). In reduced units, the onset temperature is , the mode-coupling temperature kob_testing_1995 (), and the Vogel-Fulcher temperature Sastry_2001 () (the model has intermediate fragility). We performed simulations with . As described below (see Simulation Methods and SI for detailed descriptions), we combine umbrella sampling, multi-histogram reweighting and replica exachange techniques to quantity the rare fluctuations of the global overlap that need to be studied to obtain the free energy . We find that differences between various possible estimates of MPbook () can only be distinguished in a very narrow temperature regime near which is not resolved by the present set of data, and therefore does not affect any of our conclusions.
Our central results are in Fig. 2 which displays the temperature dependence of obtained from Eq. (2) for two glass models. In both cases, we find that emerges discontinuously at a critical temperature . We obtain for harmonic spheres berthier_overlap_2013 (), and for the Lennard-Jones model. Because is very close to, or slightly below, the onset temperature , this suggests that might represent a well-defined, physically meaningful definition of the onset temperature in supercooled liquids sastry_signatures_1998 (). Note that is significantly larger than obtained from a mode-coupling analysis of the dynamics. While and are found to coincide in mean-field calculations franz_phase_1997 (), we find that remains well-defined in finite dimensions, whereas the mode-coupling singularity is replaced by a smooth crossover.
The abrupt emergence of at stands in sharp contrast with alternative methods angelani_configurational_2007 (); sciortino_inherent_1999 (), as demonstrated in Fig. 2. Therefore, the qualitative evolution of obtained in this work is in closer agreement with theoretical and physical expectations (see e.g. Ref. berthier_theoretical_2011 ()). Notice that such a discontinuous temperature dependence is of course not observed experimentally, because experimental methods (just as previous numerical methods) are not sensitive to the sharp emergence of metastable states that we are able to reveal here. Physically, our results simply suggest that a decomposition of the entropy in vibrational and configurational parts is not meaningful at high temperatures, a fact which is also hinted by the inherent structure approach sastry_signatures_1998 ().
The configurational entropy in Fig. 2 decreases steadily as temperature is lowered below . This implies that the free energy difference between localized and delocalized states in configuration space decreases as gets lower, suggesting that the thermodynamic driving force to structural relaxation also decreases. A quantitative comparison with literature data in Fig. 2 shows that the temperature evolution of below is in qualitative agreement with earlier work. However, the inherent structure formalism provides an estimate of the configurational entropy that is systematically larger than over the explored range. Despite the shortcomings mentioned above, inherent structure based approaches might still represent a valuable approximation at low .
A motivation to determine follows from Kauzmann’s study of experimentally determined suggesting the existence of an entropy crisis, , possibly close to the Vogel-Fulcher temperature kauzmann_nature_1948 (). Our data do not cover a broad enough temperature range to extrapolate an entropy crisis. However, they do support the qualitative connection between thermodynamic and dynamic fragilities found experimentally martinez_thermodynamic_2001 (), since the more fragile Lennard-Jones system also has a steeper dependence of , as implied by Fig. 2.
The Adam-Gibbs relation in Eq. (1) is a quantitative connection between thermodynamics and dynamics that can readily be tested once is known, see Fig. 3. Notice first that, by construction, this relation cannot hold above the critical temperature where is not defined. Therefore, Eq. (1) cannot be expected to work if is too large. In fact, we find that it does not work well except close to , although we would need more data to establish more firmly its validity at lower temperatures. Therefore, our results indicate that the broad range of validity of Eq. (1) reported in earlier simulations Scala_nature_2000 (); Sastry_2001 (); mossa_dynamics_2002 (); sengupta_dependence_2011 () stems from using an alternative definition of , for which Eq. (1) holds over a broader range. We emphasize that our results conform to the general physical expectations that relaxation dynamics in supercooled liquids becomes thermally activated when temperature is low enough, typically below the mode-coupling temperature. The results exposed in Fig. 3 are therefore physically welcome, as there exists no fundamental reason for Eq. (1) to be relevant in the weakly supercooled temperature regime. Additionally, experiments find clear deviations from this relation in the temperature window covered by simulations richert_dynamics_1998 (). Therefore, our results suggest a plausible resolution to the existing discrepancy between experiments and previous simulations, although more work remains to be done, especially at lower temperatures, to fully settle this issue.
Despite the above successes, we emphasize that our determination of a configurational entropy remains an approximation to the theoretical concept of a complexity counting the number of metastable states. A first approximation stems from the fact that we perform measurements of with two freely evolving copies. An alternative procedure franz_phase_1997 () consists of first drawing copy 1 from the equilibrium distribution, before studying the thermodynamics of copy 2 in the presence of the quenched disorder imposed by copy 1. This amounts to distinguishing annealed from quenched complexities cardenas_constrained_1999 (). In the mean-field limit where rigorous calculations exist, the annealed is an approximation to the quenched one, but the latter is more fundamental because it exactly counts the number of metastable states. Although the quenched potential can also be measured berthier_overlap_2013 (), the procedure is more demanding. Before quantitatively comparing the two complexities, one should first establish more firmly the existence of a critical temperature in the quenched case. Preliminary work berthier_overlap_2013 () suggests a slight depression of the critical temperature , and very close values for both , but these issues need to be examined thoroughly.
A more fundamental issue concerns the interpretation of determined from as an entropy associated to the number of metastable states. This is true at the mean-field level, where both and the complexity can be rigorously defined and computed berthier_theoretical_2011 (); cavagna_supercooled_2009 (). The situation is ambiguous in finite dimensions, where infinitely long-lived metastable states do not exist, which forbids a strict definition of a complexity associated to their number mezard__???? (). Metastability can therefore only be approximately defined, for instance using finite timescales biroli_metastable_2001 () or lengthscales bouchaud_adam-gibbs-kirkpatrick-thirumalai-wolynes_2004 (), and metastable states cannot sharply emerge at the mode-coupling temperature, as they do in mean-field approximations franz_phase_1997 (); mezard__???? (). By contrast, we note that and defined in Eq. (2) do not suffer from these ambiguities. Therefore, our estimate of a configurational entropy is well-defined in finite dimensions, even though metastable states are not. This distinction also explains that a sharp emergence of at is found in our simulations, whereas only a weak vestige of the mode-coupling transition can be observed.
We can tentatively interpret , as measured here, as the entropy related to the number of ‘metastable’ states, now defined over finite lengthscales, suggesting a possible deep connection between the emergence of found here, and the growth of a static (point-to-set) correlation length kob_non-monotonic_2012 (); bouchaud_adam-gibbs-kirkpatrick-thirumalai-wolynes_2004 (); biroli_thermodynamic_2008 (); berthier_static_2012 (); hocky_growing_2012 (). We emphasize that since in Eq. (2) quantifies the free energy cost to localize the system in configuration space, its rapid decrease upon supercooling is likely related to the slowing down of the dynamics. This scenario naturally emerges in thermodynamic theories of the glass transition, such as Adam-Gibbs and random first-order transition theories. The numerical strategy proposed herein provides a sensible measure of the configurational entropy in the temperature range currently accessible to numerical simulations and will thus allow a stringent test of theoretical approaches in which the configurational entropy plays a central role.
Iv Materials and methods
which is truncated at distance . The particle types are labeled and and the interaction parameters are , and . We consider , 108 and 256 with . We express length scales in units of , temperatures in units of , and perform simulations at constant number density .
In the Kob-Andersen mixture, particles of types interact by a Lennard-Jones potential
which is truncated and shifted at . The particle types are labelled and and the interaction parameters are and . We consider particles with and . We express length scales in units of , temperatures in units of , and perform simulations at constant number density .
Simulation methods– Both models are studied numerically using Monte-Carlo simulations berthier_monte_2007 (). Measuring from Eq. (2) is numerically challenging as it requires the determination of over a broad range of , which necessitates a quantitative analysis of atypical overlap fluctuations. This difficulty is efficiently overcome by using umbrella sampling techniques frenkel_understanding_2001 (). Briefly, is obtained by gathering the results of a series of simulations biased in such a way that distinct simulations explore distinct ranges of overlap values berthier_overlap_2013 (). Histogram reweighting techniques are then used to reconstruct over the complete range berthier_overlap_2013 (). Another challenge is the difficulty to ensure proper thermalization of each simulation, which becomes serious when is large and is low. This is not prohibitive in studies where only the vicinity of the critical temperature is explored berthier_overlap_2013 (); parisi_liquid-glass_2013 (). To access much lower temperatures, we have introduced replica-exchange Monte-Carlo moves between the biased simulations, borrowing techniques used in phase transition studies frenkel_understanding_2001 (); hukushima_exchange_1996 (); yan_hyper-parallel_1999 (). Because each simulation now performs a random walk in parameter space, thermalization is greatly enhanced, and lower temperatures can be sampled. By approaching the mode-coupling temperature, we are able to measure over a physically significant -range. The procedure can presumably be further optimized to access even lower temperatures. It can also easily be applied to different models, including hard spheres to which the inherent structure formalism does not apply angelani_configurational_2007 (). Full details about these methods as well as about finite size effects are given in the SI.
Acknowledgements.We thank G. Biroli, G. Parisi and G. Tarjus for useful exchanges. The research leading to these results has received funding from the European Research Council under the European Union’s Seventh Framework Programme (FP7/2007-2013) / ERC Grant agreement No 306845.
- (1) Ediger, M. D, Angell, C. A, and Nagel, S. R. (1996) Supercooled liquids and glasses. J. Phys. Chem. 100, 13200.
- (2) Berthier, L and Biroli, G. (2011) Theoretical perspective on the glass transition and amorphous materials. Reviews of Modern Physics 83, 587.
- (3) (year?) The Nature of the Glassy State and the Behavior of Liquids at Low Temperatures. Chem. Rev. 43.
- (4) Richert, R and Angell, C. A. (1998) Dynamics of glass-forming liquids. V. On the link between molecular dynamics and configurational entropy. J. Chem. Phys. 108, 9016.
- (5) Martinez, L.-M and Angell, C. A. (2001) A thermodynamic connection to the fragility of glass-forming liquids. Nature 410, 663.
- (6) Angell, C. A and Borick, S. (2002) Specific heats Cp, Cv, Cconf and energy landscapes of glassforming liquids. J. Non-Cryst. Solids 307, 393.
- (7) Stillinger, F. H and Weber, T. A. (1982) Hidden structure in liquids. Phys. Rev. A 25, 978.
- (8) Sciortino, F, Kob, W, and Tartaglia, P. (1999) Inherent Structure Entropy of Supercooled Liquids. Phys. Rev. Lett. 83, 3214.
- (9) Scala, A, Starr, F. W, La Nave, E, Sciortino, F, and Stanley, H. E. (2000) Configurational entropy and diffusivity of supercooled water. Nature 406, 166.
- (10) Sastry, S. (2001) The relationship between fragility, configurational entropy and the potential energy landscape of glass-forming liquids. Nature 409, 164.
- (11) Mossa, S, La Nave, E, Stanley, H. E, Donati, C, Sciortino, F, and Tartaglia, P. (2002) Dynamics and configurational entropy in the Lewis-Wahnström model for supercooled orthoterphenyl. Phys. Rev. E 65, 041205.
- (12) Angelani, L and Foffi, G. (2007) Configurational entropy of hard spheres. J. Phys.: Condens. Matter 19, 256207.
- (13) Adam, G and Gibbs, J. H. (1965) On the Temperature Dependence of Cooperative Relaxation Properties in Glass-Forming Liquids. J. Chem. Phys. 43, 139.
- (14) Biroli, G and Monasson, R. (2000) From inherent structures to pure states: Some simple remarks and examples. EPL (Europhysics Letters) 50, 155.
- (15) Wolynes, P. G. (1997) Entropy Crises in Glasses and Random Heteropolymers. J. Res. Natl. Inst. Stand. Technol. 102, 187.
- (16) Cavagna, A. (2009) Supercooled liquids for pedestrians. Phys. Rep. 476, 51.
- (17) Sengupta, S, Vasconcelos, F, Affouard, F, and Sastry, S. (2011) Dependence of the fragility of a glass former on the softness of interparticle interactions. J. Chem. Phys. 135, 194503.
- (18) Monasson, R. (1995) Structural Glass Transition and the Entropy of the Metastable States. Phys. Rev. Lett. 75, 2847.
- (19) Mézard, M and Parisi, G. (1999) Thermodynamics of Glasses: A First Principles Computation. Phys. Rev. Lett. 82, 747.
- (20) Mézard, M and Parisi, G. (2012) in Structural glasses and supercooled liquids, eds. Wolynes, P. G and Lubchenko, V. (Wiley & Sons).
- (21) Franz, S and Parisi, G. (1997) Phase Diagram of Coupled Glassy Systems: A Mean-Field Study. Phys. Rev. Lett. 79, 2486.
- (22) Cardenas, M, Franz, S, and Parisi, G. (1999) Constrained Boltzmann-Gibbs measures and effective potential for glasses in hypernetted chain approximation and numerical simulations. J. Chem. Phys. 110, 1726.
- (23) Cammarota, C, Cavagna, A, Giardina, I, Gradenigo, G, Grigera, T. S, Parisi, G, and Verrocchio, P. (2010) Phase-Separation Perspective on Dynamic Heterogeneities in Glass-Forming Liquids. Phys. Rev. Lett. 105, 055703.
- (24) Berthier, L. (2013) Overlap fluctuations in glass-forming liquids. Phys. Rev. E 88, 022313.
- (25) Parisi, G and Seoane, B. (2013) Liquid-glass transition in equilibrium, arxiv:1311.1465.
- (26) Berthier, L and Kob, W. (2007) The Monte Carlo dynamics of a binary Lennard-Jones glass-forming mixture. J. Phys.: Condens. Matter 19, 205130.
- (27) Berthier, L and Witten, T. A. (2009) Compressing nearly hard sphere fluids increases glass fragility. EPL 86, 10001.
- (28) Berthier, L and Witten, T. A. (2009) Glass transition of dense fluids of hard and compressible spheres. Phys. Rev. E 80, 021502.
- (29) Kob, W, Roldán-Vargas, S, and Berthier, L. (2012) Non-monotonic temperature evolution of dynamic correlations in glass-forming liquids. Nature Phys. 8, 164.
- (30) Flenner, E and Szamel, G. (2013) Dynamic heterogeneities above and below the mode-coupling temperature: Evidence of a dynamic crossover. J. Chem. Phys. 138, 12A523–12A523.
- (31) Kob, W and Andersen, H. C. (1995) Testing mode-coupling theory for a supercooled binary Lennard-Jones mixture I: The van Hove correlation function. Phys. Rev. E 51, 4626.
- (32) Mézard M. and Parisi G., (2012) Glasses and replicas, in Structural glasses and supercooled liquids, Eds: P. G. Wolynes and V. Lubchenko (Wiley & Sons).
- (33) Sastry, S, Debenedetti, P. G, and Stillinger, F. H. (1998) Signatures of distinct dynamical regimes in the energy landscape of a glass-forming liquid. Nature 393, 554.
- (34) Biroli, G and Kurchan, J. (2001) Metastable states in glassy systems. Phys. Rev. E 64, 016101.
- (35) Bouchaud, J.-P and Biroli, G. (2004) On the Adam-Gibbs-Kirkpatrick-Thirumalai-Wolynes scenario for the viscosity increase in glasses. J. Chem. Phys. 121, 7347.
- (36) Biroli, G, Bouchaud, J.-P, Cavagna, A, Grigera, T. S, and Verrocchio, P. (2008) Thermodynamic signature of growing amorphous order in glass-forming liquids. Nat Phys 4, 771.
- (37) Berthier, L and Kob, W. (2012) Static point-to-set correlations in glass-forming liquids. Phys. Rev. E 85, 011102.
- (38) Hocky, G. M, Markland, T. E, and Reichman, D. R. (2012) Growing Point-to-Set Length Scale Correlates with Growing Relaxation Times in Model Supercooled Liquids. Phys. Rev. Lett. 108, 225506.
- (39) Frenkel, D and Smit, B. (2001) Understanding Molecular Simulation. (Academic Press), 2 edition.
- (40) Hukushima, K and Nemoto, K. (1996) Exchange Monte Carlo method and application to spin glass simulations. J. Phys. Soc. Japan 65, 1604.
- (41) Yan, Q and de Pablo, J. J. (1999) Hyper-parallel tempering Monte Carlo: Application to the Lennard-Jones fluid and the restricted primitive model. J. Chem. Phys. 111, 9509.
Appendix A Supplementary information
In this appendix, we describe the umbrella sampling, parallel tempering and histogram reweighting techniques we used to measure numerically the configurational entropy defined by Eq. (2), and discuss finite size effects.
a.1 Umbrella sampling
To access the probability distribution function of the overlap at a given temperature, we conduct distinct simulations. In each simulation, , two copies of the system evolve according to the following Hamiltonian:
where is the Hamiltonian of the original liquid, and respectively represent the positions of the particles in copies 1 and 2, is the overlap between configurations 1 and 2, is the thermodynamic field conjugated to the overlap, and the biasing potential is taken of the form:
with parameters chosen to constrain the overlap to explore values away from its average equilibrium value.
We perform simple Monte-Carlo moves, where we attempt single particle displacements of small amplitude SIMC () (typically of size 0.05 ) which we accept using a Metropolis acceptance rate given by Hamiltonian (3). We define time steps such that one Monte-Carlo step represents attempts to displace randomly chosen particles among the two copies of the system.
Provided it is properly thermalized (see below), the main outcome of a given simulation is the probability distribution function of the overlap,
where represents a thermal average with the Hamiltonian in Eq. (3) for a given state point defined by .
The idea behind the biasing potentials in Eqs. (3, 4) is that the fluctuations of the overlap in each simulation can be tailored to explore a narrow region centered around . Therefore, each simulation explores only a small range of overlap values, and it becomes unnecessary to wait for very rare overlap fluctuations to occur. Therefore, umbrella sampling is the key method to efficiently measure atypical overlap fluctuations SIumbrella ().
In Fig. 4 we show the distributions measured in simulations of harmonic spheres at for for a given set of biasing potentials chosen to adequately cover the range of overlap between 0 and 1. We used values in the interval , and a strength of the Gaussian trap in the interval . We observe that each simulation returns a probability distribution function which is relatively narrow and approximately Gaussian, showing that the sampling of each overlap sector is no more slowed down by the emergence of non-trivial distributions of overlap fluctuations SIumbrella (). In other words, simulations are faster because the Hamiltonian has no phase transition in the plane.
While the umbrella sampling technique described above considerably accelerates the measurement of the overlap fluctuations, we have observed that when is low, is large, and/or is large, the particle dynamics slows down considerably, and it becomes difficult to perform an accurate sampling of the overlap fluctuations imposed the Hamiltonian (3), because the overlap fluctuations become slow. This sampling problem was already mentioned in Ref. SIpre (), and this prevented the exploration of temperatures much lower than .
To accurately explore the temperature regime needed for the present work, we have implemented replica-exchanges Monte-Carlo moves SIPT (). We now conduct the simulations needed for the umbrella sampling at temperature in parallel, and propose Monte-Carlo exchange moves between neighboring simulations characterized by nearby sets of parameters, say and . An exchange between simulations and is proposed with a low frequency (typically every Monte-Carlo steps) and they are accepted with a Metropolis acceptance rate given by the Hamiltonians and , ensuring that simulations satisfy detailed balance.
Because each simulation now performs a random walk in the parameter space defined by , the sampling of the overlap fluctuations is greatly enhanced, even for the ‘hard’ cases. For the method to be efficient, we need to adjust the biasing potentials such that the overlap distributions in each simulation largely overlap, as can be seen in the example shown in Fig. 4. We have used up to parallel simulations to gather our data. We have carefully checked thermalization by running simulations for very long times (up to Monte-Carlo steps per simulation), making sure that each state point was visited several times by all simulations due to the replica-exchange. This represents a significant numerical effort.
a.3 Histogram reweighting
Having obtained thermalized results from biased simulations running in parallel, we process the simulation outcome using multi-histogram reweighting methods to reconstruct the unbiased probability from the independently measured SIumbrella (); SIpre (),
where the are defined self-consistently as
Notice that the value of used in the simulations plays no conceptual role because the reweighting method allows us to directly obtain from for distinct field values and :
We have applied the combined umbrella sampling / replica-exchange technique SIyan () both with for the Lennard-Jones potential, and with adjusted to be roughly at phase coexistence for the harmonic sphere system. We found that the latter method allows for an easier convergence of the simulation parameters, without affecting qualitatively the efficiency of statistical sampling.
Two values of the field are particularly relevant for this study. First, we obtain the potential as:
Notice that is only defined up to an additive constant, which we adjust such that , where is defined as the location of the global minimum in . This additive constant is irrelevant, as we only need to determine free energy differences to determine . In Fig. 4 we show the distribution for the chosen example. The potential in Fig. 1 is simply obtained by taking (minus) the logarithm of this function which allows for a clearer view of the tail of the distribution.
A second useful quantity is the overlap distribution obtained in the presence of a field ensuring phase coexistence below . In practice, we use the strength of the reweighting method to finely explore a range of values, and define as the field for which the distribution yields a maximum amount of fluctuations (as quantified by their variance). This corresponds to the situation where the two peaks of the distribution have nearly equal height. From this bimodal distribution we define the position of the high-overlap value , as illustrated in Fig. 1.
It is useful to notice that has a simple graphical interpretation, since it represents by definition the amplitude of the field needed to ‘tilt’ the potential towards coexistence, see Fig. 1. Therefore, it represents a simple proxy to the free energy difference defined in Fig. 1, because the relation
holds to a good approximation. We find that differences between various possible estimates of SIMPbook () can only be distinguished in a very narrow temperature regime near which is not resolved by the present set of data, and therefore does not affect any of our conclusions. This estimate is also useful because it stems from a single number (the critical field ), and is less prone to statistical errors than the entire function .
We also remark that, contrary to mean-field calculations, cannot be defined from the existence of a secondary minimum in , because the potential has to be a convex function of in finite dimensions. Therefore, a secondary minimum cannot exist in the large system size limit in our simulations. This is why we have instead defined from the position of the large- peak in at coexistence.
In the same vein, the appearance in mean-field calculations of a secondary minimum in is the direct signature of the mode-coupling transition at . Because no local minimum appears in at any temperature in finite dimensions, we do not expect to detect any qualitative change in in the region of the mode-coupling temperature in our simulations. This immediately suggests that can not be expected to play any significant role in the present study.
a.4 Finite size effects
For the harmonic sphere model we have used 3 different system sizes. This allows us to discuss how finite size effects affect our determination of the configurational entropy in this system.
Because our definition of a configurational entropy results from the existence of an underlying phase transition for the constrained Hamiltonian (3), finite size effects should naturally be carefully discussed. However, while a strong system size dependence might be expected to affect the overlap fluctuations near the critical temperature , much smaller finite size effects can be expected in the low-temperature regime of interest here. Therefore a finite size system could potentially affect the value of and the free energy difference between the two phases, but these can be expected to be small, just as finite size effects would also be small in the magnetic phase of the Ising model.
To test this idea, we present in Fig. 5 the value of the critical field for different system sizes and temperatures in the harmonic sphere model, which represents a faithful estimate to the configurational entropy, see Eq. (10). The data shown in Fig. 5 clearly indicate that finite size effects are small in the temperature regime explored by the present simulations. We find similarly that the value of is not strongly affected by finite size effects. We expect that a stronger influence of the system size could be observed at much lower temperatures, when the system size competes more strongly with the growing point-to-set correlation length of the system, which is found to be relatively modest in the temperature regime above the mode-coupling transition SIpts3 ().
- (S1) L. Berthier and W. Kob, J. Phys.: Condens. Matter 19, 205130 (2007).
- (S2) D. Frenkel and B. Smit, Understanding Molecular Simulation: From Algorithms to Applications (Academic Press, San Diego, 2001).
- (S4) L. Berthier, Phys. Rev. E 88, 022313 (2013).
- (S3) K. Hukushima and K. Nemoto, J. Phys. Soc. Japan 64, 1604 (1996).
- (S5) Q. Yan and J. J. de Pablo, J. Chem. Phys. 111, 9509 (1999).
- (S6) M. Mézard and G. Parisi, in Structural glasses and supercooled liquids, Eds: P. G. Wolynes and V. Lubchenko (Wiley & Sons, 2012).
- (S7) L. Berthier and W. Kob, Phys. Rev. E 85, 011102 (2012).