Slow dynamics of Ising models with energy barriers

Using Monte Carlo simulations we study the dynamics of three-dimensional Ising models with nearest-, next-nearest-, and four-spin ~plaquette! interactions. During coarsening, such models develop growing energy barriers, which leads to very slow dynamics at low temperature. As already reported, the model with only the plaquette interaction exhibits some of the features characteristic of ordinary glasses: strong metastability of the supercooled liquid, a weak increase of the characteristic length under cooling, stretched-exponential relaxation, and aging. The addition of two-spin interactions, in general, destroys such behavior: the liquid phase loses metastability and the slow-dynamics regime terminates well below the melting transition, which is presumably related with a certain corner-rounding transition. However, for a particular choice of interaction constants, when the ground state is strongly degenerate, our simulations suggest that the slow-dynamics regime extends up to the melting transition. The analysis of these models leads us to the conjecture that in the four-spin Ising model domain walls lose their tension at the glassy transition and that they are basically tensionless in the glassy phase.


I. INTRODUCTION
A lot of effort has been devoted in the past twenty years to understanding the behavior of various glassy and disordered systems ͓1͔.Such systems, which include conventional glasses, spin glasses, amorphous semiconductors, and many others are of great importance both experimentally and theoretically.However, despite intensive research, our understanding of such systems is still limited.For example, even the very nature of the glassy phase in spin glasses is still a very controversial issue ͓2-4͔.Although they are much more abundant, conventional glasses seem to pose an even greater puzzle.Why do supercooled liquids fall out of equilibrium at a more or less well defined temperature?Why do they collapse into the glassy state when the cooling is fast enough and into the crystalline phase when the cooling is slow?These fundamental questions still await definitive answers.One of the important problems in physics of conventional glasses is the continuing lack of a satisfactory microscopic model of such systems.In this respect the situation is much better for spin glasses where it is commonly accepted that models containing quenched disorder correctly describe physics of such systems.Lattice realizations of such models are a particularly valuable source of information about spin glasses ͓5͔.The most realistic models of conventional glasses, so-called off-lattice models, still constitute an enormous computational challenge although progress in this field is also being made ͓6͔.
A model of conventional glasses should be capable of describing ͑at least͒ three phases: liquid, glass, and crystal.The actual state of the system should be determined by control parameter͑s͒ ͑e.g., temperature͒ and possibly also its history.Since the glass is regarded as a liquid trapped during the falling out of equilibrium, the model should possess such a trapping mechanism.In spin glasses the trapping mechanism is related with energy barriers generated by quenched disorder ͓3͔.On the other hand, we do not expect the quenched disorder to be a relevant factor in conventional glasses because models with strong quenched disorder are unlikely to exhibit periodic solutions ͑which are needed for the model to be in the crystal phase͒.Recently, various lattice models, which do not contain quenched disorder, were studied which have some features of conventional glasses.Some of these models are infinite-dimensional and their thermodynamical properties can be found exactly ͓7͔.There are also finite dimensional models whose dynamics exhibit some glassy behavior ͓8,9͔.
Recently, it has been shown that the three dimensional Ising model with the four-spin ͑plaquette͒ interaction also exhibits some glassy features ͓10-12͔.This model undergoes a first-order phase transition between low-temperature ͑crystal͒ and high-temperature ͑liquid͒ phases.However, when conventional simulation techniques are used, the transition is screened by a very strong metastability during heating as well as cooling.For temperatures lower than the limit of metastability of the liquid phase, the model has a very slow coarsening dynamics.In addition, the zero-temperature characteristic length increases very slowly as a function of the inverse cooling rate, which is also an expected property of glasses.Further evidence of the glassy behavior in this model has been recently reported by Swift et al. ͓13͔.They have shown that the glassy transition coincides with the divergence of a certain relaxation time and that aging properties of the model are also typical of glassy systems.They have also observed that some time correlation functions may decay as stretched exponentials.
These results strongly suggest that the model with fourspin interactions might describe important aspects of the glassy transition.It would be interesting to find which properties of this model are responsible for such a behavior.It has been already suggested ͓10-12͔ that the trapping mechanism might be related with diverging energy barriers.These bar-riers would arise in this model basically due to the same mechanism as in a model with competing nearest-neighbor and next-nearest-neighbor interactions examined by Shore et al. ͓8͔ ͑the SS model for short͒.However, the behavior of the SS model is not fully consistent with our conception of glasses since it orders too quickly under cooling ͓8͔.It was also suggested that the difference in the behavior of the SS and four-spin models might be related with the degeneracy of the ground state in the four-spin model.This degeneracy might lead to some entropy barriers, which would be responsible for the strong metastability of the liquid phase.
In the present paper, using Monte Carlo simulations, we examine a certain class of three-dimensional Ising models which generate energy barriers.These models are described by the following Hamiltonian:

͑1͒
In the above expression ͗...͘ and ͗͗•••͘͘ denote pairs of nearest and next-nearest neighbors, respectively, and ͓i, j,k,l͔ stands for summation over elementary plaquettes.In general, these models have double degenerate ground state and our simulations suggest that the dynamical properties in this case are similar to the SS model.However, when the interaction constants are such that the model has a strongly degenerate ground state ͑gonihedric case͒, the dynamical properties change.Simulations suggest that two types of dynamical behavior appear.In the first type the model behaves similarly to the already described four-spin model.In the second type, the glassy transition appears to coincide with the thermodynamic transition.
Such behavior gives rise to the following questions: why a glassy transition appears only in certain systems with slow dynamics and what is its nature.Analysis of the ground-state structure and thermodynamic properties of models studied here prompts the following conjecture, which, if confirmed, would constitute an important result of the present paper: at the glassy transition the domain walls lose their surface tension, and, as a result, the glassy phase is composed of tensionless domains.Although based on the analysis of Ising models, we hope that such an interpretation might shed some light on the nature of the glassy transition in more realistic systems also.Moreover, such an interpretation of glassy phase is in accord with some recent hypotheses concerning the nature of the glassy phase in spin glasses ͓4͔.
In Sec.II we discuss briefly the properties of the four-spin model.In Sec.III we present the results of our simulations for the gonihedric case.The case of the doubly-degenerate ground state is discussed in Sec.IV.Section V contains a summary of our results and some arguments on the nature of the glassy transition.

II. ISING MODEL WITH PLAQUETTE INTERACTIONS
This model corresponds to the case J 1 ϭJ 2 ϭ0,J 4 ϭ1, and has been already studied using cluster variational method ͓14͔ and Monte Carlo simulations ͓15,10-13͔.
Clearly, the ferromagnetic configuration is a ground-state configuration of this model.It is also easy to realize that flipping coplanar spins does not change the energy.Thus any configuration obtained from the ferromagnetic configuration by flipping coplanar spins is also a ground-state configuration.Moreover, any combination of such coplanar flippings ͑even for crossing planes͒ does not increase the energy.Simple analysis along these lines shows that for the model on the lattice of the linear size L the degeneracy of the ground state is equal to 2 3L .Although ground state of this model is strongly degenerate its ground-state entropy is zero.
The model undergoes a first-order thermodynamic transition at TϭT c ϳ3.6 which is, however, screened by very strong metastability ͓11͔.As a result, when heated or cooled, the transition observed in simulations is shifted to Tϳ3.9 or Tϳ3.4,respectively.Transitions at these spinodals are accompanied by peaks in the specific heat.
The low-temperature spinodal Tϳ3.4 seems to coincide with the glassy transition.Below this temperature the model exhibits very slow coarsening dynamics ͓11͔ as well as aging properties which are characteristic of glassy systems ͓13͔.A certain characteristic time, which governs the relaxation of energy-energy correlation functions, also seems to diverge at this temperature ͓13͔.In addition, the behavior of the model under continuous cooling supports the glassy-transition interpretation of this temperature ͓12͔.

A. Ground state and thermodynamics
It has already been suggested that the slow dynamics of the four-spin model might be related with energy barriers generated in that model ͓10,11͔.These barriers arise due to the shape dependence of the energy of excitations: it is not only the size of an excitation which determines its energy but also its shape.Such shape dependence appears also in the SS model.
Are there any other models which could have a similar property?In our opinion, the shape dependence of energy of excitations should be rather a robust feature of Ising-type models.It is only in some specific cases, like the standard nearest-neighbor case, when this energy does not depend on the shape of an excitation.In particular, energy barriers appear in model ͑1͒.The Hamiltonian of this model is quite general and it includes both the four-spin model and the Shore et al.'s model (J 4 ϭ0,J 1 Ͼ0,J 2 Ͻ0).
In the present section we examine a class of models described by this Hamiltonian, namely gonihedric models ͓16,17͔.These models correspond to the following choice of interaction constants: J 1 ϭ2k,J 2 ϭϪk/2 and J 4 ϭ 1 2 (kϪ1).Gonihedric models have a strongly degenerate ground state.In addition to the ferromagnetic ground state any configuration obtained by flipping coplanar spins also minimizes the Hamiltonian.Any combination of such flips does not increases the energy, provided that flipping planes do not cross.As a particular example of such a ground state we can mention lamellar configurations where, e.g., every second plane of spins is flipped.Although lamellar structures constitute a legitimate ground state, they do not survive at finite temperature as shown by Cirillo et al. using the cluster variational method ͓14͔.We will return to this feature in the last section.
For kϭ0 the gonihedric model is equivalent to the fourspin Ising model.In this case the model has an additional symmetry which implies a larger degeneracy of the ground state since the flipping planes can now cross.As a result we obtain that antiferromagnetic configurations belong to the ground state.Further analysis of differences between the k ϭ0 and the k 0 cases is postponed to the last section.
Gonihedric models are expected to undergo a thermodynamic transition which for kϽk tr is of first order and for k Ͼk tr is of second order.Only very rough estimations of k tr (ϳ0.5) are known ͓14͔.In this section we analyze dynamical properties of the gonihedric model for kϭ2, i.e., for a value with a continuous transition.Our results were obtained using a standard Monte Carlo method with random sequential update using Metropolis algorithm ͓18͔.Some details can be found elsewhere ͓11,12͔.
To find the thermodynamic transitions we measured the specific heat.Our simulations, which were made for various linear sizes L up to Lϭ40, locate the peak at T c ϳ2.35, which is a good agreement with the cluster variational method estimation ͓14͔.The absence of hysteresis effects confirms that the transition at TϭT c is continuous.͑Although the nature of the thermodynamic transition in gonihedric models is an interesting and still open problem, its further analysis is not an objective of the present paper.͒

B. Dynamics
An important indication of glassy dynamics is a slow evolution of a random quench.For usual models with nonconservative dynamics one expects ͓19͔ that the characteristic length l increases with time t as lϳt 1/2 .However, in glassy systems l should increase much more slowly and presumably only logarithmically with time (lϳln t).Such a behavior most likely appears in the SS model and in the four-spin model.In the following we present the results of our simulations of the evolution of quenches in the kϭ2 model.We With this identification from Fig. 1 we infer that for all the examined temperatures the asymptotic increase of l is much slower than t 1/2 .At the end of this section we will argue that the relation ͑2͒ most likely does not hold for this model and Fig. 1 actually suggest that the increase of l is slower than t 1/4 .Since there are no theoretical arguments for such a slow algebraic increase in our opinion it is quite plausible that asymptotically we have lϳln t.Such a slow increase of l is most likely due to energy barriers.An alternative technique to examine the dynamics of the model is to measure characteristic times of certain processes.For example, one can measure the average time needed for the inversion of a cubic like excitation.In this method, which was already applied to similar models ͓8,11͔, one prepares the system of the size L with fixed boundary conditions and interior spins which are opposite to the boundary spins.One expects that after some time, the system will invert the interior spins.For a two-spin Ising model or SS model above the corner rounding transition, ϳL 2 , which indicates a relatively fast dynamics ͑naive inversion of this relation gives l ϳt 1/2 ).On the other hand, in the SS model below the corner rounding transition and in the four-spin model increases much faster, presumably exponentially, with L.
We measured the time needed for magnetization of the interior spins to reach the equilibrium value at a given temperature and the results are shown in Fig. 2.These results show that even at the highest examined temperature increases faster than L 2 .This is a potential indication of an exponential increase ϳa L (aϾ1) in the entire lowtemperature phase.
Additional confirmation of such behavior is obtained from simulations of this model under the continuous cooling.Similarly to simulations of the four-spin model ͓12͔, we relax the random sample at a temperature T 0 ϾT c and then FIG. 1.The excess energy ␦E as a function of t in the log-log scale.Simulations were made for the system size Lϭ100 and T ϭ1.0, 1.3, 1.6, and 1.9 ͑from bottom to the top͒.The dotted line has a slope 0.5.FIG. 2. The characterisitc time needed for the inversion of the cubic excitation of the size L as a function of L for the gonihedric case (kϭ2).Simulations were made for Tϭ1.9 (᭺), 2.1 (ᮀ) and 2.3 (ٌ).The dotted line has a slope corresponding to ϳL 2 .continuously lower the temperature according to the formula T(t)ϭT 0 Ϫrt, where r is the cooling rate.When the temperature is reduced below the critical point T c the growth of order begins.The slower the cooling the more ordered is the system at the end of the cooling , i.e., at Tϭ0 ͑see Fig. 3͒.To quantify the zero-temperature order we measure the excess energy ␦E at Tϭ0 and the result is shown in Fig. 4.
Using the relation ͑2͒ this data suggests that asymptotically lϳr Ϫ1/2 .Such a relation, which indicates that the growth of order is relatively fast, holds for the two-spin Ising model ͓20͔ and also for the SS model ͓8͔.However, this conclusion is based on the validity of the relation ͑2͒ and, similarly to the four-spin model ͓11͔, we want to argue that this relation does not hold.Our argument refers to the following property of all gonihedric models: the energy of cubiclike excitations scales as their linear size L ͓21͔.Let us recall that in two-spin Ising model, this energy scales as the area of the excitation (ϳL 2 ).Provided that the final configuration is composed of such domains of the size L and using a simple dimensional argument ͓12͔ we obtain

͑3͒
A direct confirmation of the assumption about the structure of the configuration at the of the cooling process comes from visual inspection.In Fig. 5 we can see an example of singlelayer configuration.One can clearly see cubiclike ͑i.e., non-rough͒ domains whose energy scales linearly with their size.
Using the relation ͑3͒ Fig. 4 shows that the zero-temperature characteristic length scales as r Ϫ1/4 which is much slower than in the two-spin Ising model but faster than in the fourspin model.
In this section we used three independent techniques to probe the dynamics of the gonihedric model in the case of continuous thermodynamic transition kϭ2.Domain coarsening suggest that the model has a slow dynamics up to, at least, the temperature Tϭ1.9 (T c ϳ2. 35).An analysis of the size dependence of the characteristic time suggests that cubiclike domains remain nonrough at least up to the temperature Tϭ2.1.Thus a slow-dynamics regime is most likely extended up to this temperature.Since this is very close to the critical point is it not unlikely that a slow-dynamics regime actually covers the whole low-temperature phase.The behavior of the model under cooling confirms such a scenario: if there would be a certain temperature T 0 ϽT c such that for T 0 ϽTϽT c the dynamics would be fast then for the slow cooling the growth of order would be dominated by the time spent in this temperature interval and we would have lϳr Ϫ1/2 .Such a scenario takes place in the SS model ͓8͔.
The growth of order in our kϭ2 model is much slower l ϳr Ϫ1/4 and excludes the existence of such a temperature T 0 ͑unless it is very close to T c and our simulations are not sufficient to detect the true asymptotic behavior͒.Let us also notice that Shore et al. also analyzed certain SOS model for which ͑by necessity͒ T 0 ϭT c .Using some scaling arguments they have shown that for this model one should have l ϳr Ϫ1/4 , which is in agreement with our numerical result.

IV. OFF-GONIHEDRIC ISING MODEL
The gonihedric case corresponds to a certain choice of interaction constants in the Hamiltonian ͑1͒.As we have noted this choice has important implications: the ground state is strongly degenerate and energy of excitations scale as their linear size and not as their area.In the present section we examine what is going on when the interaction constants of model ͑1͒ deviate from the gonihedric case ͓22͔.
As a particular example we choose J 1 ϭ6, J 2 ϭϪ1, J 4 ϭ1/2, which differs from the gonihedric case kϭ2 by a modified nearest-neighbor coupling J 1 .Such a model has a double-degenerate ͑ferromagnetic͒ ground state and our rough estimation of the critical temperature is T c ϳ12.5.
To examine the dynamics of this model we used the same techniques as described in the previous section.First, we examined the coarsening behavior of this model.At low temperature ͑up to Tϳ3.0) we observed a very slow decrease of the energy toward the ground state value.Our data, which we do not present suggests that for such temperatures l most likely increases logarithmically with time.However, above this temperature the dynamics becomes much faster and presumably the characteristic length increases as lϳt 1/2 .Such behavior is confirmed by the measurements of the characteristic time defined in the same way as in the previous section.Results of simulations are shown in Fig. 6.They indicate that for temperature Tϭ6 and 8 the characteristic time increases as L 2 similarly to the SS model above the corner rounding transition.
The above results indicate that the dynamical behavior of the model is very similar to the SS model.Namely, in the low temperature regime the model has a slow dynamics and rapidly ͑faster than L 2 ) increasing characteristic time .However, within the ordered phase ͑i.e., for TϽT c ϳ12.5) there is also a high temperature regime where dynamics is much faster.Presumably, in this regime the dynamics is similar to other nonconservative systems with scalar order parameter ͓19͔.Similarity between the model examined in this section and the SS model is in our opinion related with the structure of the ground state: both models have only double degenerate ͑ferromagnetic͒ ground states.It is well known in the statistical mechanics that degeneracy of the ground state plays a very important role in determining the thermodynamic behavior of models ͑e.g., critical behavior͒.
Most likely this is an important factor for determining the dynamical behavior of models also.

V. NATURE OF THE GLASSY PHASE
In the present paper we have examined the dynamics of models described by the Hamiltonian ͑1͒.Depending on the interaction constants we can distinguish three types of behavior which are schematically shown in Fig. 7.
(a) Double-degenerate ground state.This, the most typical situation, appears in the off-gonihedric model studied in previous section and also in the SS model whose dynamics has been already examined in great details ͓8͔.The dynamics of the model a low temperatures (TϽT c ) has two regimes separated by a certain temperature T cr .For TϽT cr the model has slow dynamics with most likely logarithmically increasing characteristic length l.Such behavior is related with the fact that at such temperature the model is below the cornerroughening transition ͓23͔.As a result, an evolving quench develops complicated structures of cubiclike excitations which are very stable and effectively block further coarsening dynamics ͓24͔.At TϭT cr the model undergoes the corner rounding transition and the blocking mechanism is no longer effective.As a result the fast ͑standard͒ dynamics is restored.
(b) Gonihedric case with continuous transition.In this case the entire low-temperature phase has slow dynamics, whose origin is similar to the case ͑a͒.Namely, a quench develops cubiclike structures which block further coarsening dynamics.The degeneracy of the ground state seems to be the most important difference between this and the offgonihedric case.Thus, we relate the disappearance of the corner-rounding transition ͑or maybe its overlap with T c ) with the infinite degeneracy of the ground state.
(c) Four-spin model.As we already mentioned, the fourspin model undergoes a dynamic transition which exhibits a lot of similarities with the glassy transition.It might be interesting to examine whether such a behavior appears only at kϭ0 or persists also for some other ͑small͒ values of k.
Glassy transition: loss of surface tension.Why does this transition exist in the four-spin case and not in the other cases?It was already suggested that the difference between the four-spin model and the SS model is related to the degeneracy of the ground state ͓12͔.However, the above analysis of the gonihedric model with kϭ2 shows that the infinite degeneracy of the ground state is not sufficient for the model to have a glassy transition ͑of course, we limit our analysis to models which can generate diverging energy barriers and thus have slow coarsening dynamics͒.Why does the gonihedric case kϭ2 differ from the kϭ0 case ͑i.e., the four-spin model͒?Both models have strongly degenerate ground state.The degeneracy equals 2 3L for kϭ0 and 2 L for kϭ2.Since in both cases degeneracy increase exponentially with the linear system size this difference does not seem important.In our opinion, however, the difference in the dynamical behavior is related with the ground state structure of these models.As we already mentioned, for kϭ0 the flipping planes, which generate various ground-state configurations might cross.As a result, in addition to ferromagneticlike configurations we obtain antiferromagneticlike ones.For the k 0 such crossings are not allowed and only ferromagneticlike configurations are possible ͑we consider lamellar configurations also as ferromagnetic-like͒.
This difference has important implications.Let us note that for the four-spin model in addition to tensionless domain walls, which appear for example when a cubic ferromagnetic ''up'' domain is surrounded by ''down'' one, there are tensionful ones too.As an example of such a domain wall we can consider an antiferromagnetic domain surrounded by ferromagnetic one ͓12͔.At first sight it does not seem to be much different from the kϭ2 case.Indeed, when one considers a lamellar configuration where successive layers are of opposite sign ͑see Fig. 8͒ which is surrounded by a ferromagnetic domain than the excess energy of such a configuration scales as the area of the wall ͑i.e., the domain wall is tensionful͒.There is, however, an important question: do such configurations affect coarsening dynamics or, in other words, are they spontaneously generated in sufficient amounts?In our opinion the answer to this question is negative.Our first argument comes simply from the visual inspection of the snapshot configuration.In Fig. 5 one can see relatively large ferromagneticlike domains but there is no indication of lamellar ones.The second argument comes from the cluster variational method calculations ͓14͔ by Cirillo et al. who have shown that the lamellar structures are equivalent ͑i.e., of the same energy͒ to the ferromagnetic one but only for the ground state.At nonzero temperatures they are always metastable.These arguments show why such configurations are not spontaneously generated during the evolution of the quench.They imply that for kϭ2 the dominant domain walls which exist at the late-time evolution of the quench are tensionless.
On the other hand, for the four-spin model, antiferromagnetic structures are fully equivalent to the ferromagnetic ones.In the liquid phase both ferromagnetic and antiferromagnetic domains are intertwinned and form very complicated structures.As already noticed ͓12͔, tensionful domains usually have lowest energy barriers and the system can relatively easily remove the interior domains.On the contrary, tensionless domain walls have large energy barriers and their dynamics is much slower.It means that in the liquid phase dynamics is dominated by dynamics of tensionful domains and thus resembles the dynamics of two-spin Ising models.
It is in our opinion very likely that upon lowering the temperature the system will undergo a phase transition which will eliminate tensionful domain walls.Below that transition the energy of the system would be located mainly in tensionless domain walls.It means that at this transition the antiferromagnetic-ferromagnetic symmetry of the model would be spontaneously broken.In other words, at this transition the system selects a dominant type of domains, whether ferromagnetic or antiferromagnetic.This type of symmetry breaking is an essential ingredient of the transition which we tentatively identify as the glassy transition.
Finally, let us note that the idea that the glassy phase consists of a complicated mixture of tensionless domain walls appeared recently in the contex of spin glasses ͓4͔.It suggests that, at least at the geometrical level, spin glasses and ordinary glasses might exhibit a lot of similarities.Their further explorations is, however, left as a future problem.
measured the excess energy ␦EϭEϪE ϱ , where E ϱ is the equilibrium energy.Our results for temperatures TϽT c are shown in Fig. 1.To relate the characteristic length with the energy excess we can employ the frequently used relation ͓8,11͔ lϳ1/␦E.͑2͒

FIG. 6 .
FIG.6.The characteristic time needed for the inversion of the cubic excitation of the size L as a function of L for the offgonihedric case.Simulations were made for Tϭ2 (छ), 3 (᭝), 6 (ᮀ), and 8 (᭺).The dotted line has a slope corresponding to ϳL 2 .