This paper was converted on www.awesomepapers.org from LaTeX by an anonymous user.
Want to know more? Visit the Converter page.

Bistable Chimera Attractors on a Triangular Network of Oscillator Populations
July 26, 2025

Erik A. Martens Max Planck Research Group for Biological Physics and Evolutionary Dynamics, Max Planck Institute for Dynamics and Selforganization (MPIDS), Göttingen 37073, Germany erik.martens@ds.mpg.de
Abstract

We study a triangular network of three populations of coupled phase oscillators with identical frequencies. The populations interact nonlocally, in the sense that all oscillators are coupled to one another, but more weakly to those in neighboring populations than to those in their own population. This triangular network is the simplest discretization of a continuous ring of oscillators. Yet it displays an unexpectedly different behavior: in contrast to the lone stable chimera observed in continuous rings of oscillators, we find that this system exhibits two coexisting stable chimeras. Both chimeras are, as usual, born through a saddle node bifurcation. As the coupling becomes increasingly local in nature they lose stability through a Hopf bifurcation, giving rise to breathing chimeras, which in turn get destroyed through a homoclinic bifurcation. Remarkably, one of the chimeras reemerges by a reversal of this scenario as we further increase the locality of the coupling, until it is annihilated through another saddle node bifurcation.

chimera states, bistability, breathing chimeras, symmetry, topology, nonlocal coupling, oscillators
pacs:
05.45.Xt

I Introduction

While studying a continuum of identical oscillators on a ring with nonlocal coupling, Kuramoto et al. [13] discovered a remarkable state where the population of oscillators splits into two subpopulations, where one is synchronized and the other is desynchronized. This state was later dubbed a chimera, alluding to a monster in Greek mythology that consists of incongruous parts. Since then, several groups have explored the nonlinear dynamics of chimera states [13, 3, 2, 33, 10, 11, 1, 15, 28]. Their emergence on the ring was first investigated analytically by Abrams and Strogatz [3, 2]. They found that chimera states were born through a saddle node bifurcation, which appears to be the typical scenario for the emergence of chimeras on all network topologies investigated so far. Shima and Kuramoto [33] showed that chimeras also exist on 2D lattices with free boundaries in the shape of spiral waves: here, the center of the spiral, characterized by a topological defect, is replaced by a desynchronized core with finite positive radius [23].

A recent breakthrough in the field of coupled oscillator systems was made by Ott and Antonsen who showed that the dynamics of infinitely many oscillators could be reduced to a system of finite ordinary differential equations. The method has since then been generalized and put into a formal mathematical context by Pikovsky and Rosenblum [31] and Marvel et al. [26] independently, who show that the dynamics of each population may be reduced to a flow described by three variables plus constants of motion. Abrams et al. [1] employed the method to investigate a system of two interacting populations of phase oscillators; they were able to calculate the saddle node bifurcation for this system analytically and showed that the chimera states undergo a further change of stability to become breathing chimeras via a supercritical Hopf bifurcation.

Whereas many of these studies have been focusing on the nature of the chimera in idealized settings, others have been studying how chimeras would occur or behave in systems closer to real world situations: Omel’chenko et al. [28] show that a network of globally coupled oscillators, subjected to delayed feedback stimulation with spatially decaying profile, also induces chimera states; thereby they argue that chimera states indeed are a generic feature of coupled oscillator systems. Delay coupled systems are also investigated by Sethia et al. [32] who discover the clustering of chimera states. Makovetskiy et al. [20, 19] mention they found chimera-like states in systems of three level cellular automata related to lasing systems in combination with spiral waves. Laing investigates the robustness of chimera states against the introduction of nonidentical oscillators with heterogeneous frequencies in several systems [16]. For the case of two nonlocally coupled populations studied by Abrams et al. [1], Laing shows how chimeras may both be destabilized and stabilized as the strength of heterogeneity (i.e. the width of the frequency distribution) of the oscillators is varied[15]. Finally, we mention that in neuroscience, spatially localized “bumps” of neural activity are found in networks of spiking neurons - such states have been proposed as mechanisms for visual orientation tuning and working memory, and have been related to chimeras [17, 18].

Important questions still remain. In particular, which topologies allow for chimera states? In other terms, can we classify the network structures that allow for chimeras? Even on one-dimensional domains the situation is unclear. It has been shown that chimeras may exist on a ring [13] with a continuum of oscillators, but do chimeras also exist on a finite line segment or the infinite line, i.e. on a non-periodic domain? The observation that chimeras exist in the shape of spiral waves on 2D domains [33, 23] is a hint that this might be possible.

One may seek to answer this particular question by lowering the dimensionality of the problem and discretizing a one-dimensional domain using oscillators populations. As in Abrams et al. [1], one assumes that oscillators within a given population are coupled more strongly to each other than they are to those in neighboring populations, thereby defining a spatial structure on the network. Thus the simplest network that exhibits both chain-like and ring-like character consists of three populations, as shown in Fig. 1. By tuning a new structural parameter cc, a rotationally symmetric (ring-like or “triangular”) network (Fig. 1 (b)) is deformed into a less symmetric, chain-like structure (Fig. 1 (a)). The impact of changing the network topology in such ways is discussed in a companion article [21].

This paper is wholly devoted to the study of the purely triangular network as shown in Fig. 1 (b). While examining this simple case, we surprisingly found the coexistence of two stable chimera attractors. These attractors differ in the number of desynchronized populations (one or two). Their occurrence is unexpected because the ring with a continuum of oscillators [2, 3], sharing the same rotational symmetry, has only a single stable chimera state. The finding of multiple coexisting chimera attractors is novel and is relevant in various contexts. First, if we let the number of oscillator populations go to infinity, we would retrieve the system of a ring studied by Abrams et al.; however, this system is not known to yield any multistable chimeras, which forms an interesting discrepancy. Second, the finding of bistable chimeras is noteworthy in the context of the question raised just recently [27] about what the basins of attraction leading to chimera and globally synchronized states actually are. Third, in connection with the above mentioned ’bump states’, multistability of chimera states may be relevant in the study of mathematical neural models, where so called switching processes [4] and the competition of attractors in working memory [9] are of interest (see Discussion).

The organization of the article is as follows. We introduce the governing equations in Section II and summarize the derivation to obtain the reduced set of equations, as described in [1]. Next we derive the equations implied by special symmetries consistent with chimera states. These are analyzed in Section III, together with all their possible bifurcations. Section IV discusses our results in the context of numerical simulations. Finally, our findings are discussed in Section V. Additional results on stability of some simple symmetric states for arbitrary networks are outlined in the Appendix.

II Governing equations

The governing equations are given by

ddtθiσ\displaystyle\frac{d}{dt}\theta_{i}^{\sigma} =\displaystyle= ω+σ=13KσσNσj=1Nσsin(θjσθiσα),\displaystyle\omega+\sum_{\sigma^{\prime}=1}^{3}\frac{K_{\sigma\sigma^{\prime}}}{N_{\sigma^{\prime}}}\sum_{j=1}^{N_{\sigma^{\prime}}}\sin{(\theta_{j}^{\sigma^{\prime}}-\theta_{i}^{\sigma}-\alpha)}, (1)

where the phases of the oscillators are defined by θ\theta, ii denotes the individual oscillators belonging to the population with indices σ=1,2,3\sigma=1,2,3, each of which has NσN_{\sigma} oscillators, and where α\alpha is the phase lag parameter.

Refer to caption
Refer to caption
Figure 1: Networks of three populations of oscillators. The gray disks symbolize the populations of oscillators, populated by individual oscillators symbolized by black dots. Their bidirectional coupling is represented by black arrows. (a) Chain-like general case with structural tuning parameter cc. (b) Triangular network structure, corresponding to c=1c=1. Each population has a self-coupling of unit strength 11, and is coupled to the neighboring populations with strength 1A1-A.

The coupling kernel KσσK_{\sigma\sigma^{\prime}} describes the strength between populations σ\sigma and σ\sigma^{\prime}. The coupling strength is assumed to decay with increasing separation between the populations on the network. Within a population, the oscillators interact with strength Kσσ=1K_{\sigma\sigma^{\prime}}=1. Neighboring populations couple more weakly, with strength 1A1-A as displayed in Fig. 1 (b). We then have

Kσσ=(11A1A1A11A1A1A1).K_{\sigma\sigma^{\prime}}=\left(\begin{array}[]{lcr}1&1-A&1-A\\ 1-A&1&1-A\\ 1-A&1-A&1\end{array}\right). (2)

In the case of A=0A=0, we retrieve the case of a globally coupled network. Thus AA quantifies how ’far’ we are from global coupling. This network has the same rotational symmetry as a continuum of oscillators on a ring, studied by [13, 2, 3], and generalizes the problem with two populations discussed by Abrams et al. [1].

II.1 Reduction to low-dimensional system

We consider the limit of infinitely large populations, i.e. NσN_{\sigma}\rightarrow\infty. This allows us to reduce the problem to a finite set of equations, using the ansatz introduced by Ott and Antonsen [29], outlined below. In this limit, it is natural to describe the dynamics of the system in terms of the oscillator density distribution fσ(θ)f^{\sigma}(\theta), which evolves according to the continuity equation

fσt+θ(fσvσ)\displaystyle\frac{\partial f^{\sigma}}{\partial t}+\frac{\partial}{\partial\theta}(f^{\sigma}v^{\sigma}) =\displaystyle= 0.\displaystyle 0.\ (3)

The velocity of the oscillators is then given by

vσ\displaystyle v^{\sigma} =\displaystyle= ω+σ=13Kσσ02πsin(θθα)fσ(θ,t)𝑑θ.\displaystyle\omega+\sum_{\sigma^{\prime}=1}^{3}K_{\sigma\sigma^{\prime}}\int_{0}^{2\pi}\sin{(\theta^{\prime}-\theta-\alpha)}f^{\sigma^{\prime}}(\theta^{\prime},t)d\theta^{\prime}.

To keep the notation simple we denote θσ\theta_{\sigma} by θ\theta and θσ\theta_{\sigma^{\prime}} by θ\theta^{\prime}. It proves convenient to define a complex order parameter

zσ(t)\displaystyle z^{\sigma}(t) =\displaystyle= σ=13Kσσ02πeiθf(θ,t)𝑑θ,\displaystyle\sum_{\sigma^{\prime}=1}^{3}K_{\sigma\sigma^{\prime}}\int_{0}^{2\pi}e^{i\theta^{\prime}}f(\theta^{\prime},t)d\theta^{\prime},\ (5)

which defines a weighted average over all oscillators; we therefore refer to this as the global order parameter. Following Kuramoto’s footsteps [12, 14, 13], we rewrite the velocity in terms of the order parameter and find

vσ\displaystyle v^{\sigma} =\displaystyle= ω+Im[eiθeiασ=13Kσσ02πeiθf(θ,t)𝑑θ]\displaystyle\omega+\textrm{Im}\left[e^{-i\theta}e^{-i\alpha}\sum_{\sigma^{\prime}=1}^{3}K_{\sigma\sigma^{\prime}}\int_{0}^{2\pi}e^{i\theta^{\prime}}f(\theta^{\prime},t)d\theta^{\prime}\right] (6)
=\displaystyle= ω+12i[eiθeiαzσ(t)eiθeiαzσ(t)].\displaystyle\omega+\frac{1}{2i}\left[e^{-i\theta}e^{-i\alpha}z_{\sigma}(t)-e^{i\theta}e^{i\alpha}z_{\sigma}^{*}(t)\right].\

Following Ott and Antonsen [29], we now restrict attention to a special class of density functions in the form of a Poisson kernel

fσ(θ,t)\displaystyle f^{\sigma}(\theta,t) =\displaystyle= 12π{1+[k=1(aσ(t)eiθ)k+c.c.]},\displaystyle\frac{1}{2\pi}\left\{1+\left[\sum_{k=1}^{\infty}(a_{\sigma}(t)e^{i\theta})^{k}+c.c.\right]\right\},\ (7)

where c.c. is the complex conjugate of the expression under the sum. The implications of this special ansatz and its validity will be explained in more detail in Section V. Substitution of fσf_{\sigma} and vσv_{\sigma} into the continuity equation (3) yields an exact solution, so long as aσa_{\sigma} evolves according to

0\displaystyle 0 =\displaystyle= a˙σ+iωaσ+12aσ2zσeiα12zσeiα.\displaystyle\dot{a}_{\sigma}\ +i\omega a_{\sigma}+\frac{1}{2}a_{\sigma}^{2}\,z_{\sigma}\,e^{-i\alpha}-\frac{1}{2}\,z^{*}_{\sigma}\,e^{i\alpha}.\ (8)

It remains to express the order parameter in terms of this ansatz. We find

zσ(t)\displaystyle z_{\sigma}(t) =\displaystyle= σ=13Kσσaσ(t).\displaystyle\sum_{\sigma^{\prime}=1}^{3}K_{\sigma\sigma^{\prime}}a_{\sigma^{\prime}}^{*}(t).\ (9)

Finally, we express the amplitude aσa_{\sigma} in polar coordinates as

aσ=ρσeiϕσ.\displaystyle a_{\sigma}=\rho_{\sigma}e^{-i\phi_{\sigma}}.\

By division of (8) by eiϕσe^{-i\phi_{\sigma}}, we obtain

0\displaystyle 0 =\displaystyle= ρσ˙iϕσ˙ρσ+iωρσ\displaystyle\dot{\rho_{\sigma}}-i\dot{\phi_{\sigma}}\rho_{\sigma}+i\omega\rho_{\sigma}
+12ρσ2σ=13Kσσρσei(ϕσϕσα)\displaystyle+\frac{1}{2}\rho_{\sigma}^{2}\sum_{\sigma^{\prime}=1}^{3}K_{\sigma\sigma^{\prime}}\rho_{\sigma^{\prime}}e^{i(\phi_{\sigma^{\prime}}-\phi_{\sigma}-\alpha)}
12σ=13Kσσρσei(ϕσϕσα),\displaystyle-\frac{1}{2}\sum_{\sigma^{\prime}=1}^{3}K_{\sigma\sigma^{\prime}}\rho_{\sigma^{\prime}}e^{-i(\phi_{\sigma^{\prime}}-\phi_{\sigma}-\alpha)},\

and by separation of real and imaginary parts we find

ρ˙σ\displaystyle\dot{\rho}_{\sigma} =\displaystyle= 1ρσ22σ=13Kσσρσsin(ϕσϕσ+β),\displaystyle\frac{1-\rho_{\sigma}^{2}}{2}\sum_{\sigma^{\prime}=1}^{3}K_{\sigma\sigma^{\prime}}\rho_{\sigma^{\prime}}\sin{(\phi_{\sigma^{\prime}}-\phi_{\sigma}+\beta)},
ϕ˙σ\displaystyle\dot{\phi}_{\sigma} =\displaystyle= ω1+ρσ22ρσσ=13Kσσρσcos(ϕσϕσ+β),\displaystyle\omega-\frac{1+\rho_{\sigma}^{2}}{2\rho_{\sigma}}\sum_{\sigma^{\prime}=1}^{3}K_{\sigma\sigma^{\prime}}\rho_{\sigma^{\prime}}\cos{(\phi_{\sigma^{\prime}}-\phi_{\sigma}+\beta)},

where we introduce the definition

β=π/2α.\displaystyle\beta=\pi/2-\alpha.\ (11)

These equations describe the dynamics of our system in terms of the variables aσa_{\sigma}. Notice that in contrast to zσz_{\sigma}, they do not represent averages over all populations, and therefore, we refer to them as local order parameters. Thus any synchronized population σ\sigma of oscillators is characterized by ρσ=1\rho_{\sigma}=1. (These results are trivially generalized to the case of a network with arbitrarily many populations σ=1,2,,N\sigma=1,2,\ldots,N.)

In what follows, we will make particular assumptions about the symmetries of the solutions that we expect to find. Then we will analyze existence, stability and bifurcations of the states of interest.

II.2 Manifold of symmetric states (SDS and DSD)

To make progress on analyzing the chimera states, we will have to make certain symmetry assumptions. Perfectly synchronized populations have ρσ=1\rho_{\sigma}=1; desynchronized populations have ρσ<1\rho_{\sigma}<1, and consist of oscillators that drift relative to one another and to the synchronized populations. Let SS and DD denote synchronized and desynchronized populations, respectively. Then in a triangular network, we can distinguish only two chimera states, namely SDSSDS (sync-drift-sync) and DSDDSD (drift-sync-drift); all other permutations of SS and DD give equivalent states because of the rotational invariance inherent to the triangular network.

Another class of solutions can be described as SSSSSS, corresponding to a globally synchronized state where all oscillators are in sync, but with different synchronized mean phases ϕi\phi_{i}. We may distinguish three cases, i.e. ϕ1=ϕ2=ϕ3\phi_{1}=\phi_{2}=\phi_{3}, ϕ1=ϕ3ϕ2\phi_{1}=\phi_{3}\neq\phi_{2} and the state where all phase angles are different. These solutions are analyzed in Appendix A. These states are of less interest to us and we restrict ourselves to chimera states and their emergence in parameter space in the following.

The symmetry of our coupling kernel (2) in combination with ρ1=ρ3\rho_{1}=\rho_{3} implies that ϕ1=ϕ3\phi_{1}=\phi_{3} and hence populations 11 and 33 are phase-locked. The SDSSDS state is then defined via

ρ1\displaystyle\rho_{1} =\displaystyle= ρ3=1 and ρρ2<1,\displaystyle\rho_{3}=1\textrm{ and }\rho\equiv\rho_{2}<1,
ϕ1\displaystyle\phi_{1} =\displaystyle= ϕ3,\displaystyle\phi_{3},\

whereas the DSDDSD state is given by

ρ\displaystyle\rho \displaystyle\equiv ρ1=ρ3<1 and ρ2=1.\displaystyle\rho_{1}=\rho_{3}<1\textrm{ and }\rho_{2}=1.
ϕ1\displaystyle\phi_{1} =\displaystyle= ϕ3,\displaystyle\phi_{3},\

We define the phase difference of the angular order parameter between the synchronized and desynchronized states by

ψ=ϕ1ϕ2=ϕ3ϕ2.\displaystyle\psi=\phi_{1}-\phi_{2}=\phi_{3}-\phi_{2}.\ (12)

Applying these symmetry assumptions to (II.1) and substituting the coupling kernel defined in (2), we obtain the equations describing the SDSSDS states

ρ˙\displaystyle\dot{\rho} =\displaystyle= 1ρ22[2(1A)sin(ψ+β)+ρsinβ],\displaystyle\frac{1-\rho^{2}}{2}\left[2(1-A)\sin{(\psi+\beta)}+\rho\sin{\beta}\right],
ψ˙\displaystyle\dot{\psi} =\displaystyle= (2A)cosβ(1A)ρcos(ψ+β)\displaystyle-(2-A)\cos{\beta}-(1-A)\rho\cos{(-\psi+\beta)} (13)
+\displaystyle+ 1+ρ22ρ[2(1A)cos(ψ+β)+ρcosβ],\displaystyle\frac{1+\rho^{2}}{2\rho}\left[2(1-A)\cos{(\psi+\beta)}+\rho\cos{\beta}\right],\

and the DSDDSD states

ρ˙\displaystyle\dot{\rho} =\displaystyle= 1ρ22[(2A)ρsinβ+(1A)sin(ψ+β)],\displaystyle\frac{1-\rho^{2}}{2}\left[(2-A)\rho\sin{\beta}+(1-A)\sin{(-\psi+\beta)}\right],
ψ˙\displaystyle\dot{\psi} =\displaystyle= 1+ρ22ρ[(2A)ρcosβ+(1A)cos(ψ+β)]\displaystyle-\frac{1+\rho^{2}}{2\rho}\left[(2-A)\rho\cos{\beta}+(1-A)\cos{(-\psi+\beta)}\right] (14)
+\displaystyle+ 2(1A)ρcos(ψ+β)+cosβ.\displaystyle 2(1-A)\rho\cos{(\psi+\beta)}+\cos{\beta}.\

Note that these equations hold only if we restrict attention to symmetry-preserving perturbations. The fixed points of (II.2,II.2) correspond to phase-locked solutions of the original system (at the macroscopic level of the local order parameters).

By reduction of the full system of oscillators (1) to a low dimensional system for the local order parameters, we have cast our problem into a two dimensional system represented by Eqs. (II.2,II.2). This enables us to study the problem in the phase plane.

III Analysis

III.1 Phase portraits

Unfortunately, the equations for the SDSSDS and DSDDSD states cannot be solved in closed form. Before we get deeper into the matter of analyzing these states, let us get a quick intuition of their behavior by inspecting the phase planes of their corresponding equations, which will guide us in the subsequent analysis. Their phase portraits shown in Fig. 2 and Fig. 4 represent a sweep in parameter space with increasing values of AA while keeping the value of β=0.1\beta=0.1 constant.

Refer to caption
Refer to caption
Refer to caption
Refer to caption
Figure 2: Phase portraits for the SDSSDS chimera, with increasing values of AA at constant β\beta. Real (xx) and imaginary (yy) components of a¯(σ)\overline{a}_{(\sigma)} are shown. The unit circle is displayed in gray. Stable fixed points are shown as solid (red) and half-filled (green) circles, respectively. Limit cycles are emphasized in blue color. Stable and unstable manifolds are shown as red solid and green dashed trajectories, respectively. The point in (ρ,ψ)=(1,0)(\rho,\psi)=(1,0) is a nodal sink. The position of the nodal source depends on β\beta and moves in clockwise direction with growing values of β\beta.

Consider first the SDSSDS symmetry. For small values of AA (close to global coupling), we only observe a nodal sink and source on the unit circle; these points correspond to the in-phase SSSSSS solutions (Fig. 2(a)). Increasing AA further, a saddle-node pair is born very close to the unit circle. For larger AA, the node moves closer to the origin, implying that the order of the desynchronized population decreases and the chimera state becomes more pronounced (if we instead increase the values of β\beta, this critical point starts to move clockwise while getting closer to the origin). The node has become a stable spiral (Fig. 2(b)), and at a critical value of AA, it loses stability through a Hopf bifurcation and a limit cycle is born (Fig. 2(c)). The amplitude of the order parameter ρ\rho of the drifting population starts to oscillate, and is therefore called a breathing chimera. As we raise the value of AA more, the limit cycle gains in amplitude until it collides with the saddle: the limit cycle is destroyed in a homoclinic bifurcation (Fig. 2(d)). The resulting bifurcation diagram is shown in Fig. 3 (a). The saddle-node, Hopf and homoclinic bifurcation curves all intersect in a Takens-Bogdanov point with codimension two.

Refer to caption
Refer to caption
Figure 3: Bifurcation diagram for the SDSSDS (a) and DSDDSD (b) chimera. The displayed curves are: the saddle-node curve (solid red), the Hopf curve (dashed blue), and the homoclinic curve (dotted black). Dots mark the bifurcation points obtained by inspection of the phase plane. The homoclinc curve (dotted) is an interpolation based on these points, whereas all the solid curves were obtained analytically.

For the DSDDSD symmetry, we observe a scenario that is qualitatively similar to the previous case. However, surprisingly, the whole scenario appears a second time in the upper part of the parameter plane, but now in reversed order, as shown in Fig. 3 (b). We again increase AA while keeping β\beta constant, as shown in Fig. 4. For all parameter values, we find two synchronized SSSSSS solutions on the unit circle: one is a nodal sink in (ρ,ψ)=(1,0)(\rho,\psi)=(1,0), but what in the SDSSDS case before was a nodal source, is now a saddle. Also notice that a new fixed point has appeared in the left half of the unit circle in the form of a spiral source (Fig. 4(a)). This is the second, currently unstable, DSDDSD chimera seen in the upper half of the bifurcation diagram. As AA increases and the coupling becomes more local, a saddle-node pair is born in the right half of the unit circle (Fig. 4(b)); again, its node then becomes a spiral and loses stability as the coupling strength becomes more local, and the resulting limit cycle (Fig. 4(c-d)) gets ultimately destroyed in a homoclinic bifurcation (Fig. 4(e)). Whereas one chimera has been rendered unstable, we observe that, above a critical AA, a stable limit cycle has formed around the spiral source on the left half of the circle: the twin of the DSDDSD chimera in its breathing mode has emerged (Fig. 4(e)). From here, the bifurcations happen in reversed order, the source becomes first a spiral node (Fig. 4(f)), i.e. a stable chimera, which is eventually annihilated in a saddle-node bifurcation. The resulting bifurcation diagram is seen in Fig. 3 (b).

Refer to caption
Refer to caption
Refer to caption
Refer to caption
Refer to caption
Refer to caption
Figure 4: Phase portraits for the DSDDSD chimeras, with increasing values of AA at constant β\beta. Real (xx) and imaginary (yy) components of a¯(σ)\overline{a}_{(\sigma)} are shown. The unit circle is displayed in gray (the unstable green dashed manifold of the saddle coincides with it). Stable and unstable fixed points are shown as solid (red) and half-filled (green) circles, respectively. Limit cycles are emphasized in blue color. Stable and unstable manifolds are shown as red solid and green dashed trajectories, respectively. The point in (ρ,ψ)=(1,0)(\rho,\psi)=(1,0) is a nodal sink. The position of the saddle depends on β\beta and moves in clockwise direction with growing values of β\beta.

III.2 Calculation of bifurcation curves

In order to calculate the saddle-node and Hopf curves of the SDSSDS solutions, we must linearize (II.2) around the appropriate fixed point. This task amounts to solving the fixed point equations implied by Eqs. (II.2) and (II.2) simultaneously with the saddle node condition,

det(J)\displaystyle\det{(J)} =\displaystyle= 0,\displaystyle 0,\

or with the Hopf condition,

tr(J)\displaystyle\textrm{tr}{(J)} =\displaystyle= 0anddet(J)>0,\displaystyle 0\quad\textrm{and}\quad\det{(J)}>0,\

where JJ denotes the Jacobian of (II.2) or (II.2), respectively. For the SDSSDS symmetry, we have

J11\displaystyle J_{11} =\displaystyle= 12[(13ρ2)sinβ4(1A)ρsin(β+ψ)],\displaystyle\frac{1}{2}\left[(1-3\rho^{2})\sin{\beta}-4(1-A)\rho\sin{(\beta+\psi)}\right],
J12\displaystyle J_{12} =\displaystyle= (1A)(1ρ2)cos(β+ψ),\displaystyle(1-A)(1-\rho^{2})\cos{(\beta+\psi)},
J21\displaystyle J_{21} =\displaystyle= 1ρ2[cosβ(ρ3(1A)cosψ)\displaystyle\frac{1}{\rho^{2}}\left[\cos{\beta}(\rho^{3}-(1-A)\cos{\psi})\right.
sinβsinψ(2ρ21)(1A)],\displaystyle-\left.\sin{\beta}\sin{\psi}(2\rho^{2}-1)(1-A)\right],
J22\displaystyle J_{22} =\displaystyle= A1ρ[(1+2ρ2)cosψsinβ+cosβsinψ].\displaystyle\frac{A-1}{\rho}\left[(1+2\rho^{2})\cos{\psi}\sin{\beta}+\cos{\beta}\sin{\psi}\right].\

The fixed point condition implied by (II.2) yields the nontrivial solution

ρ\displaystyle\rho =\displaystyle= 2(A1)cscβsin(β+ψ).\displaystyle 2(A-1)\csc{\beta}\sin{(\beta+\psi)}.\ (15)

Substitution of this expression into the fixed point condition for ψ\psi results in

0\displaystyle 0 =\displaystyle= (A2)cosβ+12csc(β+ψ)sinψ\displaystyle(A-2)\cos{\beta}+\frac{1}{2}\csc{(\beta+\psi)}\sin{\psi} (16)
+\displaystyle+ (1A)2cscβ[sin2β+2(cotβsin2ψ+sin2ψ)].\displaystyle(1-A)^{2}\csc{\beta}\left[\sin{2\beta}+2(\cot{\beta}\sin^{2}{\psi}+\sin{2\psi})\right].\

Unfortunately, this equation cannot be solved in closed form for ψ\psi, which in turn would allow us to express AA in terms of β\beta. We settle therefore for a series approach in AA and ψ\psi, as follows:

ψ\displaystyle\psi =\displaystyle= k=0Npkβk+O(βN+1)\displaystyle\sum_{k=0}^{N}\,p_{k}\,\beta^{k}+{O}(\beta^{N+1}) (17)
A\displaystyle A =\displaystyle= k=0Nakβk+O(βN+1).\displaystyle\sum_{k=0}^{N}\,a_{k}\,\beta^{k}+{O}(\beta^{N+1}).\ (18)

We substitute these two expressions into fixed point equation (16) and the saddle node condition, and solve the resulting equations for each power of β\beta. This leads to the following expression for the saddle node curve:

ASN(β)\displaystyle A_{SN}(\beta) =\displaystyle= 92β634β2+1954β3235516β4\displaystyle\frac{9}{2}\,\beta-\frac{63}{4}\,\beta^{2}+\frac{195}{4}\,\beta^{3}-\frac{2355}{16}\,\beta^{4}
+3528380β5210247160β6+872617224β7\displaystyle+\frac{35283}{80}\,\beta^{5}-\frac{210247}{160}\,\beta^{6}+\frac{872617}{224}\,\beta^{7}
2949379256β8+274411626180640β9+O(β10).\displaystyle-\frac{2949379}{256}\,\beta^{8}+\frac{2744116261}{80640}\,\beta^{9}+{O}(\beta^{10}).

Using the Hopf condition and proceeding in the same way we also find the Hopf curve, approximated by

AH(β)\displaystyle A_{H}(\beta) =\displaystyle= 0.447153+1.34639β2+8.34371β4+O(β6).\displaystyle 0.447153+1.34639\,\beta^{2}+8.34371\,\beta^{4}+{O}(\beta^{6}).

The Takens-Bogdanov point is determined by numerically solving the saddle-node, Hopf and fixed point conditions simultaneously. It is located at

(β,A)SDS=(0.1974,0.5092).\displaystyle(\beta,A)_{SDS}=(0.1974,0.5092).\ (21)

The bifurcation curves for the DSDDSD chimeras are obtained in an analogous procedure. Solving the equations expanded in series is now a bit trickier due to the coexistence of two branches. For brevity, we shall only summarize our findings. The two saddle node curves are approximated by

ASN,1(β)\displaystyle A_{SN,1}(\beta) =\displaystyle= 9β72β2+10592β33855β4\displaystyle 9\,\beta-72\,\beta^{2}+\frac{1059}{2}\,\beta^{3}-3855\,\beta^{4} (22)
+113094340β510392765β6\displaystyle+\frac{1130943}{40}\,\beta^{5}-\frac{1039276}{5}\,\beta^{6}
+854234093560β7783117837β8\displaystyle+\frac{854234093}{560}\beta^{7}-\frac{78311783}{7}\,\beta^{8}
+330978868116140320β9+O(β10),\displaystyle+\frac{3309788681161}{40320}\,\beta^{9}+{O}(\beta^{10}),
ASN,2(β)\displaystyle A_{SN,2}(\beta) =\displaystyle= 1β+316β3β46421120β5+O(β6).\displaystyle 1-\beta+\frac{31}{6}\,\beta^{3}-\,\beta^{4}-\frac{6421}{120}\,\beta^{5}+{O}(\beta^{6}).

Unfortunately, the series approach for the lower second branch converges extremely slowly and doesn’t match all the way to the Takens-Bogdanov point (even going to this high order doesn’t help!). However the effort to find the series coefficients was not all in vain: a Padé approximant based on the above power series at order three does an excellent job in matching the data points retrieved from the examination of the phase portraits. We have

ASN,2(β)\displaystyle A_{SN,2}(\beta) \displaystyle\approx 9β+1521018301403β2+2874468276028060β31+2580226301403β+901238336028060β2+2002078284521045β3.\displaystyle\frac{9\,\beta+\frac{1521018}{301403}\,\beta^{2}+\frac{287446827}{6028060}\,\beta^{3}}{1+\frac{2580226}{301403}\,\beta+\frac{90123833}{6028060}\,\beta^{2}+\frac{200207828}{4521045}\,\beta^{3}}.

Finally, the Hopf curves are approximated by

AH,1(β)\displaystyle A_{H,1}(\beta) =\displaystyle= 0.593737+1.14491β2\displaystyle 0.593737+1.14491\,\beta^{2} (25)
+4.55308β4+O(β6),\displaystyle+4.55308\,\beta^{4}+{O}(\beta^{6}),
AH,2(β)\displaystyle A_{H,2}(\beta) =\displaystyle= 0.8854081.15074β2\displaystyle 0.885408-1.15074\,\beta^{2} (26)
3.96289β4+O(β6).\displaystyle-3.96289\,\beta^{4}+{O}(\beta^{6}).\

The Takens-Bogdanov points are located at

(β,A)DSD,1\displaystyle(\beta,A)_{DSD,1} =\displaystyle= (0.2132,0.6615)\displaystyle(0.2132,0.6615) (27)
(β,A)DSD,2\displaystyle(\beta,A)_{DSD,2} =\displaystyle= (0.1903,0.8359).\displaystyle(0.1903,0.8359).\ (28)

For all symmetries, we found an excellent agreement of our perturbative results with the bifurcation points obtained from inspection of the phase portraits. We did not derive an analytical expression for the homoclinic bifurcation curves; the curves shown in the Figs. 3 (a) and 3 (b) are based on data points obtained from inspecting the phase portraits while varying the parameters.

IV Numerical simulations

We have obtained analytical results describing the dynamics of our triangular network of oscillator populations, using two different reductions: firstly, we reduced the governing Eqs. (1) using the Ott Antonsen method. And secondly, we have assumed that the system attains certain symmetry states that allow for chimera states; these symmetries need not be transversely stable to perturbations off the symmetry manifold. Thus, the equations obtained from these reductions may not necessarily account for the complete dynamics of the governing equations, and we checked if our analytical results agree with numerical simulations of the governing equations. We did this with a finite, but what may be considered a sufficiently large oscillator population (N(σ)=40N^{(\sigma)=40}).

We used two different methods to generate initial conditions that would lead to the appearance of chimera states: Firstly, for the phases of the desynchronized populations we used an initial condition in the shape of a bump, specifically, a Gaussian distribution in the shape of ϕiσ,dexp(γ(i/Nσ1/2)2)\phi_{i}^{\sigma,d}\sim\exp{(-\gamma(i/N_{\sigma}-1/2)^{2})} with an appropriately chosen decay rate γ\gamma. The second method was chosen with the intention to place the system right on the Ott-Antonsen (OA) manifold. This was accomplished by generating a phase distribution that is consistent with the Poisson kernel (7). For both methods, the phases for the synchronized populations are given by ϕiσ,s=0\phi_{i}^{\sigma,s}=0. To achieve this, we solved the SDSSDS and DSDDSD Eqs. (II.2) and (II.2), respectively, for fixed points (ρ,ψ)(\rho,\psi). This in turn, enables us to compute the Poisson distribution fσ(θ,t)f^{\sigma}(\theta,t) defined by (7), by using the definition of the order parameter, aσ=ρσeiϕσa_{\sigma}=\rho_{\sigma}e^{-i\phi_{\sigma}}. Because this function defines the probability with which oscillators populate a certain phase, we may use its inverse cumulative distribution function to construct from it a set of phases that is consistent with the OA-manifold. The system should remain close to the OA-manifold, because of its invariance. Unless mentioned otherwise, we used this latter method.

We first confirmed that the unreduced system would exhibit all types of chimeras predicted by our analysis, and that they would correspond to stable states, for various points in parameter space. These states were observed with both Nσ=20N_{\sigma}=20 and Nσ=200N_{\sigma}=200 oscillators per population. The observed behavior is that the system first goes through a tiny transient and reach an attracting state which was confirmed to be stable even for long computation times. The transient may be explained by the fact that the system, due to its finite size, can only be approximately on to the OA-manifold (or a member of the OA-family, as explained in the Discussion).

Next we checked if the critical parameter values of saddle node, Hopf and homoclinic bifurcation in the full system would be in accord with the critical values obtained from our theory, see Figs. 5 (a) and (b). In order to do so, we held the value of β\beta fixed and continued a solution through 2020 increasing values of AA. To initiate the continuation we used an initial condition consistent with the OA-manifold.

In Fig. 5, we show fixed points and oscillation amplitudes of the order parameter, obtained analytically from (II.2) and (II.2), as dashed curves. Our simulation results are superposed, and were obtained as follows: To remove transient effects from our analysis, we only considered the last 2/52/5 of the computed time series. Instead of only detecting the global maximum and minimum of the series, we detected local maxima, shown as light gray dots, and minima, shown as dark gray dots. We chose to do this, because it would allows us to see the traces of new appearing periods, that potentially may occur in the unreduced, highly dimensional system. However, finite size fluctuations also cause small oscillations; this is the reason we see some small amount of blurring in the data (For similar reasons, maxima and minima on the stable branch do not coincide.)

Refer to caption
Refer to caption
Figure 5: Bifurcation diagram obtained from numerical simulation, for states SDSSDS, shown in (a), and DSDDSD in (b). Light dots and dark circles correspond to local maxima and minima, respectively, that are detected by the algorithm. The dashed curve represents the analytical result for the continuum case (NN\rightarrow\infty). The computations were performed with Nσ=40N_{\sigma}=40 oscillators per population for a simulation time of T=100T=100. The question mark indicates that we could not conclusively confirm the existence of the second DSD state. (The kink in the lower left dashed branch is an artefact from the limit cycle reaching into the left hand side quadrants, as ρ\rho is not measured relative to the limit cycle center but to the origin.)

Despite these undesired effects, we can clearly demonstrate the onset of the Hopf and the homoclinic bifurcations in the simulation. Consider first the SDSSDS chimera, shown in Fig. 5 (a). We expect that finite size effects affect the locations of all the bifurcation points. While we would have to compute at a higher resolution for the continuation to actually see this happen for the Hopf bifurcation, it is more apparent for the homoclinic bifurcation: the limit cycle oscillation brings the system periodically close to the saddle point (as seen in the phase plane); therefore, the larger finite size fluctuations are, the more likely the system is to be kicked off the limit cycle and is instead attracted to the nearby all-in-phase SSSSSS state on the unit circle. This is seen in Fig. 5 (a), where the limit cycle oscillation disappears at a smaller value of AA than it does for the analytic result (sold curve indicates SSSSSS state). The same behavior is observed for the DSDDSD chimera in Fig. 5 (b), but much more pronounced: the system snaps to the in-phase state (solid curve) much earlier than expected for a continuous system of oscillators. Similarly, we can continue the solutions in reverse direction, and check that the chimera states are annihilated via a saddle node bifurcation.

Whereas we easily managed to show DSDDSD states in the lower half of the parameter plane, our attempts to place the system on the second DSDDSD attractor had little success: for all trials the system would eventually reach the in-phase state. Increasing the number of populations didn’t seem to remedy the matter, as one might be led to think from our previous experience; rather, the system would stay obediently in the DSDDSD state for a little while and then suddenly jump into the in-phase state. However, it doesn’t seem to be entirely clear whether the second attractor of the DSDDSD state is inherently unstable, or if it is just very hard to stay on that manifold due to a combination of factors, given by finite size effects and the location and shape of the attractor (see phase plane in Fig. 4). Specifically, in the case of a breathing state, the limit cycle is always very close to the invariant field defined by the in-phase state on the unit circle. In conclusion, it is likely that this second DSDDSD attractor is unstable to symmetry breaking perturbations.

V Discussion

We have investigated a triangular network of populations of sinusoidally, nonlocally coupled oscillators with identical frequencies. Our analysis approximates the practical case of large populations by considering the limit of infinitely many oscillators per population, for which we may use the Ott-Antonsen method to reduce the governing to a finite set of equations [29]. By assuming symmetries compliant with chimera states, we further reduced these equations and studied their emergence in a phase plane analysis. Saddle node and Hopf bifurcation curves are determined using perturbation techniques, and homoclinic bifurcation curves by observation of phase portraits. The resulting stability diagram is a variation of the one found by Abrams et al. [1], but is governed by the emergence of two different types of chimeras: one with a desynchronized ’zone’ of narrow width (SDSSDS), and the other with a desynchronized zone of twice the width (DSDDSD). DSDDSD chimeras exist in two regions of parameter space, and differ in their mean (difference) phase angle ψ\psi; for the DSDDSD state near global coupling (small AA), ψ\psi remains close to zero, and for the other near local coupling (large AA) closer to π\pi. Numerical experiments demonstrate that chimeras persist for small oscillator populations, and that the SDSSDS and the DSDDSD chimera near global coupling are truly stable attractors in the unreduced system described by Eqs. (1). Finally, we observe for the first time the feature of bistable chimeras.

V.1 Reduced system and numerical experiments

Ott and Antonsen [30] showed that the long time dynamics of a large class of phase oscillator systems (including this one) is attracted to the manifold defined by (7) if the natural frequencies are drawn from peaked distributions. This result has been confirmed in several studies [5, 7, 15, 22]. The situation is somewhat different if we consider the case of identical frequencies. Pikovsky and Rosenblum [31], and Marvel et al. [24] have shown independently that the dynamics of each population can be reduced to a flow described by three variables plus constants of motion. The unifying picture is that for identical frequencies one has a whole one-parameter family of invariant manifolds, including the OA manifold. These manifolds are neutrally stable with regards to perturbations in transverse direction to themselves (such perturbation would result in placing the dynamics into a neighboring manifold). Once the frequencies are spread this family collapses and the OA manifold is left.

Because this degeneracy in the case of identical frequencies might lead to unanticipated results, and because the symmetry manifolds defined by SDSSDS and DSDDSD don’t need to be stable to perturbations in transverse directions, we checked that the chimera states are true attractors in the unreduced system, Eqs. (1), by numerical simulation. The sequence of bifurcations with saddle node, Hopf and homoclinic bifurcations is indeed reproduced in the unreduced system and the bifurcations appear close to the predicted critical values, even though the number of oscillators per population is small. The homoclinic bifurcation makes an exception in the sense that it seems to happen already for small AA; we have argued that this likely is an artefact for large limit cycles, where the trajectory gets kicked off the orbit because of finite size fluctuations in combination with the appearance of a nearby saddle point. All chimera states, with exception of the DSDDSD state near local coupling, have proved to be true attractors in the unreduced system.

V.2 Relation to heterogeneous frequency distributions

A recent study by Laing [15] generalizes the problem of a network with two oscillator populations investigated by [1] to the case of heterogeneous frequencies. Laing showed that for this and various other network topologies, the chimera state is robust – within limits – to heterogeneity in the intrinsic frequencies of the oscillators. In particular, he finds that the chimera state remains stable for populations with nearly identical oscillators, that is, with a narrow width of the distribution. The stability results obtained from our analysis should therefore be the same as the one obtained for the dynamics of oscillators with almost identical frequencies.

V.3 Aperiodicity and chaos

In our setting of the problem with identical oscillators, we were able to find solutions that do not lie on the OA-manifold (to find them, use initial conditions lying off the OA-manifold, e.g. Gaussian ’bumps’ or add noise to the OA-compliant initial conditions). Such states may be related to the quasi-periodic states shown to exist by Pikovsky et al. [31], but perhaps also to chaotic states. We haven’t gone further into this question, but have noticed irregular behavior in some of our simulations that look aperiodic, and we think it might be possible to find something like a chaotic chimera state.

V.4 All-in-phase states

In the case where all populations are fully synchronized, i.e. ρσ=1\rho_{\sigma}=1, the dynamics of the system effectively becomes the one governed by three coupled oscillators. The stability of such states is analyzed in the Appendix. This case has already been studied in the context of three coupled sinusoidal limit cycle oscillators by Mendelowitz et al. [25] and for relaxation limit cycle oscillators by Bridge et al. [6], who found that two rotating waves, clockwise and counter-clockwise rotation, are possible.

V.5 Relation to ring of oscillators and network with two oscillator populations

We mention the connection of our study to the work done by Abrams et al. [1], who studied a similar system, but with only two oscillator populations. With its two in-phase-locked populations, it is unclear whether the SDSSDS or the DSDDSD state compares best to their – using our terminology – SDSD chimera. Certainly, it can be said they both act like a two population system in disguise. In this sense, the two population system is a degenerate case of our triangular network.

The same authors also studied a continuum of oscillators on a ring [2, 3]. In our picture with oscillator populations, this system is approximated by an infinite set of populations arranged on a ring. For the continuous ring, only a single chimera is known. Both SDSSDS and DSDDSD chimeras effectively act like a system made of two populations, but differ in their ’width’ of drifting populations. In the case of the ring, this width is slightly larger than the synchronized region; from this point of view, the DSDDSD chimera seems to match their state best.

One may ask how the behavior of discrete ring-like systems is affected as we increase the number of populations. One could imagine that more and more chimera states get added as we increase the number of populations, constituting competing multistable attractors. But then, what happens as we take this continuum limit? Do they disappear, collapse such that only one of them, specifically, the one discussed in [3], wins the competition, and dominates over all others and remains stable? In the case of the triangular network, the SDSSDS and the first DSDDSD chimera are truly stable attractors in the unreduced system, and so the situation stays inconclusive. The creation and annihilation of chimera states with varying numbers of oscillator populations is also of importance if we want to characterize the basins of attractions in complex oscillator networks, an issue pointed out recently by Motter [27].

V.6 Multistability and neural networks

Our numerical experiments show that the existence of chimeras not only exist in the thermodynamic limit of NN\rightarrow\infty, but also for relatively low numbers of oscillators per population of N20N\sim 20. In view of the robustness of chimera states towards heterogeneity [15], one may therefore expect that chimeras also may emerge in biological settings such as neuronal networks. In this case, the coexistence of chimera attractors would have implications for the study of neuronal networks beyond their similarity to spatially localized ”bumps” of neuronal activity mentioned earlier. For instance, one may envision that multiple coexisting chimera attractors encode states of associative memory [25, 9] or play a role in decision making, where different initial conditions lead to different states. The characterization of basins of attractions for chimeras therefore represents an interesting issue, and represents a topic in the field that is not well explored [27]. Alternatively to stable-state dynamics, some neuronal models, based on the mechanism of switching among unstable saddles, have drawn attention lately [4]. In these models, heteroclinic connections generate characteristic patterns of neural activity and thus represent information. One may conceive oscillator networks with the possibility of switching between chimera states while varying system parameters. Provided that such chimera sequences, encoding different spiking activity patterns, are realizable, they may constitute an interesting field of research. We note that in this context oscillator populations, or rather neurons, represent redundancy in computing networks, which significantly enhances their reliability, as shown already by von Neumann [34]. In a recent experimental study Feinerman et al. [8] build neuronal computational devices and demonstrate how reliability is significantly improved by employing neuronal ensembles for each logic unit, similar to the oscillator populations.

VI Acknowledgments

This research was supported in part by NSF Grant No. DMS-0412757 and CCF 0835706. I would like to thank Steve Strogatz for advice, mentoring and helpful guidance throughout this project, and Fabio Schittler-Neves for useful conversations.

Appendix A Stability of Symmetric States

A.1 Existence of Symmetric States SSSSSS

The states we consider here are fully synchronized, i.e. ρσ=1\rho_{\sigma}=1 for all σ=1,2,3\sigma=1,2,3. We consider the symmetries where ϕ1=ϕ3\phi_{1}=\phi_{3}. In order to derive a fixed point condition for these states, it is therefore sufficient to consider the equations that are specialized to the symmetries SDSSDS and DSDDSD given by Equations (II.2) and (II.2). Applying the fixed point equation to either of them yields the condition

0\displaystyle 0 =\displaystyle= (1A)(cosβcosβcosψ+3sinψsinβ),\displaystyle(1-A)(\cos{\beta}-\cos{\beta}\cos{\psi}+3\sin{\psi}\sin{\beta}),\ (29)

which is reduced to these cases:

A\displaystyle A =\displaystyle= 1,\displaystyle 1, (30)
β\displaystyle\beta =\displaystyle= 0 with cosψ=1,sinψ0,\displaystyle 0\quad\textrm{ with }\cos{\psi}=1,\;\sin{\psi}\neq 0, (31)
ψ\displaystyle\psi =\displaystyle= 0,\displaystyle 0, (32)
cosβ\displaystyle\cos{\beta} =\displaystyle= 3sinβsinψcosψ1 with cosψ1.\displaystyle\frac{3\sin{\beta}\sin{\psi}}{\cos{\psi}-1}\textrm{ with }\cos{\psi}\neq 1.\ (33)

The first two conditions are the degenerate cases representing the AA and β\beta axis. The third condition corresponds to the fixed point (ρ,ψ)=(1,0)(\rho,\psi)=(1,0), and the last to the position of the node that is constrained to move on the unit circle in the phase portraits, see Figs. 2 and 4.

A.2 Stability of SSSSSS States

The computation of stability of these points is accomplished by computation of the Jacobian of the six dimensional system (II.1), with the coordinate system (ρ1,ρ2,ρ3,ϕ1,ϕ2,ϕ3)(\rho_{1},\rho_{2},\rho_{3},\phi_{1},\phi_{2},\phi_{3}), using a computer algebra system. Applying our symmetry assumptions, we find the eigenvalues

λk=(0sinβ+2(A1)sin(β+ψ)(A2)sinβ+(A1)sin(βψ)(A2)sinβ+(A1)sin(βψ)(A1)(2sinβ+sin(βψ))(A1)(3cosψsinβ+cosβsinψ)).\displaystyle\lambda_{k}=\left(\begin{array}[]{c}0\\ -\sin{\beta}+2(A-1)\sin{(\beta+\psi)}\\ (A-2)\sin{\beta}+(A-1)\sin{(\beta-\psi)}\\ (A-2)\sin{\beta}+(A-1)\sin{(\beta-\psi)}\\ (A-1)(2\sin{\beta}+\sin{(\beta-\psi)})\\ (A-1)(3\cos{\psi}\sin{\beta}+\cos{\beta}\sin{\psi})\end{array}\right). (40)

The first eigenvalue is a manifestation of the rotational invariance of the system (the system only depends on phase differences and is effectively five dimensional).

We first compute the stability of the trivially symmetric state defined by ψ=0\psi=0. This degenerate state has the eigenvalues

0\displaystyle 0 =\displaystyle= λ(λ(2A3)sinβ)3(λ(A1)sinβ)2\displaystyle\lambda(\lambda-(2A-3)\sin{\beta})^{3}(\lambda-(A-1)\sin{\beta})^{2}\ (41)

If we only consider parameter values A(0,1)A\in(0,1), a sufficient condition for linear stability is given by β(0,π)\beta\in(0,\pi).

The less trivial state with ψ0\psi\neq 0 must be considered in combination with the fixed point condition (33). The signs of these eigenvalues were not obtained analytically but rather by graphing the maximal eigenvalue in the (β,A)(\beta,A)-plane. It turns out that this state has at least one positive eigenvalue except for the degenerate case where β=0,π\beta=0,\pi, and is thus always a saddle. This result is consistent with the behavior observed by inspection of the phase portraits of Eqs. (II.2) and (II.2) (in the case of SDSSDS symmetry, the nontrivially symmetric state changes stability at β=π\beta=\pi, but this holds only within the SDSSDS-symmetry manifold, and has nothing to do with stability in the six (or five, for that matter) dimensional space.).

It is possible to obtain a similar stability result for the general case of a network with arbitrarily many populations NN, for the trivially symmetric point satisfying ρσ=1\rho_{\sigma}=1 and ϕ1=ϕ2==ϕN\phi_{1}=\phi_{2}=...=\phi_{N}, using Gershgorin’s circle theorem. In this general setting, the point also becomes linearly unstable as β>π\beta>\pi (provided that all row sums of the coupling matrix are strictly positive). Unfortunately, no similar result was obtained for the remaining N2N-2 fixed points that may occur on the unit-sphere related to the general problem. The calculation is tedious and will not be represented here.

References

  • [1] D. M. Abrams, R. E. Mirollo, S. H. Strogatz, and D. A. Wiley. Solvable model for chimera states of coupled oscillators. Phys. Rev. Lett., 101:084103, 2008.
  • [2] D. M. Abrams and S. H. Strogatz. Chimera states for coupled oscillators. Phys. Rev. Lett., 93:174102–174102, 2004.
  • [3] D. M. Abrams and S. H. Strogatz. Chimera states in a ring of nonlocally coupled oscillators. Intern. Journ. Bif. Chaos, 16(1):21–37, 2004.
  • [4] P. Ashwin and M. Timme. When instability makes sense. Nature, 436:36–37, 2005.
  • [5] E. Barreto, B. Hunt, E. Ott, and P. So. Synchronization in networks of networks: The onset of coherent collective behavior in systems of interacting populations of heterogeneous oscillators. Phys. Rev. E, 77(3):36107, 2008.
  • [6] J. Bridge, L. Mendelowitz, R. Rand, S. Sah, and A. Verdugo. Dynamics of a ring of three coupled relaxation oscillators. Communications in Nonlinear Science and Numerical Simulation, 14(4):1598–1608, 2009.
  • [7] L. M. Childs and S. H. Strogatz. Stability diagram for the forced kuramoto model. Chaos, 18(043128), 2008.
  • [8] O. Feinerman, A. Rotem, and E. Moses. Reliable neuronal logic devices from patterned hippocampal cultures. Nature Physics, 4, 2008.
  • [9] F. C. Hoppensteadt and E. M. Izhikevich. Oscillatory neurocomputers with dynamic connectivity. Phys. Rev. Lett., 82:2983–2986, 1999.
  • [10] Y. Kawamura. Chimera Ising walls in forced nonlocally coupled oscillators. Phys. Rev. E, 75:056204, 2007.
  • [11] Y. Kawamura and Y. Kuramoto. Hole structures in nonlocally coupled noisy phase oscillators. Phys. Rev. E, 76:047201, 2007.
  • [12] Y. Kuramoto. Chemical Oscillators, Waves and Turbulence. Springer, Berlin, 1984.
  • [13] Y. Kuramoto and D. Battogtokh. Coexistence of coherence and incoherence in nonlocally coupled phase oscillators. Nonlin. Phen. Complex Systems, 5(4):380–385, 2002.
  • [14] Y. Kuramoto, S. Shima, D. Battogtokh, and Y. Shiogai. Mean-field theory revives in self-oscillatory fields with non-local coupling. Progr. of Theor. Physics, 161:127, 2006.
  • [15] C. R. Laing. Chimera states in heterogeneous networks. Chaos: An Interdisciplinary Journal of Nonlinear Science, 19:013113, 2009.
  • [16] C. R. Laing. The dynamics of chimera states in heterogeneous Kuramoto networks. Physica D, 238:1569–1588, 2009.
  • [17] C. R. Laing and C. C. Chow. Stationary bumps in networks of spiking neurons. Neural Computation, 13(7):1473–1494, 2001.
  • [18] C. R. Laing and A. Longtin. Noise-induced stabilization of bumps in systems with long-range spatial coupling. Physica D: Nonlinear Phenomena, 160(3-4):149–172, 2001.
  • [19] S. D. Makovetskiy and D. N. Makovetskii. arXiv preprint, 0805(1319v1), 2005.
  • [20] S. D. Makovetskiy and D. N. Makovetskii. Emergence, competition and dynamical stabilization of dissipative rotating spiral waves in an excitable medium: A computational model based on cellular automata. arXiv nlin.CG, (0805.1319), 2008.
  • [21] E. A. Martens. Chimerae in a network of three oscillator populations with varying network topology. arXiv nlin.AO, 2010.
  • [22] E. A. Martens, E. Barreto, S. H. Strogatz, E. Ott, P. So, and T. M. Antonsen. Exact results for the Kuramoto model with a bimodal frequency distribution. Phys Rev E, 79:026204, 2009.
  • [23] E. A. Martens, C. R. Laing, and S. H. Strogatz. Solvable model of spiral wave chimeras. Phys. Rev. Lett., to appear, 2010.
  • [24] S. A. Marvel and S. H. Strogatz. Invariant submanifold for series arrays of Josephson junctions. Chaos, 19:013132, 2009.
  • [25] L. Mendelowitz, A. Verdugo, and R. Rand. Dynamics of three coupled limit cycle oscillators with application to artificial intelligence. Communications in Nonlinear Science and Numerical Simulation, 14(1), 2009.
  • [26] R. E. Mirollo, S. A. Marvel, and S. H. Strogatz. Identical phase oscillators with global sinusoidal coupling evolve by Mobius group action. Chaos, 19:043104, 2009.
  • [27] A. Motter. Spontaneous symmetry breaking. Nature Physics, 6, 2010.
  • [28] O. E. Omel’chenko, Y. L. Maistrenko, and P. A. Tass. Chimera states: The natural link between coherence and incoherence. Phys. Rev. Lett., 100:044105, 2008.
  • [29] E. Ott and T. M. Antonsen. Low dimensional behavior of large systems of globally coupled oscillators. Chaos, 18:037113, 2008.
  • [30] E. Ott and T. M. Antonsen. Long time evolution of phase oscillator systems. Chaos, 19:023117, 2009.
  • [31] A. Pikovsky and M. Rosenblum. Partially integrable dynamics of hierarchical populations of coupled oscillators. Phys. Rev. Lett., 101:264103, 2008.
  • [32] G.C. Sethia, A. Sen, and F.M. Atay. Clustered chimera states in delay-coupled oscillator systems. Phys. Rev. Lett., 100:144102, 2008.
  • [33] S. Shima and Y. Kuramoto. Rotating spiral waves with phase-randomized core in nonlocally coupled oscillators. Phys. Rev. E, 69:036213, 2004.
  • [34] J. von Neumann. Automata Studies, pages 43–98, 1954.