Fragile-glass behavior of a short range $p$-spin model

In this paper we propose a short range generalization of the $p$-spin interaction spin-glass model. The model is well suited to test the idea that an entropy collapse is at the bottom-line of the dynamical singularity encountered in structural glasses. The model is studied in three dimensions through Monte Carlo simulations, which put in evidence fragile glass behavior with stretched exponential relaxation and super-Arrhenius behavior of the relaxation time. Our data are in favor of a Vogel-Fulcher behavior of the relaxation time, related to an entropy collapse at the Kauzmann temperature. We however encounter difficulties analogous to those found in experimental systems when extrapolating thermodynamical data at low temperatures. We study the spin glass susceptibility investigating the behavior of the correlation length in the system. We find that the the increase of the relaxation time is not accompanied by any growth of the correlation length. We discuss the scaling properties of off-equilibrium dynamics in the glassy regime, finding qualitative agreement with the mean-field theory.


I. INTRODUCTION
The glassy state is very common in nature 1 . When it is reached from the liquid phase lowering the temperature, one finds a dramatic increase of the relaxation time, and off-equilibrium phenomena can not be avoided on experimental time scales. This leads to non analytic behavior of the thermodynamic quantities, with a "transition temperature" that depends on the cooling rate. Despite its ubiquity, the basic mechanisms underlying the common features, as well as the peculiarities of the glassy behavior in different systems are yet to be clarified.
One of the most suggestive ideas in the glass theory, proposed long time ago by Gibbs and Di Marzio 2 , relates the increase of the relaxation time, and the observed finite time singularities to the existence of a thermodynamic transition at the Kauzmann temperature where the configurational entropy collapses to zero 3 . Soon after, in a refinement of the argument, Adams and Gibbs argued in favor of a Vogel-Fulcher singularity in the relaxation time.
Disordered systems have been proposed as paradigmatic models in which glassy phenomena can be studied in a nutshell, and theoretical ideas tested on microscopic models 4 . This is due to the fact that in disordered systems, the glassy state already appears in mean-field theory. The natural separation of the variables among "quenched" and "annealed" allow for the successful use of powerful techniques as the replica method for static 5 and functional methods in dynamics 6 . In fact, a satisfactory mean-field theory of disordered systems for the static 5 as well as for the equilibrium 6 and off-equilibrium dynamics exists 7,8 . Recently, examples of mean-field deterministic models with glassy behavior very similar to the one of disordered systems have been displayed 9 . This points in the direction that common mechanisms could lead to the glassy behavior of disordered and non-disordered systems.
The simplest example in which the Gibbs-Di Marzio collapse occurs is the Random Energy Model of Derrida 10 , and is a common feature to all systems with a "discontinuous glassy transition", or technically "onestep replica symmetry breaking", where the Edwards-Anderson parameter 11 undergoes a discontinuity 4 . Examples of such models are the Potts Glass 12 , the pspin interaction model 13 and a model of manifolds in a disordered media with short-range correlated disorder 14 . This class of systems has been proposed by Kirkpatrick, Thirumalai and Wolynes as simple toy models for the structural glass transition 4 . Notably, the study of the Langevin dynamics of the spherical version of these models shows that there the Mode Coupling Theory 15 is exact, and displays a dynamical singularity of kind B in the Götze classification 15 . In fact, recent progresses in the comprehension of the dynamics of mean-field disordered systems 7,8 allowed for an extension of the Mode Coupling Theory to the broken ergodicity phase 16,17 .
Many studies 4,18,19,7,20,21 have pointed out the exis-tence of a temperature T D where, despite the fact that no singularity is observed in the free energy, there is a statical breaking of ergodicity into an exponentially large number of metastable states. A thermodynamic singularity is present at a temperature T C smaller then T D . As a genuine mean-field theory, the Mode Coupling Theory neglects the activated (droplet) processes that in finite-dimensional systems are responsible for the decay of metastable states in a finite time. Kirkpatrick and Wolynes have recently stressed how the inclusion of these processes can restore ergodicity for T C < T < T D , and give rise to a generalized Vogel-Fulcher singularity at the static transition temperature T C 4 . The argument has been confirmed and refined by Parisi with a theoretical calculation based on the potential theory in spin glasses 22 .
The aim of this paper is to test this idea in a finitedimensional disordered model where metastable states should be present, but destabilized by activated processes. In section 2 we propose a finite dimensional analogous of the p-spin interaction model in the case p = 4. We believe that, as usual, the mean-field limit is recovered for high dimensionality. In section 3 we study the thermodynamics in the high-temperature regime through Monte Carlo simulations, that demonstrate that the model behaves as a fragile glass.
A second aspect of our work concerns the offequilibrium dynamics deep in the glassy phase. There, the properties of the system depend on the thermal history, and time translation invariance does not hold 23 . The off-equilibrium Mode Coupling Theory predicts scaling relations and a definite pattern of violation of the Kubo fluctuation-dissipation relation. In section 4 we study the dynamics in this regime, showing the consistency of the aforementioned scenario. Finally, the conclusions are drawn.

II. THE MODEL
The p−spin model 10,13 is defined by the long-range Hamiltonian where the couplings J i1,i2,...,ip are independent Gaussian variables with zero mean, and variance J 2 i1,i2,...,ip = p!/(2N p−1 ). The spins σ i , i = 1, ..., N can be taken as Ising variables, or as real variables subjected to the spherical constraint N i=1 σ 2 i = N 18 . The case p = 2, that corresponds to the Sherrington Kirkpatrick model, has a glassy transition with continuous order parameter 5 and will not be considered in this paper. For p ≥ 3 both in the Ising and in the spherical case the transition is discontinuous, and the properties of the model are the ones of interest in this paper. As a finite dimensional analogous of the model (1) in the case p = 4 we take a spin system with interacting Ising spins σ i arranged on the sites of a D dimensional square lattice with periodic boundary conditions. The Hamiltonian is defined as where the sums runs over all the plaquettes ✷ of the lattice. Each spin belongs to 4 D 2 plaquettes.
where the variables J ✷ are chosen as independent Gaussian variables with zero mean and unit variance. Note that, in generic dimension, the model is not invariant under Z 2 gauge transformation as it would be if the spins were located on the links (instead that on the vertices) of the plaquettes. Thanks to this, there is no Elitzur theorem preventing non-zero global order parameters 24 . In fact there are no local symmetry operations leaving H invariant 25 . However the Hamiltonian is invariant under the contemporary inversion of all the spins that belong to any hyperplane of dimension D − 1 orthogonal to one of the Cartesian axes. It is easy to check that these operations do not change the sign of any of the plaquettes. The degeneracy due to this symmetry (plane inversion symmetry in the following) is 2 DL−(D−1) , and can be removed e.g. fixing the spins on the Cartesian axes to arbitrary values. This exponentially large degeneracy of the states is also present in the ferromagnetic version of the model (J ✷ = 1). We think that this could lead to a very interesting spinodal dynamics, above the lower critical dimension D = 2. We concentrate here to the disordered model, leaving the study of the ferromagnetic case for future work. A ferromagnetic model with four spin interactions at the vertices of plaquettes was studied in 26 in connection with random surfaces physics. In that case also pair interactions that removed the plane inversion symmatry where present in the Hamiltonian.
In the limit of infinite dimension, where the number of plaquettes to which a spin belongs tends to infinity, one can expect that, modulo the symmetry, the model is equivalent to the p−spin model for p = 4. Models with different p's could be easily constructed for other lattices, e.g. the case p = 3 would correspond to a triangular lattice 27 .
It is worth at this point to present a brief qualitative review of the results of the mean-field theory, based on the Hamiltonian (1) 4,19,20 . The study of the thermodynamics of this system leads to the following results. At high temperatures the system is paramagnetic and ergodic. At a temperature T D the ergodicity breaks down, and an exponentially large number (exp(N Σ(T ))) of pure states (ergodic components) separated by barriers of order N contribute to the partition function 4 . This transition occurs without singularities in the free energy, which is equal to the free-energy per state F in , plus an entropic contribution −T Σ(T ) coming from the multiplicity. The quantity Σ, the configurational entropy, is a decreasing function of the temperature, and at a temperature T C < T D vanishes. At T C there is a thermodynamic phase transition with singularities in the free-energy.
On the other hand, the study of the off-equilibrium dynamics after a sudden quench from a high temperature shows that the large-time limit of various dynamical quantities is non analytic and the singularity in the dynamics at T D is suppressed. The characteristic time scale τ for these processes to restore ergodicity has been recently estimated in Potts-glass models using heuristic arguments by Kirkpatrick and Wolynes 4 , and substantiated using a droplet argument in replica space by Parisi 22 . They find a generalized Adam-Gibbs relation of the kind where C is a constant and γ = D − 1. As the configurational entropy vanishes linearly near T C , eq.(3) results in a generalized Vogel-Fulcher law τ ∼ exp C (T −TC ) γ . The value γ = 2 in D = 3 is at variance with the usual value γ = 1 used to fit the experimantal data. However, in the case of the present model the value of γ = D − 1 in (3) should be lowered due to the plane inversion symmetry. We do not know if this would result in γ = D − 2, and we leave the investigation of this point for future work. As the matter of fact, simple (and trivial) results are obtained for the statics in D = 2. In that case, one can show that in the high temperature expansion only diagrams involving a number of spins proportional to L or higher, and hence irrelevant for L → ∞, are present. Accordingly the free-energy per spin is found F = −T dJ √ 2π exp(−J 2 /2) log(2 cosh(βJ)). We will see that the relaxation time follows a simple Arrhenius law in this case. The lowest dimension at which one can expect non trivial thermodynamical results is D = 3, where one can see that frustration is present. The study of the properties of the three-dimensional model through Monte Carlo simulations and the comparison with the results of the theory, will be the subjects of the rest of this paper. Some results for the two dimensional case will also be mentioned.

III. THERMODYNAMICS
In order to investigate the questions posed in the previous sections, we have performed Monte Carlo simulations of the model in three dimensions, using a standard serial single spin-flip heat-bath algorithm. The signature of glassy behaviour is easily seen in simulations of cooling experiments. In figure 1 we plot the energy as a function of the temperature for different cooling rates. We clearly see a change of behaviour corresponding to a jump of the specific heat around T ≃ 0.7, where the systems fails to reach equilibrium within the observation time. The inset shows that the "transition temperature" as well as the value of the energy at which the system freezes are dependent on the cooling rate. Since we expect the equilibrium entropy to be the relevant quantity for the transition, we integrated the high-temperature energy data in β to get the entropy, taking into account data from β = 0.01 up to β = 1.
In figure 2 we present the results of this operation, togheter with some rational function fits of the data points. The functional form that we have chosen to fit that data are: Even functions of β have been chosen coherently with the fact that the high temperature expansion of S only contains even powers of β. The fits are roughly of the same quality at high temperature and are both consistent with the entropy collapse scenario. Extrapolating at low temperature, one can estimate the point where the entropy vanishes to be in the range 0.3 ≤ T c ≤ 0.6. This range of temperatures should not be considered as more than a mere indication of what it could happen. First of all, the extrapolation at low temperatures by means of the functions (4) is highly arbitrary. Second, we expect only the configurational entropy, and not the total entropy (which take an intra-state contribution), to become negative at the transition. The impossibility to disentangle the two contributions adds incertitude to the critical temperature estimate. However, both the mean-field theory and the off-equilibrium dynamics of section 4 suggest that in all the low-temperature region the Edwards-Anderson parameter is close to 1, correspondingly the intra-state entropy is rather small, and possibly much smaller than the configurational one not too close to T c .
Coherent informations are obtained from the study of the equilibrium dynamics at high temperatures. The analysis of the autocorrelation function C(t) =< σ(t)σ(0) > shows how the relaxation follows a stretched exponential law of the type, As shown in figure 3 our data are pretty well fitted by the previous relaxation law all along the relaxation.
From that fit we are able to extract the temperature dependent relaxation time τ (T ) and the exponent b(T ). The behavior of these parameters as a function of the temperature is depicted in figure 4. The relaxation time is consistent with a Vogel-Fulcher law of the kind τ = A exp(B/(T − T 0 )) with T 0 ≃ 0.58, A = 0.5, B = 10 (continuous line in figure 4). But, analogously with what happens with experimental data, it is also well fitted by a zero temperature singularity of the type τ = A exp(B/T 2 ) with A ≃ 11.5, B ≃ 0.69. A four parameter fit of the form τ = A exp(B/(T − T 0 ) γ ) fits the data with γ = 1.01, while posing by hand γ = 2 as in (3) leads to T 0 = 0.001. The indication of a thermodynamic transition by the divergence of the relaxation time at T 0 is supported by the fact that it is in the range of temperatures where the extrapolated entropy vanishes. The exponent b is also depicted in the inset of figure 4. It appears to be linear with the temperature and of order .5 in the low temperature phase. Unfortunately we have no evidence stronger than this in favour of the finite temperature singularity. The best we can say is that our numerical data are consistent with the glass transition scenario as much as laboratory experiments on glasses give support to this singularity (with the difference that laboratory experiments can explore a larger window of time than in the numerical experiments). However the picture we get is coherent with the theoretical relation between the relaxation time and the configurational entropy. We have related the two quantities as suggested by eq.(3). In figure 5 we plot the logarithm of the relaxation time versus the inverse entropy, and we see that the data fall quite well on a straight line indicating the validity of (3) with γ = 1. This linear relation is again consitent with a small value of the intra-state entropy, that could be estimated of the order of S ∞ = 0.04, value obtained from the extrapolation of the fit in fig. 5 to τ → ∞.
It is worth at this point to mention the results of an analogous analysis for the two dimensional system, where the thermodynamics is trivial. Simulations performed in that case indicate that also in D = 2 there exists a low-temperature regime where the equilibrium correlation function behaves as a stretched exponential (we observe that stretched exponential behavior is also observed in the one-dimensional Ising model at low temperatures in the intermediate time regime 28 ). However, the relaxation time grows much slower at low temperature than in the three dimensional case. In fact, very good fits are obtained by simple Arrhenius forms τ ∼ exp(B/T ), i.e. the minimal increase to be expected in any system with discrete spin variables.
We note that the two-step relaxation characteristic of many structural glasses is not seen here.
To clarify further the picture about the thermodynamics of the three dimensional model, we have investigated the possibility of a growing static correlation length. This is done in disordered spin system through the analysis of the correlation function 29 among two different replicas σ i and τ i with same J ✷ 's, and evolving with independent thermal noise. To avoid problems with the symmetry of the system we have fixed all the spins σ i and τ i belonging to the Cartesian axes.
In the thermodynamic limit this should be an effective way to eliminate the degeneracy.
In our simulations we have not found any evidence of a growing correlation length, which appears to be independent of temperature and smaller than 1. We conclude that the divergence of the relaxation time is not accompanied by a divergent correlation length.
An alternative way to explore the existence of a divergent correlation length is to compute directly the integral of the correlation function G(x) over the whole space. This quantity yields the spin-glass susceptibility 30 defined by, where N = L 3 is the size of the system. For vanishing correlation length this quantity is equal to 1 and tends to increase whenever there is spatial ordering and the correlation length increases 30 . We have simulated several system sizes L = 4, 5, 6, 7 from T = 3.0 down to T = 0. Data is shown in figure 6. As one can expect, there are serious thermalization problems especially at low temperatures. Nevertheless, it clearly emerges form these results that, in the region of temperatures where we are in equilibrium (for instance, above T = 1.2 where the equilibrium relaxation time is smaller than the thermalization time we made the system evolve before measuring the observables) there is no evidence of a divergence of the spinglass susceptibility with the size of the system 31 . This is at variance with the results obtained in case of shortrange Edwards-Anderson spin glasses, where the growing of the correlation length has been generally observed 29 .

IV. OFF-EQUILIBRIUM DYNAMICS
In the previous section we have seen that the relaxation time increases very fast as the temperature is lowered, and it is equally well fitted by a finite-time Vogel-Fulcher law and by the law τ = A exp(B/T 2 ). Even in this second hypothesis, in the whole low temperature range T < 1, the time scales we can computationally reach for reasonably large systems are by far shorter than the equilibration time. History dependent effects and aging are then to be expected 23 .
Recent results in the study of the dynamics following a quench in the low temperature phase in spin glasses, give prediction about the scaling of the correlation function and the linear response function for large times in the off-equilibrium regime, and on the dependence of these quantities on the time t w that the system has spent at low temperature 7 . These functions are defined respectively as where h i is a local magnetic field applied to the system. For a complete exposition of the off-equilibrium theory of the glassy dynamics we refer the reader to 7,8 . Here we limit to reassume briefly some features relevant to the present discussion. For large t w the following scenario is found: there is a first regime, for small t (t << τ (t w ) see below) where the dynamics has features similar to those of an equilibrium system. The correlation function is independent of t w in this regime, and the response function is related to the correlation function by the fluctuation-dissipation theorem. In this regime the correlation function monotonically decreases from 1 to a value q EA , which defines an off-equilibrium parameter analogous of the Edwards-Anderson parameter. In addition to this equilibrium-like regime, in the class of models of interest for this paper, there is a regime in which the correlation decays from q EA to zero, and the correlation has the scaling form The "effective relaxation time" τ (t w ) is an increasing (and diverging) function of t w , that the theory is not able to predict, and seems to be rather system dependent. In some cases it is found τ (t w ) = t w 32-34 , we will see that this is not the case for the present model.
As far as the behaviour of the response function is concerned, it is found that the function sometimes called fluctuation-dissipation ratio, depends on its time arguments in a quite special way: through the dependence on t and t w of the correlation function itself, x(t, t w ) = x(C(t, t w )).
The value x = 1 corresponds to the Fluctuation-Dissipation Theorem relation and it is valid for q EA < C(t, t w ) < 1. A non zero value of x(C) is found when aging effects are present in the response function, and testifies the memory of the system about its history. In the mean-field p-spin model it is found that x(C) is equal to a constant x smaller than 1 in the whole interval 0 < C < q EA . For a discussion of the behaviour of this quantity in the Edwards-Anderson model, as well as for a qualitative discussion of its behaviour in a different glassy scenario see 35 .
The response function is measured from simulations of "zero-field cooled experiments" 36 . Starting at time zero from a random configuration, we let evolve the system at constant temperature in zero field for a time t w . At t w we switch on a small magnetic field h, and measure the relaxation of the magnetization as a function of time. In the linear response regime, the magnetization m(t, t w ) is given by a relation that, assuming the validity of (11), takes the form , t w )). (13) The mentioned behaviour of x(C) reflects in It has been noted in 37 that, while a scaling behavior of the kind (9) is common to a glassy behavior and phenomena of domain-growth in phase separation 38 , the function x(C) seems to be non zero only in the aging regime of the glassy systems. In order to discriminate glassy behavior from domain-growth like mechanisms another quantity has been recently proposed 37 . This is defined considering the evolution of two replicas of the system, {σ i } and {τ i } which follow identical evolution up to a time t w , and independent evolutions afterwards. One then considers the correlation that, by definition, is equal to one for −t w ≤ t ≤ 0. Barrat, Burioni and Mezard have discussed in detail the meaning of this variable proving that g = lim tw→∞ lim t→∞ Q(t, t w ) is different from zero in some domain-growth models, while it tends to zero in meanfield spin-glasses. This shows that, in the last case, the two typical trajectories explore different regions of the phase space. In equilibrium, it holds the relation Q(t/2, t w ) = C(t, t w ) (with both quantities independent of t w ). In has been also shown that the same relation holds in a trap model even in the aging situation.
We have investigated the behavior of the functions C, R and χ in numerical simulations for a large size L = 30 and temperatures T = 0.5, 0.7. The two temperatures give qualitatively similar results. In figure 7 we plot our data for C(t, t w ) and Q(t, t w ). As anticipated, the scaling τ (t w ) = t w is not obeyed. A best fit of the form τ (t w ) = t α w gives α = 0.77 and produce a fairly good collapse of the data on the same master curve. We did not investigate a possible dependence of α on the temperature. The behavior of Q(t, t w ) indicates that the parameter g is zero in this model as predicted by the mean-field theory. We see that the relation Q(t/2, t w ) = C(t, t w ) is obeyed at short times, and does not work so badly even at large times. However a rescaling of the kind Q(t/(1.5), t w ) = C(t, t w ) fits better the data.
As far as the function x(C) is concerned, we plot in figure 8 the rescaled magnetization χ = (T /h)m(t, t w ) versus C(t, t w ) for different t w . The figure shows clearly that for the waiting times here considered we are vary far from an asymptotic regime where χ is independent of t w . However, the slope of the curve for small C seems not to vary too much with t w , and it is roughly equal to x = 0.4. The region where χ(C) = 1 − C terminates roughly at a value q EA = 0.97. Assuming (as supported by mean-field theory) that the Edwards-Anderson parameter defined in this way and the one defined in statics take equal or comparable values, we argue that the quasi-states dominating the thermodynamics are quite "narrow", and, as announced, the intra-state entropy is quite small. To our knowledge, the function x(C), or equivalently χ(C) has never been measured in experiments. Its determination, which would involve independent measures of a time correlation function and its associated response function, would be a good test of the spin-glass scenario in glasses.

V. CONCLUSIONS
In this work we have investigated the glassy behavior of a short-range disordered p-spin interaction model. We have focused to the case p = 4 in three dimensions, reporting for comparison purpose some results for the two-dimensional system. The mean-field version of this model displays a static transition where the configurational entropy is nearly zero and replica symmetry breaks. We have found that the finite dimensional model shows interesting features, archetypical of real laboratory glasses. Both in two and in three dimensions the equilibrium autocorrelation function follows a stretched exponential form at low enough temperature. While in two dimensions the relaxation time increases according to the Arrhenius law, in dimension three our data demonstrate a faster increase. This is compatible both with a Vogel-Fulcher law, and with an exp(A/T 2 ) behavior. In favor of the finite-temperature transition scenario, we can only offer the extrapolations of the high temperature data of the entropy. A better support to the theory is furnished by the relation τ ∼ exp(A/(T S(T )) that we observe in the temperature window we explored. The search in several quantity, of a growing statical correlation length associated to the growing relaxation time gave negative answer. In all cases we found very small, and nearly temperature independent correlation lengths. Of course we can not exclude that the correlation length begins to grow at temperatures lower that the ones at which we could equilibrate the system.
The study of the off-equilibrium dynamics of the system confirmed the qualitative features predicted by Mean-Field Theory. The behavior of the auto-correlation function shows aging. As times goes by the dynamics becomes slower. The effective relaxation time, defined as the characteristic time for the correlation function to vanish, grows as a power of the waiting time. The analysis of the function Q(t, t w ) shows that typical trajectories which coincide at time t w explore uncorrelated regions of the space at later times. The behavior of the fluctuation-dissipation ratio, displays the qualitative features expected, showing long-range memory effects in the aging regime. However the residual dependence on t w indicates that at the times we have reached the eventual asymptotic behavior is still very far 39 .
The evolution of glassy systems is often described phenomenologically as rare jumps in a landscape of "traps" (metastable states) inside which most of the time is spent 32 . This picture matches quite well with the theoretical idea of the metastable states destabilized by activated processes. It has been pointed out in 37 that jumps among uncorrelated traps imply the relation Q(t, t w ) = C(2t, t w ). Our data show Q(t, t w ) is indeed close to C(2t, t w ). If this would be confirmed by more precise and systematic studies, would be an evidence in favor of the trap models. In this respect it would also be useful the analysis of the self-averaging properties of some local quantity.
We would like to conclude mentioning some of the problems left open by this work. For example it would be interesting to study the behavior of the dynamic Edwards-Anderson parameter as a function of temperature, or the large-time decay of the energy to its asymptotic value. We have also seen, that the equilibrium relaxation in high temperature does not seem to proceed in two steps (β and α) as it commonly observed in laboratory, and as it is predicted by Mode Coupling Theory. As far as we know, the β-relaxation process has never been observed in a spin system. We do not know if this is just due to a difficulties to disentangle the two processes due the high value of the Edwards-Anderson parameter in the neighborhood of the transition, or to the real absence of the β process. We think however that the understanding of this point can shed some light on the nature of the β-relaxation process. On the theoretical side, it would be very interesting to understand if the observed value of the exponent γ = 1 in three dimensions, in contrast with the replica theory prediction for the Potts-glass γ = 2 has to be imputed to the plane inversion symmetry, or if it is observed even in absence of this. In that respect simulations of the Potts-glass 40 could give important hints.    6. Spin-glass susceptibility χ for different sizes L = 4, 5, 6, 7 averaged over 100 samples versus the temperature. There is no evidence of a growing correlation length. Fig. 7. The functions C and Q at T = 0.5. We plot on the same graph C(t, t w ) and Q(t/2, t w ). We see that the relation C(t, t w ) = Q(t/2, t w ) is quite well obeid even for small values of C. The curves for t w = 10, 100, 1000, 10000 are depicted.