Rahul N Chacko

François P Landes

Giulio Biroli

Olivier Dauchot

Andrea J Liu

David R Reichman

Dynamical Facilitation Governs the Equilibration Dynamics of Glasses

Convincing evidence of domain growth in the heating of ultrastable glasses suggests that the equilibration dynamics of super-cooled liquids could be driven by a nucleation and growth mechanism. We investigate this possibility by simulating the equilibration dynamics of a model glass during both heating and cooling between poorly and well-annealed states. Though we do observe the growth of domains during heating, we find that domains are absent during cooling. This absence is inconsistent with classical nucleation theory. By comparing the equilibration dynamics of our glass with that of two models with kinetic constraints, we demonstrate that dynamical facilitation generically leads to heating driven by domain growth and cooling without domains. Our results provide strong evidence that dynamical facilitation, not nucleation and interfacial-tension-driven domain growth, is the driving mechanism for the equilibration dynamics of glass-formers.

I. INTRODUCTION

The phenomenology of glasses and supercooled liquids (i.e., glasses in metastable equilibrium such that detailed balance is satisfied on all relevant time scales) bears a striking resemblance to that of crystalline solids [START_REF] Kirkpatrick | Scaling concepts for the dynamics of viscous liquids near an ideal glassy state[END_REF][START_REF] Douglass | Can a stable glass be superheated? Modelling the kinetic stability of coated glassy films[END_REF][START_REF] Berthier | Evidence for a Disordered Critical Point in a Glass-Forming Liquid[END_REF][START_REF] Jack | The melting of stable glasses is governed by nucleation-and-growth dynamics[END_REF][START_REF] Cubeta | Glass softening kinetics in the limit of high heating rates[END_REF][START_REF] Vila-Costa | Nucleation and Growth of the Supercooled Liquid Phase Control Glass Transition in Bulk Ultrastable Glasses[END_REF][START_REF] Herrero | Two-step devitrification of ultrastable glasses[END_REF][START_REF] Herrero | Front propagation in ultrastable glasses is dynamically heterogeneous[END_REF]. As with transformations of Eshelby disclinations in crystals [START_REF] Eshelby | The determination of the elastic field of an ellipsoidal inclusion, and related problems[END_REF], local plastic events in these systems perturb their neighborhoods with long-ranged quadrupolar displacement fields [START_REF] Picard | Elastic consequences of a single plastic event: A step towards the microscopic modeling of the flow of yield stress fluids[END_REF][START_REF] Chacko | Elastoplasticity Mediates Dynamical Heterogeneity Below the Mode Coupling Temperature[END_REF]. A large body of experimental evidence also finds that ultrastable glasses evolve upon heating in a manner consistent with the Avrami framework of classical nucleation theory [START_REF] Fanfoni | The Johnson-Mehl-Avrami-Kohnogorov Model: A Brief Review[END_REF]: fronts of high-temperature supercooled liquid advance at constant speed into the bulk ultrastable glass [START_REF] Cubeta | Glass softening kinetics in the limit of high heating rates[END_REF][START_REF] Swallen | Stable Glass Transformation to Supercooled Liquid via Surface-Initiated Growth Front[END_REF][START_REF] Kearns | Observation of low heat capacities for vapor-deposited glasses of indomethacin as determined by AC nanocalorimetry[END_REF][START_REF] Kearns | One Micrometer Length Scale Controls Kinetic Stability of Low-Energy Glasses[END_REF][START_REF] Dawson | Anisotropic Structure and Transformation Kinetics of Vapor-Deposited Indomethacin Glasses[END_REF][START_REF] Whitaker | Vapor-deposited α,α,β-tris-naphthylbenzene glasses with low heat capacity and high kinetic stability[END_REF][START_REF] Sepúlveda | Anomalous Transformation of Vapor-Deposited Highly Stable Glasses of Toluene into Mixed Glassy States by Annealing Above Tg[END_REF][START_REF] Chen | Dynamics of glass-forming liquids. XVI. Observation of ultrastable glass transformation via dielectric spectroscopy[END_REF][START_REF] Sepúlveda | Manipulating the properties of stable organic glasses using kinetic facilitation[END_REF][START_REF] Sepúlveda | Role of Fragility in the Formation of Highly Stable Organic Glasses[END_REF][START_REF] Rodríguez-Tinoco | Evaluation of Growth Front Velocity in Ultrastable Glasses of Indomethacin over a Wide Temperature Interval[END_REF][START_REF] Bhattacharya | Enthalpy and high temperature relaxation kinetics of stable vapor-deposited glasses of toluene[END_REF][START_REF] Rodríguez-Tinoco | Transformation kinetics of vapor-deposited thin film organic glasses: The role of stability and molecular packing anisotropy[END_REF][START_REF] Dalal | Influence of Substrate Temperature on the Transformation Front Velocities That Determine Thermal Stability of Vapor-Deposited Glasses[END_REF][START_REF] Tylinski | Vapor-deposited glasses of methylm -toluate: How uniform is stable glass transformation?[END_REF][START_REF] Walters | Thermal stability of vapor-deposited stable glasses of an organic semiconductor[END_REF][START_REF] Ràfols-Ribé | The role of thermodynamic stability in the characteristics of the devitrification front of vapour-deposited glasses of toluene[END_REF][START_REF] Cubeta | Communication: Surface-facilitated softening of ordinary and vapor-deposited glasses[END_REF], with the fraction of material retaining ultrastable glass structure decaying with time according to a compressed exponential law [START_REF] Kearns | One Micrometer Length Scale Controls Kinetic Stability of Low-Energy Glasses[END_REF][START_REF] Dawson | Anisotropic Structure and Transformation Kinetics of Vapor-Deposited Indomethacin Glasses[END_REF][START_REF] Sepúlveda | Anomalous Transformation of Vapor-Deposited Highly Stable Glasses of Toluene into Mixed Glassy States by Annealing Above Tg[END_REF][START_REF] Fullerton | Density controls the kinetic stability of ultrastable glasses[END_REF]. By contrast, poorlyannealed ("ordinary") glasses are found to transform homogeneously into the supercooled liquid upon heating above the glass transition temperature [START_REF] Sepúlveda | Manipulating the properties of stable organic glasses using kinetic facilitation[END_REF][START_REF] Walters | Thermal stability of vapor-deposited stable glasses of an organic semiconductor[END_REF]. This difference of behavior, also observed in pressure-controlled molecular dynamics (MD) simulations [START_REF] Herrero | Two-step devitrification of ultrastable glasses[END_REF][START_REF] Herrero | Front propagation in ultrastable glasses is dynamically heterogeneous[END_REF][START_REF] Fullerton | Density controls the kinetic stability of ultrastable glasses[END_REF][START_REF] Flenner | Front-Mediated Melting of Isotropic Ultrastable Glasses[END_REF], raises the possibility that a thermodynamic first-order transition mediated by nucleation and domain growth separates poorly-annealed glass from well-annealed glass [6-8, 32, 33], as in the case of crystal freezing and melting.

Indeed, as with crystals, a notion of structural order also applies to supercooled liquids. In these amorphous systems, local equilibrium structure constrains [START_REF] Das | Soft pinning: Experimental validation of static correlations in supercooled molecular glass-forming liquids[END_REF] and is constrained by [START_REF] Berthier | Zero-temperature glass transition in two dimensions[END_REF] neighboring structure up to a finite distance ξ PTS (the "point-to-set" length) away [START_REF] Biroli | The Random First-Order Transition Theory of Glasses: A Critical Assessment[END_REF]. Random First Order Transition Theory (RFOT), known to be exact in infinite dimension [START_REF] Charbonneau | Exact theory of dense amorphous hard spheres in high dimension. III. The full replica symmetry breaking solution[END_REF], predicts systems below a critical temperature T MCT to be separated into a mosaic of ξ PTS -sized domains of mutually-incompatible amorphous order [START_REF] Kirkpatrick | Scaling concepts for the dynamics of viscous liquids near an ideal glassy state[END_REF][START_REF] Biroli | The Random First-Order Transition Theory of Glasses: A Critical Assessment[END_REF][START_REF] Bouchaud | On the Adam-Gibbs-Kirkpatrick-Thirumalai-Wolynes scenario for the viscosity increase in glasses[END_REF]. As demonstrated by randomly pinned glasses [START_REF] Hocky | Equilibrium ultrastable glasses produced by random pinning[END_REF] and supercooled liquids coupled to a reference configuration [START_REF] Berthier | Evidence for a Disordered Critical Point in a Glass-Forming Liquid[END_REF], amorphous order can in principle lead to first order transitions in amorphous systems, supporting the possibility of a first-order transition in supercooled liquids [START_REF] Cammarota | A general approach to systems with randomly pinned particles: Unfolding and clarifying the Random Pinning Glass Transition[END_REF].

Nucleation and domain growth is not the only available explanation for the equilibration dynamics of glassformers, however. Particles in supercooled liquids are trapped into "cages" by their nearest neighbors, rattling around for long periods between cage-escape events in which they hop out of one cage and into another [START_REF] Weeks | Three-Dimensional Direct Imaging of Structural Relaxation Near the Colloidal Glass Transition[END_REF][START_REF] Candelier | Spatiotemporal Hierarchy of Relaxation Events, Dynamical Heterogeneities, and Structural Reorganization in a Supercooled Liquid[END_REF]. As a local plastic event within an otherwise elastic material, a cage escape weakens its neighborhood, creating the conditions for more cage escapes [START_REF] Chacko | Elastoplasticity Mediates Dynamical Heterogeneity Below the Mode Coupling Temperature[END_REF]. Such a process, in which rearrangements propagate mobility to their local neighborhoods, is known as dynamical facilitation [START_REF] Keys | Excitations Are Localized and Relaxation Is Hierarchical in Glass-Forming Liquids[END_REF][START_REF] Tarjus | An overview of the theories of the glass transition, in Dynamical Heterogeneities in Glasses, Colloids, and Granular Media[END_REF] and plays an important role in glassy liquids [START_REF] Dauchot | Dynamical heterogeneities in grains and foams[END_REF]. This self-propagation of mobility results in avalanches of plastic events correlated in space and time [START_REF] Candelier | Spatiotemporal Hierarchy of Relaxation Events, Dynamical Heterogeneities, and Structural Reorganization in a Supercooled Liquid[END_REF]. On long time scales, these avalanches appear as large clusters of spatially-correlated mobility that grow as a function of time scale [START_REF] Keys | Excitations Are Localized and Relaxation Is Hierarchical in Glass-Forming Liquids[END_REF][START_REF] Scalliet | Thirty Milliseconds in the Life of a Supercooled Liquid[END_REF], resembling domain growth. Dynamical facilitation plays an important role in the equilibrium dynamics of supercooled liquids, and it could also be the key mechanism at play in equilibration dynamics, i.e. relaxation towards equilibrium. Indeed, many of the phenomena observed in glasses upon heating or cooling are also seen in Kinetically Constrained Models (KCMs), spin lattice models in which spins interact only via facilitation [START_REF] Garrahan | Kinetically constrained models[END_REF], upon heating [START_REF] Douglass | Can a stable glass be superheated? Modelling the kinetic stability of coated glassy films[END_REF][START_REF] Butler | Glassy relaxation at surfaces: The correlation length of cooperative motion in the facilitated kinetic Ising model[END_REF][START_REF] Gutiérrez | Front propagation versus bulk relaxation in the annealing dynamics of a kinetically constrained model of ultrastable glasses[END_REF] or cooling [START_REF] Garrahan | Glassiness and constrained dynamics of a short-range nondisordered spin model[END_REF]. Dynamical facilitation therefore offers a possible alternative explanation for the behavior of glasses during heating and cooling. In this paper, we show that it is indeed dynam-ical facilitation that governs the equilibration dynamics of glasses, not nucleation followed by domain growth.

We first address the question of nucleation, interpreted in the classical sense [START_REF] Kashchiev | Nucleation: Basic Theory with Applications[END_REF] of a competition between bulk and interfacial free energies which sets a critical size above which randomly-nucleated domains will grow and below which they will shrink. We distinguish nucleation, which involves a thermodynamic drive towards domain growth as small domains (of linear size smaller than the nucleation length) are suppressed, from the evolution of structure via random, thermally-activated rearrangements, with no suppression of small domains. In the latter case, an additional mechanism is required to explain domain growth. In the second part of this paper, we show that dynamical facilitation fills this role.

We conduct MD simulations of a polydisperse glass at fixed temperature T ∞ and pressure P for two types of initial condition. In experiments using the first type of initial condition, which we refer to as homogeneous experiments, the initial state corresponds to equilibrium at an initial temperature T 0 different from the thermostat temperature T ∞ . In experiments using the second type of initial condition, which we refer to as slab experiments, we introduce a system-percolating slab of T ∞ equilibrium structure lying in the x ≡ 0 plane that interrupts the T 0 structure of the surrounding system (see Figs. [START_REF] Jack | The melting of stable glasses is governed by nucleation-and-growth dynamics[END_REF](a) and 4(b)). In this geometry, the slab can expand while keeping the size and shape of its boundary fixed, avoiding any putative interfacial free energy cost as it expands. In a classical nucleation theory scenario, the slab represents a nucleus above the critical nucleus size. We fix T low = min {T 0 , T ∞ } and T high = max {T 0 , T ∞ } so that any critical point separating T 0 and T ∞ during heating will also do so during cooling. These simulations provide a number of tests for nucleation and domain growth:

1. We can characterize the local structure and dynamics, and simply look (by eye) for growing domains of T ∞ structure in individual realizations of the homogeneous experiment.

2. We can also look for evidence of a nucleation time in these experiments.

3. We can quantitatively track domain growth via the displacement of the slab boundary in slab experiments.

4. We can assess the roughness of domains, both qualitatively in individual realizations of the homogeneous experiment and quantitatively (via the width of the slab boundary) in the slab experiment.

As we will show, all four tests argue against nucleation, with the most striking evidence provided by the absence of domain growth during cooling, as we shall discuss. Before applying this methodology to our polydisperse model glass, we validate these tests by applying them to a crystallizing monodisperse system with temperatures T low and T high straddling its melting point, verifying our expectations for the case of a genuine first order transition.

After testing the nucleation and growth mechanism and obtaining results inconsistent with it, we support our claim that dynamical facilitation is responsible for the equilibration dynamics of glasses by applying the same four tests to two plaquette models, the Triangular Plaquette Model (TPM) [START_REF] Garrahan | Glassiness and constrained dynamics of a short-range nondisordered spin model[END_REF][START_REF] Garrahan | Glassiness through the emergence of effective dynamical constraints in interacting systems[END_REF][START_REF] Garrahan | Transition in coupled replicas may not imply a finite-temperature ideal glass transition in glassforming systems[END_REF] and the Square Pyramidal Plaquette Model (SPPM) [START_REF] Turner | Overlap and activity glass transitions in plaquette spin models with hierarchical dynamics[END_REF]. Plaquette models are spin lattice models with a dual representation in terms of plaquettes with trivial thermodynamics and kinetically-constrained dynamics [START_REF] Ritort | Glassy dynamics of kinetically constrained models[END_REF], corresponding to KCMs in the plaquette representation. The TPM and SPPM are similar by design, and the mechanisms underlying dynamical facilitation in these two systems differ greatly from the elastoplastic mechanism for facilitation in supercooled liquids and glasses [START_REF] Chacko | Elastoplasticity Mediates Dynamical Heterogeneity Below the Mode Coupling Temperature[END_REF]. We find that just like our polydisperse model glass, these plaquette models exhibit domain growth during heating and none during cooling. For all three systems, domain growth during heating proceeds at a constant speed that persists even as the system relaxes. The only notable difference between the systems is the evolution of domain roughness, which is different for all three systems and therefore not a test for facilitation-driven dynamics.

II. SYSTEM

Our three-dimensional monodisperse crystal and twodimensional polydisperse glass systems comprise particles of uniform mass m interacting via the pair potential

V (r) = V 0 r-12 + c 0 + c 2 r2 + c 4 r4 r ≤ 1.25, 0 r > 1.25, (1)
where

r = 2r/ [(σ i + σ j ) (1 -0.2 |σ i -σ j |)] .
Here, r is the separation between particles i and j with sizes σ i and σ j , and the constants c 0 , c 2 , and c 4 are chosen such that

V (1.25) = V ′ (1.25) = V ′′ (1.25) = 0 ( §A 1).
Following [START_REF] Ninarello | Models and Algorithms for the Next Generation of Glass Transition Studies[END_REF], for the polydisperse system, σ i for each particle is drawn independently from a distribution with probability density ∝ σ -3 in the interval [σ min , σ max ], with σ min and σ max chosen such that the coefficient of variation of σ is c σ = 0.23. This system, introduced in [START_REF] Ninarello | Models and Algorithms for the Next Generation of Glass Transition Studies[END_REF], is designed to allow equilibration to temperatures well below T MCT using swap Monte Carlo simulations. We take V 0 , σ, m, and V 0 /k B to be our energy, length, mass, and temperature units, where k B is the Boltzmann constant and where σ is the expected value of the particle diameter distribution. The system coordinates are centered at the origin and periodic boundary conditions are imposed along each axis of the system. We prepare our initial states using constant-volume Monte Carlo simulations implemented as described in §A 2. For the monodisperse case, the homogeneous experiments have N = 27648 particles and initial lengths 48, 12 √ 3, and 8 √ 6 along the x, y and z axes respectively, corresponding to number density ρ = √ 2. These system dimensions are compatible with a perfect FCC crystal, the ground state for this system ( §A 7). For melting, we heat the system from an initial configuration equilibrated to T low = 1.90 using a thermostat at T high = 2.35, while for freezing, we cool the system from an initial configuration equilibrated to T high = 2.10 using a thermostat at T low = 1.48.

For the slab experiments in the monodisperse system, we fix T low = 2.0 and T high = 2.1, which straddle the melting point. We preserve the compatibility of the system dimensions with a perfect FCC crystal by not attempting to eliminate the initial pressure gradient at the slab boundary in slab experiments, other than by our choice of narrowly-separated T low and T high . Instead, we double the number of particles and the system size along the x-axis, as compared to the systems used for homogeneous experiments. This allows for time between the early-time pressure wave dissipating and the left and right slab fronts colliding at the periodic boundary of the system. For freezing, we define our frozen slab to be the region with x-coordinate satisfying |x| < 6. For melting, our preparation protocol ( §A 2) requires us to specify the size of the frozen non-slab region instead, leaving the size of the melted slab region to be an emergent property of our preparation protocol dependent on the size of the interface. We define the non-slab region to be that with x-coordinate satisfying |x| ≥ 8.

For the polydisperse case, we simulate N = 20000 particles and fix T low = 0.06 and T high = 0.14 for both homogeneous and slab experiments. The homogeneous heating experiments have initial lengths 200.00 and 100.00 along the x and y axes respectively, corresponding to number density ρ = 1. At this density, T MCT ≈ 0.12 [START_REF] Scalliet | Thirty Milliseconds in the Life of a Supercooled Liquid[END_REF], such that (accounting for density changes [START_REF] Sengupta | Densitytemperature scaling of the fragility in a model glassformer[END_REF]) T low and T high straddle T MCT . The homogeneous heating experiments follow these systems to steady state at T high , corresponding to (average) system lengths 204.70 and 102.35 (ρ ≈ 0.955) along the x and y axes. These T high equilibrium configurations are then used as the initial T 0 = T high states for the homogeneous cooling experiments.

Our choice of T low = 0.06 and T high = 0.14 comes from practical considerations. We require T low < T MCT to ensure that T low is low enough for facilitation to play an important role in the equilibrium dynamics [START_REF] Keys | Excitations Are Localized and Relaxation Is Hierarchical in Glass-Forming Liquids[END_REF] and for amorphous order to set in [START_REF] Kirkpatrick | Scaling concepts for the dynamics of viscous liquids near an ideal glassy state[END_REF][START_REF] Bouchaud | On the Adam-Gibbs-Kirkpatrick-Thirumalai-Wolynes scenario for the viscosity increase in glasses[END_REF]. The lower we set T low , the stronger we expect dynamical facilitation to be and the slower we expect the dynamics to be. We choose T low = 0.06 as a compromise between strong facilitation and dynamics fast enough that we are able to see appreciable structural evolution on the time scale of a week of run time. We have checked that increasing T high , even above the onset temperature, does not affect our results. However, choosing T high = 0.14, below the onset temperature, ensures that we are probing a potential ultrastable-to-poorly-annealed glass transition, rather than a transition between an ultrastable glass and a high temperature liquid.

For the slab experiments in the polydisperse system, we define the slab to be the region satisfying |x| < b, where b = 25.00 for heating and 25.56 for cooling. In slab experiments, we greatly reduce the initial pressure gradient at the slab boundary by first equilibrating the system in its entirety to T ∞ , with system dimensions 200.00 and 100.00 if T ∞ = T low or 204.70 and 102.35 if T ∞ = T high , then shrinking or stretching the non-slab region along the x-axis such that ρ ≈ 0.955 in this region if T 0 = T high or ρ = 1 if T 0 = T low , before finally evolving this non-slab region with Monte Carlo moves at temperature T 0 until it reaches steady state.

Our MD simulations are performed in LAMMPS [START_REF] Plimpton | Fast Parallel Algorithms for Short-Range Molecular Dynamics[END_REF] (see §A 3), adopting a Nosé-Hoover thermostat and a Berendsen barostat that keeps the relative sizes of the system along each axis constant. In order to minimize the perturbation to the system due to the start-up of the barostat, we fix the pressure such that the virial contribution to the pressure at the beginning of the MD simulation matches that of the target T ∞ equilibrium state. For any given MD experiment, we average our results over 100 realizations, except when explicitly considering individual realizations or when performing the homogeneous experiments on the polydisperse glass, in which case we average over only 10 realizations.

We define time t = 0 as the time at which the MD simulations begin. For our analysis, it will be necessary to measure the distributions of structural and dynamical quantities in the initial, untransformed state during homogeneous experiments, as well as at equilibrium at the target temperature T ∞ . Our methodology for obtaining these distributions can be found in §A 4.

III. CRYSTAL CASE STUDY

We first validate our tests by applying them to the first-order melting and freezing transitions of a monodisperse system of particles, in which the growth of nuclei is driven by bulk free energy differences [START_REF] Peng | In situ observation of coalescence of nuclei in colloidal crystalcrystal transitions[END_REF], with structural frustration at the boundary of the nucleus suppressing roughness and promoting domain growth [START_REF] Penn | Imperfect Oriented Attachment: Dislocation Generation in Defect-Free Nanocrystals[END_REF][START_REF] Penn | Oriented attachment and growth, twinning, polytypism, and formation of metastable phases; insights from nanocrystalline TiO 2[END_REF]. We track the evolution of structure in this system using Tong and Tanaka's Θ order parameter [START_REF] Tong | Revealing Hidden Structural Order Controlling Both Fast and Slow Glassy Dynamics in Supercooled Liquids[END_REF], a local measure of packing inefficiency that vanishes in the limit of perfect FCC structure. We also track the corresponding evolution of the mobility µ (t), defined here for each particle as the magnitude of its inherent state displacement between times t and t + 0.1.

Test 1: Here we simply look for the emergence and growth of T ∞ -structure domains in individual realizations of the homogeneous experiment. Figs. 1(a) for melting and 1(e) for freezing (see also §B 4) both show a domain of T ∞ structure and dynamics which grows to span the system along the narrow horizontal axis. At this point, the domain can expand without gaining any more interfacial area, making this width an upper bound for the critical nucleus size, and we see the domain quickly grow to envelop the entire system. Test 2: Here we look for evidence of a waiting time for nucleation. In Fig. 1(a), we see that the time between forming a system-spanning nucleus and melting the entire system is of order t ∼ 10 0 , while it takes a time of order t ∼ 10 2 for this system-spanning nucleus to appear. For cooling (Fig. 1(e)) these time scales are t ∼ 10 2 and t ∼ 10 3 , respectively. This scale separation between the nucleation and domain growth times holds true across realizations ( §VI).

Two-state scenario: Our third and fourth tests involve slab experiments. To ensure that results from these are relevant to systems without an artificial slab, we first demonstrate that the structural evolution observed in the homogeneous experiments corresponds to a two-state scenario, such that at any time t, a fraction f of the system has local structure distributed according to equilibrium at the thermostat temperature T ∞ (the target state), while the remaining fraction 1-f retains its initial distribution of structure, corresponding to equilibrium at temperature T 0 (the initial state). In this scenario, the distribution P (ϕ, t) of a local structural indicator ϕ with initial-state distribution P 0 (ϕ) and target-state distribution P ∞ (ϕ) will evolve according to

P (ϕ, t) = (1 -f (t)) P 0 (ϕ) + f (t)P ∞ (ϕ) , (2)
interpolating between P 0 (ϕ) and P ∞ (ϕ) as the transformed fraction f grows with time. Given distributions P 0 (ϕ), P ∞ (ϕ), and P (ϕ, t) for a structural quantity ϕ, we can extract the transformed fraction f from P (ϕ, t) under the assumption that Eq. 2 holds (see §A 8). Alternatively, Eq. 2 implies that

1 -f (t) = ϕ ∞ -ϕ (t) ϕ ∞ -ϕ 0 , (3)
where ϕ 0 , ϕ ∞ , and ϕ (t) are the mean values of the quantity ϕ distributed according to P 0 (ϕ), P ∞ (ϕ), and P (ϕ, t), respectively, allowing us to calculate f from the evolving mean of ϕ.

If the two-state scenario holds, we should obtain the same value for the transformed fraction f , regardless of whether we use the distribution P (ϕ, t) or the mean ϕ (t) to calculate it, and regardless of which local structural variable ϕ we choose. We can therefore test the validity of the two-state scenario by calculating f from both P (ϕ, t) and ϕ for a variety of different structural quantities ϕ. In practice, we only calculate f from P (ϕ, t) for a single structural quantity ϕ, since this is more complicated to do ( §A 8). Though it is not a requirement of the two-state equilibration scenario, we also check if we can obtain f from P (ϕ, t) and ϕ when ϕ is a dynamical, rather than structural, quantity, since the coupling of structure and dynamics is a question of physical interest.

In Figs. 1(b), 1(c), 1(f), and 1(g), we see that the T 0 and T ∞ distributions of log 10 Θ and log 10 µ are wellseparated [63], making the validity of the two-state scenario immediately apparent. We calculate 1f from these distributions and plot the resulting curves against 1-f calculated from ϕ for ϕ corresponding to Θ, log 10 Θ, µ, log 10 µ, the potential energy V associated with each particle, and the inherent state potential energy V IS associated with each particle. As shown in Figs. 1(d) and 1(h), these different methods of calculating 1f agree well with one another, confirming validity of the twostate scenario and therefore the soundness of the slab experiment as a test for domain growth, and showing that structure and dynamics are coupled together in crystal freezing and melting.

Having validated the slab experiments, we calculate the (ensemble-averaged) Θ profile along the x-axis, shown as transparent curves in Figs. 2(c) and 2(d). We identify the interface at the boundary between a T ∞ structure slab and the surrounding T 0 structure region by fitting this profile with a logistic function

Θ (x) = Θ ∞ + Θ 0 2 + Θ ∞ -Θ 0 2 tanh x -x * 2λ , (4)
where the argument x = 2 |x| /L x accounts for changes in L x (the system size along x) due to the barostat and takes advantage of the x → -x symmetry of our slab geometry [START_REF]Symmetry-breaking sample-to-sample fluctuations are suppressed upon ensemble averaging in practice[END_REF]. The parameters Θ ∞ , Θ 0 , x * , and λ are fit at each snapshot of time t, representing, respectively, the plateau value of Θ within the slab of T ∞ structure, the plateau value of Θ outside this slab, the interface center, and the interface width. The order parameter Θ is invariant with respect to isotropic expansion, and is therefore insensitive to the transient effect of the pressure gradient at time t = 0, making it an especially suitable choice for our analysis.

Test 3: Here we check for domain growth in the slab experiment. As we see in Figs. 2(c) and 2(d), the slab of T ∞ structure indeed grows in both melting and freezing experiments. Figs. 2(e) and 2(f) show that after an initial transient of around 150 time units, which we associate with the dissipation of the pressure gradient present in the initial configuration, the displacement ∆x * = x * (t)x * (0) of the interface center increases at constant speed. Both this constant front speed and the two-state equilibration scenario accord with the description of domain growth within the Avrami formalism [START_REF] Fanfoni | The Johnson-Mehl-Avrami-Kohnogorov Model: A Brief Review[END_REF] for crystallization and melting.

Test 4 : In our final test we look for evidence that interfacial area is being minimized. In Fig. 1(a), we see a transition from the domain boundary being preferentially aligned along a crystal lattice vector to being aligned horizontally, while in Fig. 1(e), the domain boundary has no particular initial orientation (the surrounding structure is amorphous), but also quickly becomes horizontal. This is evidence of an interfacial energy penalty promoting a reduction of interfacial area. For the slab experiment, the initial width λ 0 of the interface has to do with the abil- ity of T 0 and T ∞ structure to constrain other structure nearby, but the two-state equilibration scenario means we can associate further growth in λ with roughening of the interface. In Figs. 2(g) and 2(h), after the initial pressure-driven transient (t ≥ 150), fluctuations in ∆λ = λ (t)λ 0 are comparable in magnitude to the size of the noise floor. This suppression of roughening is further evidence of an interfacial energy penalty promoting a reduction of interfacial area. In summary, the results of our four tests match our expectations for a nucleation and domain growth process, validating our methodology.

IV. GLASSES

By adjusting our methodology for calculating Θ ( §A 5) and redefining the mobility µ (t) to be the inherent state displacement of a given particle between times t and t + 100 (accounting for Mermin-Wagner fluctuations as described in §A 6), we can use Θ and µ to characterize the local structure and dynamics of our two-dimensional polydisperse glass.

Test 1: Here we look for growing T ∞ -structure domains in an individual realization of the homogeneous experiment ( §B 4). In Fig. 3(a), we see growing high-Θ and matching high-µ domains, representing clear evidence that heating is mediated by domain growth. We do not, however, see any visible low-Θ or low-µ domains during cooling in Fig. 3(e), suggesting a lack of domain growth in this case. This is a first and important difference from the crystal case.

Two-state scenario: Before applying our second test, which will benefit from the results of the slab experiments, we confirm the validity of the two-state scenario. We find good agreement between 1f calculated from the different structural quantities P (Θ), V , V IS , Θ, and log 10 Θ in Figs. 3(d) and 3(h), establishing the validity of the two-state scenario for both heating and cooling. These figures also highlight a major difference between heating and cooling: the coupling of structure and dynamics. For heating, this coupling is evident from the matching high-Θ and high-µ domains, and the good agreement between 1f calculated using any of our structural quantities with 1f calculated from P (log 10 µ) and µ (Fig. 3(d)) demonstrates this yet further. (We do see a deviation from this trend of 1f calculated from log 10 µ due to the rightwards shift of the P 0 (log 10 µ) mode with time in Fig. 3(c); see also §B 1.) For cooling, however, 1f calculated from µ or log 10 µ decays much faster than 1f calculated from structural quantities (Fig. 3(h)), and while P (Θ) evolves in accordance with the two-state scenario of Eq. 2 (Fig. 3(f); see also §A 8), the evolution of P (log 10 µ, t), comprising a continuous shift of the initial distribution P 0 (log 10 µ) towards lower mobility (Fig. 3(g)), does not. We therefore see that structure and dynamics are decoupled during cooling.

Test 2: Our second test is to look for evidence of a nucleation time. With no evidence of domain growth during cooling, we focus here on heating. We see in Fig. 3(a) that T high -structure domains first appear at time t ≲ 10 4 . In a classical nucleation theory scenario, small domains are suppressed, so the nucleation time should be much larger than the time scale of individual rearrangements (successful nucleation should require many attempts), as was the case for the crystal (Fig. 1). We can obtain an estimate of the rearrangement time scale from the time scale for the growth of T ∞ -structure domains, which is a lower bound for the rearrangement time scale within domains of T low -structure. In our slab experiments, we find (Fig. 4(e)) that T high -structure fronts advance with speed 7.85 × 10 -4 into the bulk (neglecting deviations of the system size from L x = 100) during heating, corresponding to 1.27 × 10 3 time units per mean particle diameter. Fronts advancing at this rate into the bulk in both the y and -y directions would take 6.36 × 10 4 time units to cover 100 distance units to span the system, a result consistent with Figs. T ∞ = T high -structure slabs indeed grow in heating experiments, with the center of the slab interface, x * , advancing with constant speed up to the system size limit (Figs. 4(c) and 4(e)). In the case of cooling, however, we see that slabs of T ∞ = T low structure do not grow at all (Figs. 4(d) and 4(f)). This shows that the reason domains of T low structure are not visible in Fig. 3(e) is that such domains do not grow in the first place. In a classical nucleation theory scenario, the competition between bulk and interfacial free energies would apply to both heating and cooling, with the free energy difference between T low and T high structure in fact being larger during cooling given the lower thermostat temperature T ∞ = T low . This asymmetry and lack of domain growth during cooling is therefore a strong piece of evidence against this scenario.

Test 4: Here we investigate whether the size of the interface is being minimized. The high-mobility domains visible in Fig. 3(a) clearly have rough edges, and we see no alignment along the horizontal axis to minimize the length of the interface. Indeed, spanning the system appears to have no impact on the growth of domains in Fig. 3(a), unlike in the crystal case, where spanning the system was associated with an acceleration of domain growth and alignment of the domain boundaries with the short axis of the system (Figs. 1(a) and 1(e)). The growth of the interface width λ in slab experiments provides further evidence against a thermodynamic drive to minimize the size of the interface. After an early-time transient (t ≲ 10 3 ) of diffusive growth (∆λ ∼ t 1/2 ), we see that that the growth in λ during heating (Fig. 4(g)) is super-diffusive until the time t ≈ 12000 at which the opposite sides of the slab meet at the periodic boundary. This indicates significant interfacial roughening, faster even than the naïve diffusive expectation for an unsuppressed, freely-diffusing interface.

In summary, the glass behaves very differently from the crystal upon equilibration. All four tests reveal discrepancies between observations and expectations based on a nucleation and growth picture.

V. PLAQUETTE MODELS

Here we show that the very same findings that are inconsistent with a nucleation and growth picture are explained naturally with dynamical facilitation. We mean dynamical facilitation in its broadest sense as the selfpropagation of mobility, independent of the mechanisms through which this propagation is achieved. To this end, we apply our four tests to the TPM and SPPM. These models have been shown to display a diverging (at zero temperature) point-to-set length, but the interfacial tension between amorphous metastable states in these systems is zero and their equilibrium dynamics is completely driven by dynamical facilitation, not by RFOT dynamics [START_REF] Jack | Caging and mosaic length scales in plaquette spin models of glasses[END_REF]. Here, the initial interface width λ0 has values λ0 = 0.00428 and 0.00414 in (e) and (f) respectively. The dashed line in (e) is a linear fit, while the dashed line in (g) represents diffusive growth ∼ t 1 2 , showing that λ is growing super-diffusively. The dotted curve in (c), and vertical lines in (e) and (g) correspond to time t = 12000, above which we judge the left and right boundaries of the slab to interact with one another at the x = ±1 periodic boundary of the system.

A. Model details

The TPM and SPPM are, respectively, two-and threedimensional lattice models comprising spins σ i ∈ {-1, 1} (i being the lattice coordinate (i, j) or (i, j, k)) interacting via Hamiltonians

H = - J 2 i p i , (5)
where J is the interaction energy and the plaquette p i ∈ {-1, 1} is the product of the triplet of spins s i-1,j , s i-1,j+1 , and s i,j in the TPM and the quintuplet of spins s i-1,j-1,k+1 , s i-1,j,k+1 , s i,j-1,k+1 , s i,j,k+1 , and s i,j,k in the SPPM. We evolve these systems according to Glauber dynamics [START_REF] Glauber | Time-Dependent Statistics of the Ising Model[END_REF]: nodes are randomly selected, one at a time, to undergo a trial flip σ i → -σ i with probability

P (flip) = 1 1 + e β∆E (6)
of acceptance, where ∆E is the change in energy that would result from this spin flip and β = 1/k B T for temperature T . Given periodic boundary conditions and a number of particles along at least one axis that is a power of two, there is a one-to-one correspondence between spins and plaquettes [START_REF] Garrahan | Glassiness and constrained dynamics of a short-range nondisordered spin model[END_REF], and it is in the plaquette representation that the TPM and SPPM are seen to be KCMs, with non-interacting Hamiltonians (Eq. 5) and facilitation resulting from the kinetic constraint that trial moves correspond to a triplet of plaquettes p i,j , p i+1,j-1 , and p i+1,j flipping simultaneously in the TPM and a quintuplet of simultaneously flipping plaquettes p i,j,k , p i,j,k-1 , p i,j+1,k-1 , p i+1,j,k-1 , and p i+1,j+1,k-1 in the SPPM. We take the indices i, j, and k to vary along the x, y, and z axes, respectively, with bond length b separating adjacent nodes along an axis, and adopt b, J, and J/k B as our length, energy, and temperature units. We assume a constant rate of N trial flips per unit of time, where N is the number of lattice points in the system. For both the TPM and SPPM, we simulate systems with 128 lattice points along each axis. We once again conduct homogeneous and slab experiments. For the slab experiments, we define our slab to be the set of points with axis coordinate α satisfying |α| < 12, where α is the x-coordinate in the TPM case and the z-coordinate in the SPPM case. Since plaquettes do not interact in these systems (Eq. 5), we equilibrate the slab and non-slab regions by randomly sampling the excitation state of each plaquette site with equilibrium probability 1/ 1 + e 1/T of being excited, where T = T ∞ in the slab and T 0 outside the slab. Because the TPM is not invariant under x → -x and the SPPM is not invariant under z → -z, we treat the interfaces on the left (α < 0) and right (α > 0) sides of the slab separately.

B. Test results

Test 1: Here we look for growing domains in an individual realization of the slab experiment. As we see in Fig. 5, there is domain growth in both the TPM and SPPM during heating but not during cooling. When present, the domains grow with a characteristic triangular shape set by the kinetic constraint, highlighting the fact that this domain growth is due to dynamical facilitation. This asymmetry of having domain growth during heating but not during cooling is consistent with the results of §IV for the polydisperse glass, and different from the results of §III for a genuine first-order phase transition in the monodisperse system.

Test 2: Here we look for evidence of a nucleation time. Ironically, given that there is by design no nucleation in these systems (the plaquettes do not interact, so there can be no interfacial tension), the high energy barrier to thermal activation in an excitation-free system means that we do see a significant waiting time of order t ∼ 10 4 before the emergence of domains, which grow to envelop the system on a time of order t ∼ 10 3 . This waiting time arises due to the need to escape the large potential energy basin of the defect-free initial system in order for domain growth to begin.

Two-state scenario: The TPM and SPPM comprise non-interacting plaquettes in one of two states, so we can safely assume that the results of slab experiments will relate to the behavior seen in the homogeneous experiment. There is therefore no need to test the two-state scenario as we did for the glass and crystal.

Test 3: Here we track domain growth quantitatively in the slab experiment. As we show in Figs. 6(e), 6(f), 6(i), and 6(j), the result of the third test is that slabs in the TPM and SPPM indeed grow upon heating, with the displacement of the interface center from the slab center increasing at a constant rate. The speed at which the center of the left interface advances is different from the corresponding speed for the right interface, as expected given the asymmetric kinetic constraints. By contrast, we see in Figs. 6(g), 6(h), 6(k), and 6(l) that T ∞structure slabs in the TPM and SPPM do not grow during cooling, but sometimes even shrink slightly. These results are, once again, consistent with the results of the polydisperse glass.

Test 4: Here we look for evidence of constraints on the roughness of growing domains. During heating, the growing TPM domains in Fig. 5(a) have rough boundaries, but we see little roughness in the growing SPPM domains in Fig. 5

(b).

There is also no evidence of accelerated domain growth once the system spans the box width, as was seen in Figs. 1(a) and 1(e) for the monodisperse crystal. This is confirmed by the growth of the interface width λ in the slab experiments. As shown in Fig. 6(m), λ grows approximately diffusively during heating in the TPM, while, as seen in Fig. 6(n), the growth of λ in the heated SPPM is strongly sub-diffusive. Despite the similar designs of the TPM and SPPM, the roughness of growing domains during heating has very different dynamics in the two cases. We note that even within a given model, the power law for the growth of λ for the left interface differs from the power law for the right interface. This shows that the asymmetric constraint is affecting the rate of roughening, and highlights the sensitivity of roughening to the mechanism underlying facilitation. The fact that neither the TPM nor the SPPM see the super-diffusive growth of λ seen in the heated polydisperse glass is consistent with this.

We see that the main phenomena associated with the equilibration dynamics of glasses-Avrami-like domain growth during heating and an absence of domains during cooling-can be obtained purely via dynamical facilitation, even with trivial thermodynamics, as in the KCMs studied here, suggesting that these two-state scenario be- haviors are generic to dynamically-facilitated systems.

The one phenomenon seen in the glass but not recovered in the KCMs is super-diffusive roughening during heating, but this roughening behavior is highly sensitive to the specific details of the system, varying not only between KCMs, but even between the left and right boundaries of a single slab.

C. Phenomenological model

From the shared phenomenology of the polydisperse glass and the KCMs, we can build a minimal, phenomenological model of slab growth, of excited plaquettes within bins of constant x (TPM) or z (SPPM) at times (e) t = 0, 10000, . . . , 100000, (f) t = 0, 1000, . . . , 10000, (g) t = 0, 10 3 , 10 4 , 10 5 , and 10 6 , and (h) t = 0, 10 1 , 10 2 , 10 3 , 10 4 and 10 5 (larger time, lighter color). Here, c is the unnormalized concentration, and c0 and c∞ are the equilibrium concentrations 1/ 1 + e 1/T for T = T0 and T∞, respectively. Transparent curves show raw data, solid curves are fits to Eq. 4. The solid curves fully overlap with the transparent curves, apart from in (g). Bottom rows: change in the fit parameters (i)-(l) x * and (m)-(p) λ relative to their values in the fit to the t = 0 profile. The orange curves correspond to the slab interface initially at negative x, the blue curves to the slab interface initially at positive x. The dashed lines in (i)-(l) are linear fits, the dashed lines in (m)-(p) are power-law fits. There is no blue curve in (o) because it is below the noise floor. The dotted curve in (e), and vertical lines in (i) and (m) correspond to time t = 50000, above which we judge the left and right boundaries of the slab to interact with one another.

∂ t c = (c ∞ -c) /τ + s |∂ x c| + D∂ 2 x c, (7)
where c is the excitation state of the system (e.g. Θ or the plaquette concentration c), τ is the local relaxation time of the system, s is the speed of mobility propagation (the absolute value ensuring that transport is always directed from high to low mobility), and D is a diffusion constant accounting for the roughening of the system. This allows us to discuss the observed dynamics of these systems in terms of specific dynamical ingredients. For simplicity, we focus on slab growth in the TPM, identifying the contributions to the rate of change ∂ t c of the plaquette concentration c, and noting system-dependent considerations when they arise.

The first contribution we consider is the mobility propagation s |∂ x c|. In the ultrastable glass heating literature, experiments [START_REF] Rodríguez-Tinoco | Transformation kinetics of vapor-deposited thin film organic glasses: The role of stability and molecular packing anisotropy[END_REF][START_REF] Ràfols-Ribé | The role of thermodynamic stability in the characteristics of the devitrification front of vapour-deposited glasses of toluene[END_REF] and a KCM study [START_REF] Léonard | Macroscopic facilitation of glassy relaxation kinetics: Ultrastable glass films with frontlike thermal response[END_REF]

find front speeds s = l s /τ ∞ (T ∞ )
γ , where γ ≈ 1. Compar-ing Figs. 4(c) and 6(f) with Figs. 4(i) and 6(j), we see that in both the polydisperse glass and the TPM, the velocity of the interface center remains constant even after significant relaxation of the non-slab region towards equilibrium, so we take l s to be a constant. The specific form of the relaxation time τ is systemdependent, but any sufficiently fragile system would yield the speeds s ≈ 0 for low T ∞ and s > 0 for high T ∞ seen in Figs. 4 and6. From the equilibrium behavior τ ∼ e 1/ log(3)T 2 [START_REF] Jack | Caging and mosaic length scales in plaquette spin models of glasses[END_REF] of the TPM, we obtain

τ = τ 0 exp 1/ log 3 ((1 -a) T ∞ + aT f ) 2 (8)
by modifying the Tool, Narayanaswami, and Moynihan (TNM) model [START_REF] Moynihan | Dependence of the Fictive Temperature of Glass on Cooling Rate[END_REF] to yield the TPM relaxation in equilibrium and to work in the limit of vanishing fictive temperature

T f = 1/ log c -1 -1 .
Here, τ 0 is a constant time scale and a is a fit parameter (see §B 3). Finally, the rate of roughening is system-dependent, as we have discussed. In the TPM, this roughening yields an approximately diffusive growth of the slab width λ ∼ t 1/2 , which can be obtained using a diffusive term D∂ 2

x c. From Figs. 4 and6, we see that when the interface speed is small, so is the diffusion rate, so we assume proportionality between the two, D = l 2 D /τ ∞ , where l D is a constant.

Comparing Figs. 4 and7, we see that this phenomenological model of relaxation dynamics based on dynamic facilitation captures the key features of slab evolution, both during heating and during cooling, justifying the separation of dynamics into distinct local relaxation, constant-speed drive, and roughening pieces.

VI. DISCUSSION

We have investigated whether the relaxation of a glass towards equilibrium during heating and cooling proceeds via nucleation and domain growth in the sense of classical nucleation theory. By applying the same four tests to the heating and cooling of a three-dimensional monodisperse crystal, a two-dimensional polydisperse glass, and the TPM and SPPM, we find that the equilibration dynamics of our glass does not resemble that of the melting or freezing crystal, but instead resembles that of the dynamically-facilitated plaquette models in every discriminating aspect. This is made immediately apparent in Fig. 8, which shows the evolution of 1-f calculated from Eq. 3 applied to Θ (8(a)-8(d)) or the plaquette concentration c (8(e)-8(h)) for ten individual homogeneous experiment trajectories (c.f. the ensemble-averaged data in Figs. 1, 3, and5). For crystal melting and freezing in our monodisperse system, we see a nearly discontinuous change of 1f in Figs. 8(a Dynamics of the phenomenological model governed by Eq. 7 for constants ls = 0.18 and lD = 0.25 given an artificial target-state slab at time t = 0. Heating (left column) is performed between temperatures T0 = 0.0 and 0.3 using τ0 = 0.549 and a = 0.253 and cooling (right column) is performed between temperatures T0 = 0.3 and 0.2 using τ0 = 0.705 and a = 1.00. Parameters are chosen so as to yield behavior similar to the TPM, establishing the reasonableness of this model, but no explicit fit to the TPM is performed. (a, b): normalized excitation state c defined as in Fig. 6(e) at times (a) t = 0, 1000, . . . , 10000 and (b) t = 0, 10 3 , 10 4 , 10 5 , and 10 6 (larger time, lighter color). Transparent curves correspond to numerically-solved Eq. 7, solid curves to fits to the x >= 0 part of this solution, noting that c (-x) = c (x). (c-f): change, relative to its initial value, of the interface center position x * (c, d), and width λ (e, f) obtained when fitting the profiles to Eq. 4. The dashed lines in (c) and (d) are linear fits, while those dashed line in (e) and (f) represent diffusive growth ∼ t of metastable states in Fig. 8(b). By contrast, the evolution of 1f in both our polydisperse glass and the two plaquette models (Figs. 8(c)-8(h)) is much smoother and more gradual, reflecting the absence of a nucleation time in these systems. For cooling (Figs. 8(d), 8(f), and 8(h)), the lack of variation between individual trajectories reflects the lack of domain growth in these systems, which results in better averaging within individual realizations.

Fig. 8 highlights the fact that while the absence of a nucleation time and the asymmetry between heating being driven by domain growth and cooling being a homogeneous process is not consistent with a classical nucleation theory scenario, the equilibration dynamics of our glass is in fact generic to dynamically-facilitated systems. Individual trajectories of the glass resemble those of the (e-h): Evolution of 1f calculated from Eq. 3 using the concentration c of excited plaquettes for ten individual realizations in the same regimes as those shown in Fig. 5. plaquette models and not the crystal because dynamical facilitation, not nucleation and domain growth, governs the equilibration dynamics of glasses.

A. Nucleation and Domain growth

While we do find that the equilibration of our glass during heating is mediated by domain growth (Figs. 3,4, and 8), we do not find evidence of nucleation for either heating or cooling. We find no evidence of a nucleation time for the emergence of T high domains during heating, which appear on a similar time scale to the time scale over which the front of a growing T high domain advances by one mean particle diameter. We also find no evidence of any thermodynamic drive to minimize the total interface between the domains of T low structure and T high struc-ture during this heating process. Unlike in the crystal case (Figs. 1(a) and 1(e)), there is no alignment of the growing T high domain in Fig. 3(a) along the horizontal axis to minimize its size, and in slab experiments, the slab boundary exhibits super-diffusive roughening. Cooling ultimately provides the most striking evidence against a drive to minimize the interface between domains of T high and T low structure, since the absence of domain growth combined with the validity of the two-state equilibration scenario implies that small domains of T low structure surrounded by T high structure proliferate during the cooling process.

This result concerns only nucleation in the thermodynamic sense of classical nucleation theory. First-order transitions in systems coupled to a reference configuration [START_REF] Berthier | Evidence for a Disordered Critical Point in a Glass-Forming Liquid[END_REF] and dynamical phase transitions [START_REF] Hedges | Dynamic Order-Disorder in Atomistic Models of Structural Glass Formers[END_REF] are found in KCMs, and an explanation of relaxation dynamics based on those transitions may be possible [START_REF] Jack | The melting of stable glasses is governed by nucleation-and-growth dynamics[END_REF][START_REF] Turner | Overlap and activity glass transitions in plaquette spin models with hierarchical dynamics[END_REF]. As we have shown, however, this explanation cannot involve traditional nucleation and growth, and we believe that it would ultimately boil down to a complementary view on how dynamical facilitation influences relaxation dynamics. Note, moreover, that our results also do not contradict the existence of a finite point-to-set length ξ PTS , nor the possibility of an RFOT mosaic. In fact, the asymmetry between heating and cooling, including the decoupling of structure and mobility during the latter process, can be derived via the microscopic RFOT theory [START_REF] Lubchenko | Theory of Structural Glasses and Supercooled Liquids[END_REF][START_REF] Wolynes | Spatiotemporal structures in aging and rejuvenating glasses[END_REF]. As shown in §B 2, individual rearrangements perturb the surrounding structure over a length scale ξ perturb ∼ 10 0 similar to the point-to-set length ξ PTS ∼ 10 0 for our temperature range [START_REF] Berthier | Zero-temperature glass transition in two dimensions[END_REF], with the typical size of perturbations to Θ comparable to the size of Θ itself. This raises the possibility that individual rearrangements create individual domains of T low structure within a mosaic as the system cools.

B. Dynamical facilitation

We have established that dynamical facilitation drives the domain growth observed when heating the polydisperse glass. Despite very different mechanisms for dynamical facilitation, the interface center x * evolves in the same manner in both the polydisperse glass and our two plaquette models. In all cases its behavior is well-captured by a constant-speed driving term s |∂ x c| in Eq. 7. In recent work, Herrero, Ediger and Berthier [START_REF] Herrero | Front propagation in ultrastable glasses is dynamically heterogeneous[END_REF] conducted very similar slab experiments to ours using the same polydisperse system and find, suggestively, that it is individual rearrangements near the slab boundary that advance the transformation front. They propose a separation of time scales between the mobile T high -structure slab and the immobile T low -structure background that would explain the constant speed of the transformation front, given a constant probability of rearrangement for particles at the leading edge of the front. The validity of the two-state equilibration scenario allows us to relax the requirement of a separation of time scales, which in turn allows us to explain the constant front speed seen in the slab heating experiments at times t ≥ 10 3 despite significant relaxation of the non-slab region (Fig. 4(e)). As this region relaxes, domains of T high structure grow which contribute to raising the height of the T low -region plateau in Fig. 4(c) (parameter Θ 0 in Eq. 4), rather than to advancing the front. Instead, the front advances via the transformation of T low -structure adjacent to the front facilitated by rearrangements near the interface. The same scenario explains observations of the same behavior in the TPM (Fig. 6(i)).

The interface width λ, for its part, provides systemspecific information. The initial interface widths λ 0 relate to the ability of local structure to constrain structure nearby, hence λ 0 in the crystal (Fig. 2) being larger than in the glass (Fig. 4), while λ 0 ≈ 0 in the non-interacting plaquette systems (Fig. 4). Given the two-state equilibration scenario, the evolution of roughness is captured by the evolution of λ in the ensemble-and spatiallyaveraged profile Θ (x) (Eq. 4) in slab experiments. As we have seen, this is strongly system-dependent, to the point where the widths of the interfaces on the left and right sides of TPM and SPPM slabs grow as power laws λ ∼ t ν with different power law exponents ν. Among our glass and two plaquette models, we find exponents ν spanning values well above (Fig. 3(g)) to well below (Fig. 5(n)) the naïve diffusive exponent ν = 1/2 for roughness arising from random fluctuations in the local front speed. It was shown in [START_REF] Herrero | Front propagation in ultrastable glasses is dynamically heterogeneous[END_REF] that slab roughness in the slab-experiment heating of very cold systems at very long times reaches a steady state with characteristic length scale set by the dynamic correlation length ξ 4 [START_REF] Lačević | Spatially heterogeneous dynamics investigated via a time-dependent four-point density correlation function[END_REF]. Future work may similarly uncover the physical origin of the power-law exponent ν.

C. Technical advances

Our investigations of nucleation and domain growth and dynamical facilitation were aided by two key technical advances with implications of their own. The first was the use of Eqs. 2 and 3 to confirm the two-state equilibration scenario. While, on the basis of heat capacity measurements, this scenario has long been expected to hold for the heating of ultrastable glasses [START_REF] Vila-Costa | Nucleation and Growth of the Supercooled Liquid Phase Control Glass Transition in Bulk Ultrastable Glasses[END_REF][START_REF] Kearns | Observation of low heat capacities for vapor-deposited glasses of indomethacin as determined by AC nanocalorimetry[END_REF][START_REF] Rodríguez-Tinoco | Evaluation of Growth Front Velocity in Ultrastable Glasses of Indomethacin over a Wide Temperature Interval[END_REF][START_REF] Rodríguez-Tinoco | Transformation kinetics of vapor-deposited thin film organic glasses: The role of stability and molecular packing anisotropy[END_REF][START_REF] Ràfols-Ribé | The role of thermodynamic stability in the characteristics of the devitrification front of vapour-deposited glasses of toluene[END_REF][START_REF] Cubeta | Communication: Surface-facilitated softening of ordinary and vapor-deposited glasses[END_REF][START_REF] Rodríguez-Tinoco | Surface-Bulk Interplay in Vapor-Deposited Glasses: Crossover Length and the Origin of Front Transformation[END_REF], our work here contains the most direct confirmation of the two-state scenario to date, and confirms it for cooling as well as heating, despite the lack of domain growth in the former. This has encouraging implications for the modelling of glasses out of equilibrium, showing that local structure after a rearrangement can simply be sampled directly and independently from the equilibrium distribution of structure at the target equilibrium temperature, as for instance assumed in trap models [START_REF] Monthus | Models of traps and glass phenomenology[END_REF].

A second technical advance of our work is the introduction of the phenomenological model of slab growth, Eq. 7. This model shows how slab growth during heat-ing can be reconciled with homogeneous relaxation during cooling given the same set of dynamical ingredients. By identifying the distinct components of the dynamics of equilibrating glasses, this model should also aid the development of future models for the equilibration of glasses in non-slab geometries. ultimately yielding

A ≈ 1.3145271918 (A8) σ min ≈ 0.7244461244 (A9) σ max ≈ 1.6138530488 (A10)
for σ = 1 and c σ = 0.23.

Preparation protocols

To homogeneously equilibrate to low temperatures, we conduct Monte Carlo simulations at constant volume and temperature until steady state, randomly selecting particles and conducting trial displacements with components randomly sampled from the interval (-0.05, 0.05). For polydisperse systems, we implement a swap Monte Carlo methodology originally developed by Misaki Ozawa for [START_REF] Berthier | Zero-temperature glass transition in two dimensions[END_REF]. With probability 0.2, a trial move corresponds to a trial swap instead of a trial displacement. If a particle with size σ i is selected for a swap move, a second particle with size σ j satisfying |σ jσ i | < 0.2 is randomly selected as its swap partner. We start our Monte Carlo simulations from a configuration with random initial positions, except when initializing a crystal in the monodisperse system, in which case we start from a perfect FCC crystal in order to accelerate the equilibration process.

When preparing polydisperse systems containing an artificial target-state slab, we use swap Monte Carlo to equilibrate the system to the target temperature T ∞ . We then designate a slab-shaped region, corresponding to points with x-coordinate satisfying |x| > a for some a > 0, to be outside the target-state slab. We choose a new, off-target temperature T 0 for this region, and find the number density for a homogeneous system at equilibrium at temperature T 0 and with virial pressure matching that of the target equilibrium system. We stretch or compress the off-slab region along the x-axis such that it matches this number density, thus reducing the pressure gradient that would otherwise exist between the slab and the non-slab region. MD simulations are then run at the temperature and pressure corresponding to the target equilibrium system.

When preparing monodisperse systems containing an artificial target-state slab, we choose crystal and liquid temperatures that are similar in magnitude, such that the pressure gradient is less pronounced. It is not possible to stretch or compress the off-slab region along the x-axis to reduce the pressure gradient, as we do in the polydisperse case, as this would prevent us from guaranteeing system dimensions compatible with a perfect FCC crystal. We first homogeneously equilibrate the system to the target crystal state, then conduct Monte Carlo simulations at high temperature, sampling only from the liquid region. Once the liquid region has melted, we conduct Monte Carlo simulations at the target, lower temperature for the liquid region, sampling only from this region. MD simulations in this case are run at constant temperature and a constant pressure chosen such that the virial pressure at equilibrium matches the mean virial pressure in the initial configuration.

Molecular dynamics

We conduct our MD simulations in LAMMPS [START_REF] Plimpton | Fast Parallel Algorithms for Short-Range Molecular Dynamics[END_REF] (23 Jun 2022 -Update 2) at constant temperature and pressure implemented via LAMMPS' Nosé-Hoover thermostat (fix nvt) with damping parameter Tdamp = 1 and Berendsen barostat (fix press/berendsen) with damping parameter Pdamp = 1. If a system has been initialized to temperature T 0 using Monte Carlo, we find the virial contribution to the pressure at this temperature, and choose a pressure such that this contribution is maintained once the temperature changed to the target temperature T ∞ . We do this so as to minimize changes in the system size upon changing the temperature.

Initial and steady-state distributions

Because the start-up of the thermostat and barostat perturbs the system at early times, we must wait until a time t 0 > 0 before we can measure initial distributions for the homogeneous experiments. We use t 0 = 55 for simulations of the monodisperse system and t 0 = 100 for the polydisperse system. For homogeneous melting experiments in the monodisperse system, equilibrium distributions are obtained from t = 600 data, restricting to realizations with bulk average Θ (t = 600) > 0.091. For the homogeneous freezing experiments, "equilibrium" distributions are obtained from t = 10000 data. This does not correspond to equilibrium at the target temperature T low , but rather a collection of metastable states that each realization is trapped in at that time (see §VI). For the polydisperse system, equilibrium distributions at temperature T ∞ = T high are obtained from homogeneous heating experiment data for t ≥ 80000, while equilibrium distributions at temperature T ∞ = T low are obtained from simulations at that temperature in the micro-canonical ensemble.

Θ

The main structural indicator used in this work is Tong and Tanaka's Θ order parameter [START_REF] Tong | Revealing Hidden Structural Order Controlling Both Fast and Slow Glassy Dynamics in Supercooled Liquids[END_REF], which measures how inefficiently packed the structure local to each particle is. For each particle i in the system, we iterate over pairs of nearest neighbors j and k such that j and k are also nearest neighbors of each other. We consider them perfectly efficiently packed if the pair potentials V ij , V jk , and V ik are all equal. If these particles have corresponding sizes σ i , σ j and σ k , this means their centers form a triangle with sides of length l ij = 1 2 (σ i + σ j ), FIG. 9. Subplot (a): pair distribution g (r) for the polydisperse system at temperature T = 0.14, calculated separately for particles at the origin of size σ within each of ten logarithmically-spaced bins, plotted as separate curves (larger σ, darker curve). This pair distribution is calculated using the same normalization as the usual pair distribution g (r), but replacing the raw interparticle separation r with r, hence why g (r) ↛ 1 as r → ∞. Subplot (b): r-position rNN of the trough following the first nearest neighbor peak of g (r) as a function of the diameter σ of the particle at the origin for temperatures T = 0.06 (blue crosses) and T = 0.14 (orange circles). The dashed line is the linear fit yielding Eq. A14.

l ik = 1 2 (σ i + σ k ) and l jk = 1 2 (σ j + σ k ). At i, this trian- gle makes an angle θ * ijk such that cos θ * ijk = l 2 ij + l 2 ik -l 2 jk 2l ij l ik . (A11)
Letting θ ijk be defined such that

cos θ ijk = r 2 ij + r 2 ik -r 2 jk 2r ij r ik , (A12)
where r ij is the distance separating particles i and j, with r jk and r ik defined similarly, we define Θ at i to be

Θ i = θ ijk -θ * ijk j,k , (A13)
where the average is over nearest-neighbor pairs j and k that are nearest neighbors of i as well as each other. The difference in angle θ ijkθ * ijk is calculated modulo π. For the purposes of calculating Θ i in our twodimensional polydisperse system, we consider i and j to be nearest neighbors if they have normalized separation rij < rNN (σ i ), where rNN (σ i ) = 2.025 -0.515σ i .

(A14)

This corresponds to the first minimum after the nearestneighbor peak of g (r), where g (r) is proportional to the probability of finding a particle at normalized separation r given a particle at the origin with diameter σ, as shown in Fig. 9. For our three-dimensional monodisperse system, the corresponding quantity is r = 1.3. When a particle i does not have any pair of nearest neighbors j and k that are nearest neighbors of each other, we leave its Θ i undefined and exclude it when calculating averaged quantities. V FIG. 10. Evolution in the canonical ensemble of the mean potential energy of an initially perfect HCP crystal with 22 particles along each basis vector. Curves correspond to temperature T = 1.5, 1.6, . . . , 5.0 (lighter color, higher temperature). Time is in units of number of trial Monte Carlo steps per particle.

Mermin-Wagner fluctuations

Our two-dimensional glass is influenced by Mermin-Wagner fluctuations [START_REF] Flenner | Fundamental differences between glassy dynamics in two and three dimensions[END_REF][START_REF] Vivek | Long-wavelength fluctuations and the glass transition in two dimensions and three dimensions[END_REF][START_REF] Illing | Mermin-Wagner fluctuations in 2D amorphous solids[END_REF]. We account for this by subtracting from any displacement vector ∆p i of a particle i the mean displacement vector ⟨∆p j ⟩ j averaged over neighboring particles j ̸ = i within a distance ξ of particle i. We choose the value of ξ to be that appropriate to equilibrium at temperature T ∞ , determined as in [START_REF] Chacko | Elastoplasticity Mediates Dynamical Heterogeneity Below the Mode Coupling Temperature[END_REF].

FCC Ground State

To determine that the ground state of the monodisperse system is an FCC crystal, we reason that the ground state of an isotropic repulsive pair interaction decreasing monotonically with distance must correspond to the distance-maximizing FCC or HCP (hexagonal close packing) structures. As shown in Fig. 10, HCP packings at temperature T < 2.7 (thermostat implemented using Monte Carlo trial steps at constant volume) transiently increase in potential energy before relaxing into what we find to be an FCC crystal, implying that the FCC state is more stable. (For T ≥ 2.7, the packings melt into a liquid instead.)

Interpolating distributions

In the main text, we calculate the transformed fraction f from the evolving distribution P (ϕ, t) for different quantities ϕ characterizing the local structure or mobility of the system, plotting them as dashed curves in Figs. 1(d), 1(h), 3(d), and 3(h). We achieve this using one of two approaches, depending on the behavior of the growing P ∞ (ϕ) mode. When fluctuations in the distribution near the P ∞ (ϕ) mode are not too severe, we use where ϕ * = arg max ϕ P ∞ (ϕ). This is the case for 1f calculated from P (log 10 Θ) and P (log 10 µ) in Fig. 1(d) (dashed red curve and dashed green curve, respectively), and for 1-f calculated from P (Θ) in Figs. 3(d) and 3(h) (dashed red curves). However, where the growing P ∞ (ϕ) mode is sensitive to noise, or where it shifts slightly as it grows (Fig. 3 We also calculate f from Eq. 3, using the averages ϕ of various quantities ϕ and plotting them as solid or dotted curves in Figs. 1(d), 1(h), 3(d), and 3(h). To demonstrate that our different methods of calculating f , (c): untransformed fraction 1f for f calculated from (dashed red) P (log 10 Θ) and (dashed green) P (log 10 µ), plotted against the mean (solid orange) potential energy V , (solid red) Θ, (dotted red) log 10 Θ, (solid green) µ, and (dotted green) log 10 µ, affinely transformed according to Eq. 3 with t0 = 2.0. On these short time scales, start-up of the barostat leads to erratic behavior of the inherent state energy, which is therefore not shown.

f (t) = P (ϕ * , t) -P 0 (ϕ * ) P ∞ (ϕ * ) -P 0 (ϕ * ) , (A15)
which collapse onto one another, successfully yield the interpolation parameter f in Eq. 2, we compare in Fig. 11 the evolutions P (ϕ, t) shown in Figs. 1(d), 1(h), 3(d), and 3(h) to the evolution obtained from Eq. 2 where P 0 (ϕ) and P ∞ (ϕ) are taken as input from simulation data and f is obtained from Eq. 3 using ϕ = log 10 Θ (for the crystal) or Θ (for the glass).

tion of the decay of the mean of log 10 µ from that of the fraction f of un-transformed structure in Fig. 3(d).

Here, we show in Fig. 12 that slightly increasing the temperature at which we heat our crystal (T high = 1.40, as compared to T high = 1.35 for melting in Fig. 1) leads to large rightwards shifts of the P 0 (log 10 Θ) and P 0 (log 10 µ) peaks as the melting proceeds (Figs. 12(a) and 12(b)), with corresponding deviations of 1f calculated from Eq. 3 using the mean values of log 10 Θ and log 10 µ from the correct value of 1-f obtained from the heights of the P ∞ (log 10 Θ) and P ∞ (log 10 µ) peaks (Fig. 12(c)). Inspection of individual trajectories (not shown) suggests that this deviation is due to a large defect density, with pairs of defects imposing strain on crystalline domains trapped between them (note the visible interface of length ≈ 4 particle diameters at the edge of the growing domain in Fig. 1(a), top row).

By contrast, we find that increasing T high to 0.17 in the homogeneous heating of the polydisperse glass, far above the onset temperature [START_REF] Scalliet | Thirty Milliseconds in the Life of a Supercooled Liquid[END_REF], does not lead to a perturbation of the un-transformed structure and resultant violation of the two-state equilibration scenario beyond that seen in Fig. 3 for log 10 µ. This is likely the result of the highly localized nature of facilitation at low temperatures, as seen in Fig. 3(d) of [START_REF] Chacko | Elastoplasticity Mediates Dynamical Heterogeneity Below the Mode Coupling Temperature[END_REF], and the fact that domains of T high and T low do not noticeably statically interact, as we establish in the main text.

Structural perturbation due to rearrangements

As the glass evolves during cooling in Fig. 3(e), it becomes possible to distinguish individual rearrangement events from the characteristic quadrupolar displacement fields they induce, as seen in the snapshots at times t ≥ 10 3 . In Fig. 13, we leverage this to identify the effect of such individual rearrangements on the local structure, identifying rearrangements with particles that have mobility µ > 0.1 and have a larger mobility than any other particle within 8 mean particle diameters of distance. We see that changes in Θ near the rearrangement are of the same magnitude ∼ 10 -1 as Θ itself (Figs. 3(b) and 3(b)) up to a distance r ≈ 5 from the rearrangement, decaying exponentially with decay length ξ perturb = 2.34.

Local relaxation in the TPM

Appropriate values for parameters τ 0 and a in Eq. 8 are obtained by fitting Eq. 7 for a homogeneous system to the decay of c (x = ±64) towards equilibrium in Figs. 6(e) and 6(g). We show these fits in Fig. 14.

Supplemental Movies

In [79], we include supplemental movies for the homogeneous experiments found in the main text, including In these plots, c is transformed according to Eq. 3 to make the decay clear.

many more snapshots than could be put in the main text. We show multiple realizations of each experiment, with the realization in the main text included as the left-most column. We pair the snapshots with the corresponding evolution of 1f for the given realization, calculated as in Fig.

FIG. 1 .

 1 FIG. 1. Crystal melting (a-d) and freezing (e-h) transitions in a monodisperse system. (a, e): snapshots in the z ≡ 0 plane of (top row) the Θ order parameter and (bottom row) the mobility µ as a nucleus appears and grows in an individual trajectory. (b, c, f, g): distributions (b, f) P (log 10 Θ) and (c, g) P (log 10 µ) at times (b, c) t = 100, 200, . . . , 600 and (f, g) t = 1000, 2000, . . . , 7000 (later time, lighter color). The red curves in (b) and (c) are equilibrium distributions.(d, h): untransformed fraction 1f for f calculated from P (log 10 Θ) (dashed red), P (log 10 µ) (dashed green), potential energy V (solid orange), inherent state potential energy VIS (solid blue), Θ (solid red), log 10 Θ (dotted red), µ (solid green), and log 10 µ (dotted green).

FIG. 2 .

 2 FIG. 2. Dynamics of melting and freezing for temperaturesT low = 2.0 and T high = 2.1 tightly straddling the critical point, given an artificial target-state slab at time t = 0. (a, b): slice along z = 0 of the initial configuration, colored according to Θ on a logarithmic scale from 10 -3 to 10 -1 (higher Θ, lighter color). (c, d): average Θ (x) at times t = 0, 50, . . . , 300 for melting, and t = 0, 40, . . . , 400 for freezing (later time, lighter color). Transparent curves correspond to raw data, solid curves to fits. Here, x = 2 |x| /Lx, where Lx is the time-dependent system size along the x-axis. (e-h): change, relative to their initial value, of the (e, f) interface center x * and (g, h) relative width λ/λ0 obtained when fitting the profiles to Eq. 4. Here, the initial interface width λ0 has values λ0 = 0.0213 and 0.0159 in (e) and (f) respectively. The dashed lines in (e) and (f) are linear fits.

TFIG. 3 .

 3 FIG.3. Heating (a-d) and cooling (e-h) of a two-dimensional polydisperse supercooled liquid between temperatures 0.06 and 0.14. (a, d): snapshots of the mobility µ as the system evolves, exhibiting domain growth during heating and quadrupoles during cooling, once rearrangements are sufficiently rare. (b, c, f, g): evolution of the (b, f) Θ and (c, g) log 10 µ distributions at times (b, c) t = 6000, 12000, . . . , 60000 and (f, g) t = 10 2.0 , 10 2.5 , . . . , 10 6.0 (later time, lighter color). The red curves in (f) and (g) are equilibrium distributions, obtained from equilibrium trajectories at the target temperature and pressure. (d, h): untransformed fraction 1f for f calculated from P (log 10 Θ) (dashed red), P (log 10 µ) (dashed green), V (solid orange), VIS (solid blue), Θ (solid red), log 10 Θ (dotted red), µ (solid green), and log 10 µ (dotted green).

FIG. 4 .

 4 FIG.4. Dynamics of (left column) heating and (right column) cooling supercooled liquids between temperatures T low = 0.06 and T high = 0.14 straddling TMCT, given an artificial targetstate slab at time t = 0. (a, b): initial configuration for (a) heating and (b) cooling, colored according to Θ on a linear scale from 0.0 to 0.1 (higher Θ, lighter color). (c, d): average Θ (x) at times (c) t = 0, 2000, . . . , 20000 and (d) t = 10 2.0 , 10 2.5 , . . . , 10 6.0 (larger time, lighter color). Transparent curves correspond to raw data, solid curves to fits. As before, x = 2 |x| /Lx, where Lx is the time-dependent size of the system along the x-axis. (e-h): change, relative to their initial value, of the (e, f) interface center x * and (g, h) relative width λ/λ0 obtained when fitting the profiles to Eq. 4.Here, the initial interface width λ0 has values λ0 = 0.00428 and 0.00414 in (e) and (f) respectively. The dashed line in (e) is a linear fit, while the dashed line in (g) represents diffusive growth ∼ t

6 FIG. 5 .

 65 FIG. 5. Structural evolution of the triangular plaquette model (TPM) and square pyramidal plaquette model (SPPM) during heating and cooling. Systems are initialized at equilibrium at temperature T0 = 0.0 for heating (a, b) and 1.0 for cooling (c, d), and the system is then evolved with Glauber dynamics at temperature T∞ = 0.3 (a), 0.5 (b), or 0.2 (c, d). Markers correspond to the locations of excited plaquettes in the given snapshot of the system. In (b), we show a slice along the y-plane containing the largest number of excited plaquettes at time t = 15000, so as to capture the birth of the domain. In (d), we instead show a slice along the y-plane containing the smallest number of excited plaquettes at time t = 10 6 , such that if a low-excitation domain were present at this time, we would capture it.

FIG. 6 .

 6 FIG.[START_REF] Vila-Costa | Nucleation and Growth of the Supercooled Liquid Phase Control Glass Transition in Bulk Ultrastable Glasses[END_REF]. Dynamics (left columns) heating and (right columns) cooling the TPM and SPPM from temperature T0 to T∞ for (T0, T∞) = (0.0, 0.3) (first column), (0.0, 0.5) (second column), (0.3, 0.2) (third column), and (0.5, 0.4) (last column), given an artificial target-state slab at time t = 0. (a-d): initial configuration for (a) TPM heating, (b) SPPM heating, (c) TPM cooling, and (d) SPPM cooling, with markers corresponding to excited plaquettes (in the y = 0 plane for the SPPM). The hotter region is shaded red, while the colder region is shaded blue. (e-h): normalized concentration c = 1 2 + c-(c∞+c 0 )/2 |c∞-c 0 |

  FIG.7. Dynamics of the phenomenological model governed by Eq. 7 for constants ls = 0.18 and lD = 0.25 given an artificial target-state slab at time t = 0. Heating (left column) is performed between temperatures T0 = 0.0 and 0.3 using τ0 = 0.549 and a = 0.253 and cooling (right column) is performed between temperatures T0 = 0.3 and 0.2 using τ0 = 0.705 and a = 1.00. Parameters are chosen so as to yield behavior similar to the TPM, establishing the reasonableness of this model, but no explicit fit to the TPM is performed. (a, b): normalized excitation state c defined as in Fig.6(e) at times (a) t = 0, 1000, . . . , 10000 and (b) t = 0, 10 3 , 10 4 , 10 5 , and 10 6 (larger time, lighter color). Transparent curves correspond to numerically-solved Eq. 7, solid curves to fits to the x >= 0 part of this solution, noting that c (-x) = c (x). (c-f): change, relative to its initial value, of the interface center position x * (c, d), and width λ (e, f) obtained when fitting the profiles to Eq. 4. The dashed lines in (c) and (d) are linear fits, while those dashed line in (e) and (f) represent diffusive growth ∼ t

FIG. 8 .

 8 FIG. 8. (a-d): Evolution of 1-f in homogeneous experiments calculated from Eq. 3 using the mean value of Θ for ten individual realizations contributing to the averaged results shown in (a) Fig. 1(d), (b) Fig. 1(h), (c) Fig. 3(d), and (d) Fig. 3(h).(e-h): Evolution of 1f calculated from Eq. 3 using the concentration c of excited plaquettes for ten individual realizations in the same regimes as those shown in Fig.5.

10 1 10 2 10 3 10 4 10 5 10

 10

FIG. 11 .

 11 FIG. 11. Validation of Eq. 2. Each row corresponds to data from the main text: (a, b): Fig. 1(b), (c, d): Fig. 1(f), (e, f): Fig. 3(b), (g, h): Fig. 3(f). Left hand side: same data as in the corresponding figure from the main text. Right hand side: output of Eq. 2 (where ϕ = log 10 Θ (a-d) or Θ (e-h)) with f calculated from Eq. 3 applied to the quantity Θ.

  to these fluctuations. Here, I is an appropriately-chosen interval. This approach is used to obtain the red and green dashed curves in Fig.1(d) (I = (-∞, -1.75) (red) and (-∞, -2) (green)), as well as the green dashed curve in Fig.3(d) (I = (-0.5, ∞)).

FIG. 12 .

 12 FIG.12. Melting transition in a monodisperse system heated from T0 = 1.90 to T∞ = 2.40 (c.f. T∞ = 2.35 in Figs.1(a)-1(d)). Distributions (a) P (log 10 Θ) and (b) P (log 10 µ) at times t = 0.0, 3.2, . . . , 48.0 (later time, lighter curve). (c): untransformed fraction 1f for f calculated from (dashed red) P (log 10 Θ) and (dashed green) P (log 10 µ), plotted against the mean (solid orange) potential energy V , (solid red) Θ, (dotted red) log 10 Θ, (solid green) µ, and (dotted green) log 10 µ, affinely transformed according to Eq. 3 with t0 = 2.0. On these short time scales, start-up of the barostat leads to erratic behavior of the inherent state energy, which is therefore not shown.

1 FIG. 13 .FIG. 14 .

 11314 FIG. 13. Perturbation to the local Θ field due to rearrangements in the homogeneously cooled glass studied in Figs. 3(e)-3(h) at time t = 10 6 . (a): Mean absolute change in Θ between times 10 6 and 10 6 + 10 2 for particles a distance r from a rearranging particle, averaged over ten independent realizations. The inset shows this data on a log-scale after subtracting the global mean absolute change in Θ. The dashed line corresponds to a fit to an exponential decay with decay length ξ perturb = 2.43. (b): Change in Θ across the same time interval as (a) for the realization shown in Fig. 3(e).

8 .

 8 Supplemental Movie 1 corresponds to Figs. 1(a) and 8(a), Supplemental Movie 2 corresponds to Figs. 1(e) and 8(b), Supplemental Movie 3 corresponds to Figs. 3(a) and 8(c), and Supplemental Movie 4 corresponds to Figs. 3(e) and 8(d). Supplemental Movies 5 and 6 are copies of Supplemental Movies 3 and 4, respectively, but with the Θ plot replaced with the change ∆Θ in Θ between times t and t + 100.

ACKNOWLEDGEMENTS

The authors thank C. Scalliet for sharing data sets from [46] with them. R.N.C. is grateful to L. Berthier, H.-H. Boltz, R.C. Dennis, Z. Fakhraai, M.A. Galvani Cunha, C. Herrero, P. Luo, P.D. Olmsted, M. Ozawa, and A.G. Yodh for helpful discussions concerning this project. The design of Figs. 1 and 3 in this paper was heavily inspired by the figures in [75]. This work was supported by the Simons Foundation via the "Cracking the glass problem" collaboration (#454945, RNC and AJL). AJL also thanks the Simons Foundation for support via #327939 as well as for the hospitality of the Center for Computational Biology at the Flatiron Institute, as well as the Isaac Newton Institute for Mathematical Sciences at Cambridge University (EPSRC grant EP/R014601/1), for support and hospitality.

Appendix A: Methodological Details

Swap potential

For convenience, we note here that in Eq. 1,

where rcut = 1.25. We also note that A, σ min , and σ max in the probability density function

for the particle size σ can be found by solving

numerically for A and noting that

and

Appendix B: Supporting Details

Perturbed Crystal

In §III and §IV, we considered systems well-described by the simple interpolation between initial and target states expressed in Eq. 2. We stressed that this implies a clean two-state equilibration scenario, as the systems evolve directly from the initial state into the final state. A good test of this idea is provided by systems in which the initial state is perturbed during the heating process, violating Eq. 2. In Fig. 3, for instance, there is a visible rightwards shift of the P 0 (log 10 µ) peak in Fig. 3(c) as target-state domains grow and a corresponding devia-


Arthur Adjedj

Meven Lennon-Bertrand

Kenji Maillard

Loïc Pujet

Engineering logical relations for MLTT in Coq

.

MLTT and its metatheory Establishing meta-theoretic properties such as the existence of canonical forms or decidability of the derivability of judgements of dependent type systems is a notoriously complex endeavour. For instance, the MetaCoq project [START_REF] Sozeau | The MetaCoq Project[END_REF][START_REF] Sozeau | Coq Coq correct! verification of type checking and erasure for Coq[END_REF] aims at entirely formalizing the Predicative Calculus of Cumulative Inductive Constructions (PCUIC) underlying Coq, and showing the correctness of an implementation of a typechecker. However, to do so it assumes an axiom stating that the theory is normalizing.

Indeed, for these meta-theoretical properties, type dependency precludes most proof strategies, which ultimately rely on a stratification between types and terms. To go beyond these limitations, Abel et al. [START_REF] Abel | Decidability of conversion for type theory in type theory[END_REF] formalize an inductive-recursive [START_REF] Dybjer | Induction-recursion and initial algebras[END_REF] definition of a logical relation for a representative fragment of MLTT in Agda, to show normalization and decidability of conversion for this theory. This technique was further extended to more complex theories [START_REF] Gilbert | Definitional Proof-Irrelevance without K[END_REF][START_REF] Pujet | Impredicative Observational Equality[END_REF].

The use of the induction-recursion scheme however introduces a new gap between the object theory being formalized (which only supports a handful of simple inductive types), and the meta-theory used to formalize the result. While exploring normalization proofs for complex inductive schema is a very valuable endeavour, we wish to go the other way around, and narrow this gap by using only regular indexed inductive types. This is both a requirement and a benefit of working in Coq, which only handles this class of inductive types. Thus, we reformulate the logical relation using small induction-recursion, which can in turn be encoded using simple indexed inductive types [START_REF] Hancock | Small Induction Recursion[END_REF]. This strategy requires definitions that are replicated across several universe levels, for which the universe polymorphic features of Coq come in handy.

A bidirectional presentation of MLTT In [START_REF] Abel | Decidability of conversion for type theory in type theory[END_REF], only decidability of conversion is shown. While this is definitely the most intricate part of showing decidability of type-checking, going from the former to the latter is non-trivial. Indeed, type-checking for MLTT as defined in [START_REF] Abel | Decidability of conversion for type theory in type theory[END_REF] is not, in general, decidable, for lack of annotations [START_REF] Dowek | The undecidability of typability in the Lambda-Pi-calculus[END_REF][START_REF] Heine | Lectures on the Curry-Howard isomorphism[END_REF].

In our development, we show decidability of typing, by extending algorithmic conversionchecking to a full account of algorithmic typing, described in a bidirectional fashion [START_REF] Lennon-Bertrand | Complete Bidirectional Typing for the Calculus of Inductive Constructions[END_REF][START_REF] Lennon-Bertrand | Bidirectional Typing for the Calculus of Inductive Constructions[END_REF]. Following the strategy implemented for instance by Coq, we use annotated (Church-style) abstractions, so that inference is complete, i.e., every well-typed term infers a type.

More generally, ideas from bidirectional typing help greatly in guiding the definition and handling of the algorithmic parts of the system. For instance, for most proofs on the algorithmic system we rely on a custom induction principle which threads the invariants maintained by a bidirectional algorithm, giving us extra hypotheses for each induction step. This lets us handle once and for all these invariants that are required for most proofs, rather than bundling them in the predicate proven by induction, which would mean showing their preservation again during each proof by induction.

Three logical relations in one We show that a bidirectional presentation of our type theory is equivalent to its standard declarative presentation given in the Martin-Löf Logical Framework. To do so, we parametrize the definition of the logical relation by a generic typing interface. The interface is instantiated 3 times: once with the declarative typing, once with a mixed system with declarative typing and algorithmic conversion, and finally with a fully algorithmic system using bidirectional typing. These different instances are used to gradually show properties of the system: the declarative instance lets us show enough good properties of the declarative system to be able to show that the mixed system fits the generic typing interface; the mixed instance proves that declarative and algorithmic conversion coincide, which is used to show that bidirectional typing fits in the generic interface; finally the fully algorithmic instance gets us the desired equivalence, showing that the type-checking algorithm is sound and complete.

Engineering aspects We rely on Autosubst [START_REF] Stark | Autosubst 2: reasoning with multi-sorted de Bruijn terms and vector substitutions[END_REF][START_REF] Dapprich | Generating Infrastructural Code for Terms with Binders using MetaCoq and OCaml[END_REF] to deal with all the aspect of the raw syntax, defining untyped renamings and substitutions, generating boilerplate lemmas on these, as well as providing tactical support to discharge equational obligations.

To support working with multiple different notions of conversion and typing, sometimes simultaneously, we devised a generic notation system based on type classes in the Math Classes style [START_REF] Spitters | Type classes for mathematics in type theory[END_REF], with a system of tags that lets us disambiguate between the different notions when needed, but can be ignored safely when working on a single notion or parametrically.

To ease the development, we use tactics to provide for judgement-independent notions, e.g. we use a single irrelevance tactic to discharge goals requiring a lemma stating that some form of the logical relation is irrelevant in some of its parameters. An important part of the work achieved by the definition of the logical relation consist in its generalization of typing contexts through Kripke-style quantifications over renamings and substitutions, and we use instantiation tactics to automatically apply lemmas to the relevant hypotheses.

Future work

We could add more universes to obtain a hierarchy of arbitrary finite length, and see no theoretical obstacle in doing so. We plan to extend the formalization to a scheme of indexed inductive types used in Coq, hence narrowing the gap between the object theory and the metatheory used to prove its normalization: this would lead to a formalization of MLTT with n universes into MLTT with n + k universes for a (small but strictly positive) constant k.

On the bidirectional side, we do not cover the common pattern, used for instance in the kernel of Agda, of having some terms that only check (typically, unannotated abstraction). It would be interesting to obtain a decidability result that covers these as well.

Finally, there is a large space left to improve automation, taking inspiration from the rich Coq ecosystem. Indeed, the main difficulty for a proof by logical relations is in the setup of the relation, but most proof obligations are rather repetitive and unsurprising. While our tactics already relieve us from quite a bit of this tedious work, they are far from making it all disappear.