An infinite-period phase transition versus nucleation in a stochastic model of collective oscillations
A lattice model of three-state stochastic phase-coupled oscillators has been shown by Wood et al. (2006 Phys. Rev. Lett. 96 145701) to exhibit a phase transition at a critical value of the coupling parameter, leading to stable global oscillations. We show that, in the complete graph version of the model, upon further increase in the coupling, the average frequency of collective oscillations decreases until an infinite-period (IP) phase transition occurs, at which point collective oscillations cease. Above this second critical point, a macroscopic fraction of the oscillators spend most of the time in one of the three states, yielding a prototypical nonequilibrium example (without an equilibrium counterpart) in which discrete rotational () symmetry is spontaneously broken, in the absence of any absorbing state. Simulation results and nucleation arguments strongly suggest that the IP phase transition does not occur on finite-dimensional lattices with short-range interactions.
Models of coupled oscillators exhibit a surprising range of symmetry-breaking transitions to a globally synchronized state. In the paradigmatic Kuramoto model, for instance, oscillators with intrinsic frequencies coupled via their continuous phases can exhibit stable collective oscillations, with breaking of time translation symmetry [1, 2, 3, 4, 5]. When , the Kuramoto model acquires an additional reflection symmetry, which can also be broken at the onset of synchronization, yielding surprising response properties .
Within the class of discrete-phase models, the paper-scissors-stone game is an example of a system with three absorbing states which can present either global oscillations or spontaneous breaking of discrete rotational () symmetry [7, 8, 9, 10, 11]. More recently, Wood, Van den Broeck, Kawai and Lindenberg proposed a family of models of phase-coupled three-state stochastic oscillators which can undergo a phase transition to a synchronized state [12, 13, 14, 15] for sufficiently strong coupling. From here on each of them is called Wood et al.’s cyclic model (WCM). Although the WCM also has symmetry, it has no absorbing state.
As noted previously , the first WCM  has an unusual behavior after the phase transition to a synchronized state. As the phase coupling is further increased, collective oscillations slow down while oscillators remain partially synchronized. Here we have a closer look at this phenomenon, in particular addressing whether or not there is a second phase transition with breaking of symmetry.
In the WCM, the phase at site () can take one of three values: , where (for simplicity, from here on we employ to denote sites and to denote states — always modulo 3). The only allowed transitions are those from to (see Fig. 1), with transition rates
where is a constant (which can be set to unity without loss of generality), is the coupling parameter, is the number of nearest neighbors in state , and is the number of nearest neighbors. The rates are invariant under cyclic permutation of the indices, i.e., belong to the symmetry group of discrete rotations, and strongly violate the detailed balance condition.
The mean-field approximation is obtained by replacing in the argument of the exponential of Eq. 1 with the corresponding probability, , leading to the following set of equations:
where now . These also represent the equations for a complete graph in the limit : since , we can replace with in the infinite-size limit. Normalization reduces these equations to a two-dimensional flow in the plane (Fig. 2).
Starting from the complex order parameter,
one can define quantities that characterize the collective behavior of the system. A configuration with has unequal numbers of sites in the three states, which led Kuramoto to propose
as the order parameter for synchronization [1, 3, 12], where denotes a time average over a single realization (in the stationary state), and an average over independent realizations. Note that is consistent with, but does not necessarily imply, globally synchronized oscillations. The latter is characterized by a periodically varying phase of [16, 17, 18, 19].
In the mean-field equations (2), the transition to the synchronized regime is associated with a supercritical Hopf bifurcation at . The trivial fixed point loses stability at , and a limit cycle encircling this point appears (Fig. 2a-b) . For , sustained oscillations in characterize synchronization among the oscillators (Fig. 3a). Correspondingly, grows continuously at the transition (Fig. 3c), with a mean-field exponent . The scaled variance diverges with the system size at criticality, as shown for simulations of complete graphs using in Fig. 3d.
Note that the phase transition to stable collective oscillations is associated with the breaking of the continuous time-translation symmetry: the time series , and become periodic for . But since they are statistically identical, except for a phase, symmetry still holds. Increasing above enhances synchronization among the oscillators, leading to increasing oscillation amplitudes, as shown in Figs. 2d and 3b.
3 The Second Phase Transition
3.1 Mean-field theory and the complete graph
Wood et al. found that the increasing amplitude of oscillation is accompanied by a decreasing angular frequency , where is the mean time between peaks in (Figs. 3a-c). This can be understood qualitatively from the exponential dependence of the rates of Eq. 1 on the neighbor fractions: when a state is highly populated, the rate at which oscillators leave this state becomes very small. In addition, the rate at which a highly populated state attracts oscillators from its predecessor is also very high.
We study the effect of a further increase in the coupling , beyond the range of values investigated in Ref. , using mean-field theory and analysis on a complete graph of sites. Note that the state of the latter is completely specified by two integer random variables, and . We study the complete graph (for up to 1000) via exact (numerical) stationary solution of the master equation  and (for ranging from 10 to 10) via event-driven Monte Carlo simulations. The results of the two approaches are fully consistent.
In mean-field theory, when reaches an upper critical value , three symmetric saddle-node bifurcations occur simultaneously, and the period of the collective oscillations diverges. Above , there are three symmetric attractors in the system, and 3-fold rotational () symmetry is spontaneously broken: Figure 2f illustrates how the final fate of three different trajectories starting from the same initial conditions is determined by fluctuations. Moreover, the transition can be considered reentrant, in the sense that time-translation symmetry, which had been broken at , is restored at (although, unlike other cases seen in the literature [21, 22], the phase above is not exactly the same as that below , because of symmetry breaking).
Analogously to what occurs in a condensation or ferromagnetic phase transition , the freezing of the relative occupancy of each state does not imply that individual oscillators freeze as well. The frequencies of individual oscillators do decrease with increasing , but only vanish in the limit , when one of the states is fully occupied.
It is convenient to define an order parameter that is identically zero (in the infinite-size limit) for . Since the angular frequency does not fulfill the usual requirements for an order parameter , we propose
where is the Kronecker delta and is the transition rate at site (see Eq. 1). Thus reflects not only the configuration, but the rate at which the latter evolves.
On the complete graph, is the same for all sites in the same state . Denoting this rate by , the order parameter can be rewritten for mean-field (MF) analysis as
Note that both the disordered phase () as well as the frozen phase () have stable fixed points, so from Eq. 2. This in turn renders in Eq. 6 for both cases [a similar line of reasoning can be applied directly to Eq. (5)].
Figure 4 shows the variation of () with in mean-field theory, and on the complete graph, confirming that functions as an order parameter to detect both phase transitions. The MF critical behavior is for and for . In the pair approximation for a lattice with coordination number we find and ; for the corresponding values are and . In the pair approximation, no transition is observed for .
On the complete graph, the order parameter decays with system size as at and as at , as shown in Fig. 5. The first result is typical of mean-field-like scaling with system size at a continuous phase transition, and can be understood as follows. Our numerical results show that with increasing , as expected, the stationary probability distribution in the plane becomes (for ) increasingly concentrated around the symmetric point . The maximum of the distribution, , scales as and the number of points at which the probability is relatively large (i.e., ) grows . In terms of the densities , the probability distribution is concentrated on a region of radius . Note that is a quadratic form of the variables and that for . Thus for the purposes of scaling analysis we may write the stationary average of as
with . At the first transition () the same argument applies to the order parameter , and we indeed find numerically .
At the second transition, the probability distribution is peaked near the three stable fixed points of the mean-field analysis, and so has a very different form than that found at . It is therefore not surprising that decays with a different exponent. The value of 0.4203(3) (which does not correspond to any simple ratio of integers) nevertheless remains something of a puzzle. We suspect that understanding this result will require an asymptotic (large-) analysis of the Fokker-Planck equation for the probability density in the plane, a task we defer to future work.
3.2 Hypercubic lattices
Motivated by the appearance of a novel phase transition in mean-field theory and on the complete graph, we search for such a transition on finite-dimensional structures, i.e., the simple cubic lattice and its four-dimensional (hypercubic) analog, in systems of sites. On these lattices, we do observe a marked tendency for to fall rapidly with for values well above the first transition at (, and for , 4 and 5, respectively ). We also observe a broad maximum in at some point in this range. To infer the existence of a phase transition, rather than a rapid but smooth variation of and other properties with , we require evidence of emergent singularities as the linear system size is increased.
The global dynamics may be characterized in terms of a process that begins with all sites in state , and ends when the fraction of sites remaining this state falls to . Let the “escape time” denote the mean duration this process, averaged over realizations of the stochastic dynamics. If there is an infinite-period phase transition at some critical value , then should grow at this point, where denotes the dynamic critical exponent. Although simulations reveal an exponential growth of with at fixed , they show that for large, fixed , the escape time decreases with increasing system size! The mean period (see Fig. 6) exhibits two regimes for : for smaller , relatively slow, exponential growth with essentially independent of system size, and, for larger , rapid exponential growth with , with . We identify the latter regime as that of single-droplet nucleation . In the former regime, the density of “advanced” sites (i.e., in state in a configuration with the majority in state ) is large enough that essentially no barrier to the cyclic dynamics exists. The (smooth) crossover between these regimes occurs at a larger value for larger system sizes.
The results for the mean transition rate on the four-dimensional hypercubic lattice, , clearly suggest a nucleation mechanism. The exponential factor represents a barrier to nucleation, that is, the mean time to formation of a critical cluster. (By a “critical” cluster we mean one that is equally likely to grow as to diminish in size; larger clusters tend to grow while smaller ones tend to shrink.) Consider a region in which all spins are in the same state. The rate at which a given spin flips to the next state (e.g., 1 2) is leading to the exponential growth in . The prefactor corresponds to essentially independent contributions due to each site in the system. (Given the extremely low density of clusters in this regime, we can treat the nucleation and initial growth of a cluster as occurring in isolation, without interference from other clusters.)
The constant in the exponential reflects the probability of formation of a critical cluster, which in turn depends on the critical cluster size, . Intuitively, the nucleation barrier exists because an isolated advanced site has a higher rate to advance once again, and then rapidly return to the majority state, than it does of inducing a neighbor to advance. Given the difficulty of visualizing four-dimensional space, we begin by studying nucleation on the square lattice. Of course, the model does not exhibit an ordered phase in two dimensions, but this issue can be bypassed by using small systems (a lattice of 1010 sites, in the present case). Then the system remains near one of the three attractors (almost all sites in the same state) except for the rapid transitional periods following a nucleation event. (In larger systems global synchronization is quickly lost and the system breaks into a number of domains.) On such a small lattice, one may easily identify transitional regimes between, for example, states with majority-1 and majority-2. We define the transition (in the macroscopic sense) as the moment when 20% of sites are in state 1, with nearly all the rest (70% or more) in state 2. To estimate we determine the time required for a large number (4000, in practice) of transitions to occur, in event-driven Monte Carlo simulations. A typical evolution is shown in Fig. 7. We find that the mean transition time grows exponentially with ; for we find . This suggests that the critical cluster size is approximately 3 - 4.
Several other pieces of evidence support this estimate. First, in the pre-transition regime, the number of sites in state 2 is seen to fluctuate between zero and three or so. In the example shown in Fig. 7, the second time grows to four, a transition occurs. Typical formats of a growing cluster are shown in Fig. 8; the cluster is seen to be compact, as expected for a growth probability that is highest for the sites outside the cluster having the largest number of neighbors within it.
A semiquantitative analysis of the cluster size distribution in the pretransition regime can be developed as follows. Consider first the density of isolated sites in state 2, in a background of 1s. The rate for a given 1 to change to 2 is , while an isolated 2 (with all neighbors in state 1) flips to state 3 at unit rate. There are also transitions of the form 21 22. The total rate for an isolated 2 to become a dimer is . Thus we have a rate equation of the form
leading to in the quasi-stationary (QS) regime. By the QS regime we mean the long-lived metastable state in which the great majority of sites are in state 1. Similarly, let denote the probability that a given site is the leftmost (or lowermost) site of an isolated dimer (the restriction is to avoid double counting). Then we have
so that in the QS state. Consider next the density of L-shaped trimers. The gain term in the equation for is while the dominant loss term (due to formation of a four-site, square cluster) is simply (the transition rate is unity). Thus . The exponential factor, , is quite close to the factor found in simulations, suggesting that this is essentially the critical cluster. (Note that this is the smallest cluster with a larger probability of growing than of shrinking.) Applied to the simple cubic lattice, this line of reasoning suggests that the critical cluster consists of seven sites (in the form of a cube with one corner removed) and yields . [A crude estimate for the value in four dimensions can be obtained by supposing that ; this yields , reasonably close to the simulation value, .]
The nucleation scenario outlined above essentially rules out the second transition on any finite-dimensional lattice with short-range interactions. There is always a cluster size such that growth is more likely than shrinkage, and once such a cluster forms, a transition of the global state is likely. From the perspective of equilibrium phase transitions, this is not surprising: the ordered state of any (finite) Ising or Potts model is always subject to a reversal of magnetization at a finite rate. (This rate, however, decays with system size in equilibrium, as a finite cluster is always more likely to shrink than to grow, in the absence of an external field.) In the present cyclic dynamics, by contrast, there is no impediment to growth (no “surface tension”, as it were), once a cluster has attained the critical size.
We have identified a second phase transition in the cyclic model proposed by Wood et al. [12, 13, 14, 15] This continuous transition is of a fundamentally different nature than the first, as it corresponds to an infinite-period transition at which macroscopic properties freeze and a discrete rotational () symmetry is spontaneously broken, in the absence of an absorbing state. We find that an appropriate order parameter for the second transition involves the local transition rates. While the phase transition is observed on the complete graph, we argue that it does not exist on finite-dimensional structures with short-range interactions. The question naturally arises whether such an infinite-period transition exists in systems with long-range interactions. This indeed appears likely, for interactions that decay sufficiently slowly with distance. In the case of graphs with interactions between many pairs of sites, independent of distance, we expect to observe an infinite-period transition if (in the infinite-size limit) a given site interacts with a finite fraction of the others. We plan to investigate these questions, as well as the behavior on scale-free networks, in future work.
VRVA and MC acknowledge financial support from CNPq, FACEPE, CAPES,
FAPERJ and special programs PRONEX, PRONEM and INCEMAQ. RD
acknowledges financial support from CNPq.
-  Y. Kuramoto, editor. Chemical Oscillations, Waves and Turbulence. Dover, Berlin, 1984.
-  SH Strogatz and I Stewart. Coupled oscillators and biological synchronization. Sci. Am., 269(6):102–109, 1993.
-  S. H. Strogatz. From Kuramoto to Crawford: Exploring the onset of synchronization in populations of coupled oscillators. Physica D, 143:1–20, 2000.
-  S. H. Strogatz. Sync: How Order Emerges From Chaos In the Universe, Nature, and Daily Life. Hyperion, New York, 2004.
-  A. Pikovsky, M. Rosenblum, and J. Kurths. Synchronization: A Universal Concept in Nonlinear Sciences. Cambridge University Press, Cambridge, UK, 2001.
-  P. Reimann, C. Van den Broeck, and R. Kawai. Nonequilibrium noise in coupled phase oscillators. Phys. Rev. E, 60(6):6402–6406, 1999.
-  K. Tainaka. Lattice model for the Lotka-Volterra system. J. Phys. Soc. Japan, 57(8):2588–2590, 1988.
-  K. Tainaka. Stationary pattern of vortices or strings in biological systems: Lattice version of the Lotka-Volterra model. Phys. Rev. Lett., 63(24):2688–2691, 1989.
-  K. Tainaka and Itoh Y. Topological phase transition in biological ecosystems. Europhys. Lett., 15(4):399–404, 1991.
-  Y. Itoh and K. Tainaka. Stochastic limit cycle with power-law spectrum. Phys. Lett. A, 189:37–42, 1994.
-  K. Tainaka. Vortices and strings in a model ecosystem. Phys. Rev. E, 50(5):3401, 1994.
-  K. Wood, C. Van den Broeck, R. Kawai, and K. Lindenberg. Universality of synchrony: Critical behavior in a discrete model of stochastic phase-coupled oscillators. Phys. Rev. Lett., 96:145701, 2006.
-  K. Wood, C. Van den Broeck, R. Kawai, and K. Lindenberg. Critical behavior and synchronization of discrete stochastic phase-coupled oscillators. Phys. Rev. E, 74(3):031113, 2006.
-  K. Wood, C. Van den Broeck, R. Kawai, and K. Lindenberg. Effects of disorder on synchronization of discrete phase-coupled oscillators. Phys. Rev. E, 75(4):041107, 2007.
-  K. Wood, C. Van den Broeck, R. Kawai, and K. Lindenberg. Continuous and discontinuous phase transitions and partial synchronization in stochastic three-state oscillators. Phys. Rev. E, 76(4):041132, 2007.
-  Y. Kuramoto, T. Aoyagi, I. Nishikawa, T. Chawanya, and K. Okuda. Neural Network Model Carrying Phase Information with Application to Collective Dynamics. Prog. Theor. Phys., 87:1119–1126, 1992.
-  H. Ohta and S. Sasa. Critical phenomena in globally coupled excitable elements. Phys. Rev. E, 78(6):065101, 2008.
-  S. Shinomoto and Y. Kuramoto. Phase transitions in active rotator systems. Progr. Theoret. Phys., 75(5):1105–1110, 1986.
-  F. Rozenblit and M. Copelli. Collective oscillations of excitable elements: order parameters, bistability and the role of stochasticity. J. Stat. Mech., 2011:P01012, 2011.
-  R. Dickman. Numerical analysis of the master equation. Phys. Rev. E, 65(4):047701, 2002.
-  C. Van den Broeck, J. M. R. Parrondo, R. Toral, and R. Kawai. Nonequilibrium phase transitions induced by multiplicative noise. Phys. Rev. E, 55(4):4084–4094, 1997.
-  R. Kawai, X. Sailer, L. Schimansky-Geier, and C. Van den Broeck. Macroscopic limit cycle via pure noise-induced phase transitions. Phys. Rev. E, 69(5):051104, 2004.
-  K. Huang. Statistical Mechanics. John Wiley & Sons, 1987.
-  M. Plischke and B. Bergersen. Equilibrium Statistical Physics. World Scientific Publishing Company, 2nd edition, 1994.
-  P. A. Rikvold, H. Tomita, S. Miyashita, and S. W. Sides. Metastable lifetimes in a kinetic ising model: Dependence on field and system size. Phys. Rev. E, 49(6):5080–5090, 1994.