Phase transition for quenched coupled replicas in a plaquette spin model of glasses

We study a three-dimensional plaquette spin model whose low temperature dynamics is glassy, due to localised defects and effective kinetic constraints. While the thermodynamics of this system is smooth at all temperatures, we show that coupling it to a second system with a fixed (quenched) configuration can lead to a phase transition, at finite coupling. The order parameter is the overlap between the copies, and the transition is between phases of low and high overlap. We find critical points whose properties are consistent with random-field Ising universality. We analyse the interfacial free energy cost between the high- and low-overlap states that coexist at (and below) the critical point, and we use this cost as the basis for a finite-size scaling analysis. We discuss these results in the context of mean-field and dynamical facilitation theories of the glass transition.

Introduction.-There are several theories that aim to explain the rapid increase in the viscosity of supercooled liquids close to their glass transitions [1][2][3]. Some theories predict that phase transitions should occur, either on cooling or in response to some kind of external field [4][5][6][7][8]. Such phase transitions are relevant because their associated order parameter fluctuations may explain the dynamical fluctuations characteristic of supercooled liquids [9]. The theory of dynamical faciliation (DF) [5,10] is based on a dynamical order parameter, and the associated dynamical phase transitions occur in ensembles of trajectories [7,11,12]. This theory is encapsulated by a class of simple model systems that includes kinetically constrained models [13] and plaquette spin models [14]. A different theoretical approach is based on mean-field calculations within a replica formalism [6,15]. It proposes a static order parameter-the overlap-which measures the similarity between different configurations of the system. This approach predicts phase transitions [16], as found in certain spin-glass models, that occur when equilibrium configurations of the system are biased to be similar to fixed reference configurations.
Here, we consider a three-dimensional plaquette spin model-the square pyramid model (SPyM) [17]-whose relaxation behavior is well described by DF theory. Using computer simulations, we show that this model also exhibits phase transitions when biased by its overlap to (quenched) reference configurations. The properties of the thermodynamic phase transitions that we find are consistent with the universality class of the random-field Ising model (RFIM), as predicted by field theories for generic glassy systems [16,18]. Evidence for this scenario has also been found in atomistic simulations [19,20]. Thus, at low temperatures, the SPyM exhibits both the dynamic fluctuations characteristic of DF and the overlap fluctuations characteristic of mean-field theory. These results are important because (i) they strongly support the predictions of RFIM criticality [16,18] in glassy systems, and (ii) the fact that a single model displays features of both mean-field theories and DF theory hints at a connection between these apparently contradictory scenarios [3].
Model.-The SPyM [17] consists of Ising spins s i ¼ AE1 on the vertices of a body-centered cubic lattice, with periodic boundaries. We denote a configuration of the system by C ¼ ðs 1 ; …; s N Þ; its energy is E 0 ¼ −ðJ=2Þ P μ s i μ s j μ s k μ s l μ s m μ , where the sum runs over upward-pointing square pyramids within the lattice, and the five spins s i μ Á Á Á s m μ are on the five vertices of pyramid μ. The parameter J sets the energy scale, and the linear system size is L, with a total of N ¼ L 3 spins. Pyramids μ for which the product of spins is equal to −1 are defects that carry an energy J: from a thermodynamic viewpoint their statistics are those of an ideal gas. Hence, for a single system there are no thermodynamic phase transitions at any finite temperature. However, the dynamical behavior of the defects is complex and cooperative: by analogy with a similar model in two dimensions [14], it is believed that for low temperatures T, the relaxation time τ diverges as ln τ ∼ T −2 [17], and the range of certain multipoint spin correlations also diverges [21]. The overlap is QðC; C 0 Þ ¼ ð1=NÞ P i s i s i 0 , which measures the similarity between spin configurations. We draw a configuration C 0 at random from an equilibrium distribution at temperature T 0 . Holding C 0 fixed, we calculate expectation values with respect to the distribution pðCjC 0 Þ ¼ ZðC 0 Þ −1 e ½εNQðC;C 0 Þ−E 0 ðCÞ=T , where ZðC 0 Þ is a partition function. Thus, for ε ¼ 0 we consider equilibrium behavior of an isolated SPyM at temperature T, while increasing the coupling ε biases the configuration C to increase its overlap with the reference C 0 . Finally, we perform an average over the reference configuration C 0 . We fix the energy scale J ¼ 1, so the dimensionless parameters of the system are (T, T 0 , ε). We use a Monte Carlo scheme to study these coupled systems, as described in the Supplemental Material [22].
Results.- Figure 1 shows results for T 0 ¼ T, in which case the reference configuration is representative of thermal equilibrium. The mean overlap hQi in Fig. 1(a) increases with ε. The increase is gradual for small ε, but becomes steeper at larger ε. In spatial dimensions d ≥ 3, theory predicts [6,16,18] that the gradient dhQi=dε should diverge at a RFIM critical point at some (T c , ε c ), and that there is a first-order transition for T < T c , with a discontinuous jump in hQi. Numerical calculations for atomistic systems are consistent with these predictions of static (thermodynamic) phase transitions [19,20]. Similar critical points are predicted [8,25] in systems with pinned particles, and are also consistent with numerical simulations [26,27]. However, in numerical calculations on finite systems, divergent features in hQi are not observed: here, we use finite-size scaling methods to provide strong evidence for the existence of a RFIM critical point in the SPyM.
The reference configuration C 0 is a source of quenched disorder in this problem, and averages are calculated in two stages: first a thermal average at fixed C 0 , denoted by h·i C 0, and then a disorder average, ð·Þ. The notation h·i ≡ h·i C 0 indicates the double average. Figure 1(b) shows the average hQi as well as the behavior of hQi C 0 for eight representative configurations C 0 . For each C 0 , one sees a very sharp jump in the overlap at a sample-dependent field ε Ã C 0 . However, the jump in the average overlap is broadened out due to fluctuations in ε Ã C 0 between samples. This indicates that the disorder is the dominant source of fluctuations in the problem, as expected at a RFIM critical point [28]. Figure 1(c) shows the (total) susceptibility χ tot ¼ As the temperature is reduced, the fluctuations grow rapidly and depend increasingly strongly on the system size, as expected in the vicinity of a phase transition. Finally, Fig. 1(d) shows the distribution of the overlap Q, evaluated at the field ε Ã that maximizes χ tot . For the lower temperatures, a two-peaked structure is clearly visible. Considering Figs. 1(b) and 1(d) together, one sees that typical reference configurations contribute to either the low-Q or high-Q peak, as expected for RFIM criticality: the proportion of configurations C 0 that contribute simultaneously to both peaks scales as L −θ [28], where θ is a critical exponent for the RFIM (for d ¼ 3 then θ ≈ 1.5 [29,30]). The numerical uncertainties in PðQÞ are considerable: where the probability is small our estimate of PðQÞ may be dominated by just one or two out of the 64 reference configurations considered. We return to this point below.
Taken together, the data in Fig. 1 are consistent with a RFIM critical point for T 0 ¼ T in the range 0.3-0.4. However, while the Monte Carlo methods used here enable efficient simulation close to phase transitions [22], the slow (glassy) relaxation of the model limits the ranges of accessible temperatures and system sizes. If a RFIM critical point exists for T 0 ¼ T, we expect similar critical points at other reference temperatures T 0 . We now consider lower T 0 , which provides further evidence for RFIM criticality.
For T 0 ¼ 0.25, Fig. 2(a) shows a sharp jump in hQi as ε increases, this time at a higher temperature T ≈ 0.8 for which relaxation is faster, allowing investigations of larger systems, up to N ¼ 32 3 . Figure 2(b) shows the distribution PðQÞ, whose form changes from unimodal to bimodal as T decreases.
To identify the critical point, we use finite-size scaling. However, the RFIM has unusual scaling properties: the transition is almost first order and the distribution PðQÞ has two well-separated peaks even at the critical point [29]. The order parameter exponent β is also very close to zero [30]. Together with the significant numerical uncertainties in PðQÞ, these two features lead to difficulties with classical finite-size scaling based on universal cumulant ratios or order parameter distributions. To address this problem, we follow Ref. [29]: for each reference configuration C 0 , we calculate the coupling ε Ã C 0 that maximizes the variance hδQ 2 i C 0 . At this coupling, the distribution of Q is typically bimodal, even above T c [29]. We define a free energy for this reference configuration F Ã ðQÞ ¼ −T ln P C 0 ðQÞ. This free energy has two minima, and the height of the barrier between them is the interfacial free energy cost between high-and low-overlap states. To obtain an average free energy cost, we calculate the average free energy F Ã ðQÞ from which we extract the barrier height ΔF . (For large systems, F Ã ðQÞ=N → VðQÞ − ε Ã Q, where VðQÞ is the Franz-Parisi potential [6].) Figure 2(c) shows the free energy (scaled by β ¼ 1=T to allow interpretation as a log probability). In contrast to PðQÞ, the numerical uncertainties in F Ã ðQÞ are straightforward to estimate, and relatively small. RFIM critical points are characterized by unusual scaling behavior for the interfacial free energy cost ΔF . For T < T c one finds a surface cost ΔF ∼ L d−1 , as in conventional critical behavior. For T > T c one finds ΔF ¼ 0 for sufficiently large N, but ΔF may be significant in finite systems. For T ¼ T c , one has ΔF ∼ L θ [28,29], in contrast to conventional critical points where ΔF tends to a constant value. In Fig. 2(d), we plot βΔF =L θ for two system sizes (we take θ ¼ 1.5, which is known from detailed numerical analyses of the RFIM [29,30]). The interfacial costs ΔF are large and they grow with system size [ Fig. 2(b)]; they grow increasingly rapidly at low temperatures [ Fig. 2(c)]. This represents strong evidence for the existence of a critical point, and is consistent with RFIM universality. From the crossing point in Fig. 2(d), we estimate T c ≈ 0.83. Since we have data for only two system sizes, this estimate is subject to some uncertainty, but we argue that Fig. 2 presents strong evidence for a RFIM critical point in this system.
Interfacial costs for T 0 ¼ T.-Armed with these results for T 0 ¼ 0.25, we now return to the important case T 0 ¼ T. Figure 3(a) shows the free energy F Ã ðQÞ, at T ¼ T 0 ¼ 0.333. The interfacial cost increases strongly with system size; its temperature dependence is shown in Fig. 3(b). For T 0 ¼ 0.357 and T ≥ T 0 , the scaled free energy cost always decreases as the system size is increased from N ¼ 8 3 to N ¼ 16 3 , indicating that these systems are all above T c . The free energy costs are also nonmonotonic in temperature-our interpretation of this unusual feature is that the high-overlap and low-overlap states have increasingly similar structures as T approaches T 0 , which tends to reduce the surface tension between them. Moreover, for T 0 ¼ 0.333, the scaled free energy costs for the two system sizes are very close to each other over a wide range of temperature T. [This is in contrast to the clear crossing of the curves in Fig. 2(d).] Our interpretation of this last result is that all these systems are close to criticality.
Summary of phase behavior.- Figure 4(a) shows the estimated positions of the critical points that we have found in the SPyM. For T 0 ¼ 0 the reference configuration C 0 is in its ground state (s i ¼ 1 for all i), and the system reduces to a SPyM in a magnetic field, for which there is an Ising critical point at T c ¼ 0.98 [17,31]. For T 0 > 0, we expect a line of critical points, all of RFIM type, as in Fig. 4(a). They separate a region of parameter space in which phase coexistence is possible ("two phase") from a one-phase region where the response to the coupling ε is smooth. The Franz-Parisi potential VðQÞ is strictly convex in the onephase region but includes a linear segment (Maxwell construction) in the two-phase region. (Evidence for the critical point with T 0 ¼ 0.286 is shown in the Supplemental Material [22].) The dashed line in Fig. 4(a) indicates T ¼ T 0 . It intersects the line of critical points at T 0 ≈ 0.33. When coexistence takes place at T ¼ T 0 , the nature of the interfacial free energy cost has important implications for theories of the glass transition. Since C 0 is an equilibrium configuration of the isolated system, one expects the high-Q state to have configurations C that are also close to equilibrium. Similarly, the low-Q states are weakly perturbed from equilibrium, so their structure is also close to equilibrium. Thus, the phase coexistence illustrated in Fig. 3 is taking place between different states that have (statistically) very similar structures. The interfacial cost associated with such states plays a central role in random first-order transition (RFOT) theory [32], because C 0 may be viewed as a form of self-generated disorder. However, the nature of these interfaces is not fully understood [33,34]. As noted above, the reduction in interfacial cost as T approaches T 0 in Fig. 3(b) reflects the increasingly similar structures of the high-and low-overlap states. It also indicates that the cost is primarily entropic and not energetic (otherwise βΔF ¼ ΔE=T would likely increase on cooling).
The existence of RFIM phase transitions in the SPyM is consistent with the general picture proposed by RFOT theory, but we emphasize that a single SPyM (ε ¼ 0) has no finite-temperature phase transitions. In general glassy systems one expects ε Ã =T ≈ s c ΔQ, where s c is the configurational entropy density and ΔQ is the jump in the overlap at ε Ã [35]. For plaquette models (including the SPyM) we expect s c to be comparable to the total entropy s ¼ e −1=T ð1=T − 1Þ, and ΔQ ∼ 1. Hence, for low temperatures, ε Ã ∼ e −1=T . Figure 4(b) shows the resulting scenario for the SPyM, including three data points that indicate the locus of ε Ã as T is varied, holding T 0 ¼ T [recall Fig. 1(a)]. Given the limited data, the fitted line ε Ã ∼ e −1=T should be regarded only as a qualitative prediction for the phase boundary, but we emphasise that the first-order line can meet the ε ¼ 0 axis only at T ¼ 0. (Similar scaling for ε Ã is found for annealed coupling between plaquette models [17,36].) Within mean-field theory, T c is the highest temperature for which metastable states (and configurational entropy) are well defined, and is sometimes identified as an onset temperature T o [8,35]. However, this temperature has no clear signature for a single SPyM (ε ¼ 0). Other definitions of T o are also possible and may yield different (usually larger) numerical values [3,37,38].
Discussion.-The behavior of the SPyM is characteristic of the DF perspective of the glass transition (absence of a thermodynamic transition, effective kinetic constraints, facilitated relaxation, and heterogeneous dynamics). We have shown, nevertheless, that when coupled to a quenched reference state, the SPyM displays thermodynamic overlap transitions in the RFIM class, as expected from mean-field calculations and RFOT theory. It follows that quenched overlap transitions in atomistic simulations [20] are not incompatible with DF theory; also transitions occurring at T, ε > 0 do not imply any finite-temperature singularities at ε ¼ 0 [ Fig. 4(b)]. More generally (see also Refs. [39][40][41]) the presence of overlap transitions in the SPyM shows how these models can be consistent at the same time with predictions of both DF and RFOT theory, offering a link between these scenarios.
A crucial aspect of such a link is whether (and how) phase transitions for ε, T > 0 are related to the unbiased (ε ¼ 0) properties of the system. Near the first-order transition line in Fig. 4(b), we expect a nucleationlike mechanism for dynamical relaxation, with an initially highoverlap state decaying to a low-overlap one via nucleation and growth of a droplet of the low-overlap phase. This situation is close to the RFOT picture described in Ref. [32]: for ε < ε Ã a nucleation argument yields a free energy barrier βΔF ∼ ðε Ã − εÞR d c ∼ ðε Ã − εÞ 1−d , where R c ∼ 1=ðε Ã − εÞ is the size of the critical nucleus, which diverges at the first-order transition line. However, for ε ¼ 0 the natural picture is that of DF theory: pointlike defects facilitate cooperative rearrangements over a length scale ξ ∼ e 1=ðTd f Þ , where d f , is a scaling exponent [14]. The free energy barrier for these processes scales as βΔF ∼ ð1=TÞ ln ξ ∼ ð1=TÞ 2 . Bridging between the two regimes ε ∼ ε Ã and ε ≃ 0 remains a challenge, but it does represent a possible route for connecting the defectmediated dynamics of DF theory to the nucleation picture predicted by RFOT and mean-field theories.
The supporting data for this article will be openly available from the University of Bath data archive [42], shortly after publication.