Chimera states through invariant manifold theory
Jaap Eldering, Jeroen S.W. Lamb, Tiago Pereira, Edmilson Roque dos Santos
CChimera states through invariant manifold theory
Jaap Eldering , Jeroen S.W. Lamb , Tiago Pereira , , andEdmilson Roque dos Santos Institute of Mathematical and Computer Sciences, University of S˜ao Paulo, Brazil Department of Mathematics, Imperial College London, London SW7 2AZ, UKE-mail: [email protected]
Abstract.
We establish the existence of chimera states, simultaneously supportingsynchronous and asynchronous dynamics, in a network consisting of two symmetricallylinked star subnetworks consisting of identical oscillators with shear and Kuramoto–Sakaguchicoupling. We show that the chimera states may be metastable or asymptotically stable. If theintra-star coupling strength is of order ε , the chimera states persist on time scales at least oforder 1 /ε in general, and on time-scales at least of order 1 /ε if the intra-star coupling is ofKuramoto–Sakaguchi type. If the intra-star coupling configuration is sparse, the chimeras areasymptotically stable. The analysis relies on a combination of dimensional reduction using aM¨obius symmetry group and techniques from averaging theory and normal hyperbolicity. a r X i v : . [ m a t h . D S ] F e b ONTENTS Contents1 Motivation and main results 22 The coupled stars network 5
In 2002, Kuramoto and Battogtokh [15] observed the coexistence of spatiotemporalsynchronous and asynchronous oscillations in a ring of identical coupled oscillators. Thisphenomenon of partial synchronisation in networks of identical coupled oscillators, wassubsequently branded chimera and henceforth observed in a wide range of experimentalsettings, including lasers [16], photoelectrochemical oscillators [9, 24] and coupledmetronomes [17]. For comprehensive reviews of the prolific occurrence of chimeras inphysical and numerical experiments, see [21, 11].
ONTENTS r + and r − for twothe star subnetworks + and − . Equations of motions are given by (1) with N = h ( ϕ ± j , ϕ ∓ i ) = sin( ϕ ± j − ϕ ∓ i + δ ) for i , j ∈ { , . . . N } , δ = . ω = λ = c is the shear parameter and ε the intra-star coupling strength. After approximately 10 cycles of natural frequency the chimera collapses to asynchrounous state. Choosing di ff erentparameter values, the chimera may collapse to a synchronous or asynchronous state. For moredetails, see Section 2.3.On the theoretical side, various mechanisms for the occurrence of chimeras have beenproposed [21, 19, 29, 1]. In the limit of infinitely many oscillators, a spectral theory forchimera states was developed [20], Notably, the comprehension of chimera states in large butfinite size networks is more challenging [28]. Carefully designed coupling functions havebeen shown to generate chimeras in finite size systems [3, 5]. Overall, despite a good twodecades of fascination with chimeras, we remain far from a comprehensive mathematicalunderstanding of their occurrence and behaviour.The main objective of this paper is to obtain mathematical results on chimera behaviourin a large but finite size network of identical oscillators, with a specific modular structureand a relatively generic type of coupling. To our best knowledge, our results are the firstexample of rigorous results in this context. The example concerns the union of two identical,symmetrically coupled star subnetworks. The key ingredients that underlie our analysis area dimensional reduction due to a M¨obius group symmetry [27], averaging theory [22] andnormally hyperbolic invariant manifold theory [8, 10].The network we consider is sketched in Figure 1. It consists of two symmetricallycoupled star subnetworks. Each such star subnetwork consists of a central hub connectedto a large number ( N ) of leaves , where all nodes are occupied by identical phase oscillatorsand interactions arise through a di ff usive Kuramoto–Sakaguchi coupling with shear. Thestrength of the coupling between nodes in each of the star motif subnetworks is uniform,while the strength of the interaction between the star motif components is represented byanother parameter ε . For more details on the equations of motion, see Section 2.1 below.An important feature of the system is that when isolated, i.e., when ε =
0, each of the
ONTENTS synchrony and another stable statein which the motion is asynchronous . In our setting the bi-stability is induced by the shear.We express the level of synchrony in terms of a complex order parameter z , see Eq. (3), relatedto the average phase di ff erence between corresponding phase oscillators in the two-star motifsubnetworks, where r = r = ε >
0) with shear parameterchosen in the bi-stability regime of the isolated dynamics and initial conditions for eachsubnetwork are chosen near the alternative stable states, so that the resulting dynamics in theabsence of coupling ( ε =
0) would display a combination of synchronous and asynchronousdynamics. In the presence of coupling ( ε > cycles of naturalfrequency) when the coupling strength ε between the star motif is small, or collapse to acompletely synchronized or asynchronous state when this coupling strength is not so small.Representative examples from our simulations are presented in Figure 1.In this paper we address mathematically the explanation of these observations forsu ffi ciently large N and su ffi ciently small ε . We obtain the following main results: (Short metastability) Chimera states exist for a time O ( ε − ) ‡ , see Theorem 4.2 [i]. (Long metastability) If the intra-star coupling is of Kuramoto–Sakaguchi type then chimerastates exist up to a time O ( ε − ), see Theorem 4.2 [ii]. (Asymptotically stable) Chimera states, which persist for all time, exist if the intra-starcoupling configuration is su ffi ciently sparse. (Theorem 4.4)The choice of our specific model was inspired by earlier observations of Ko andErmentrout [13, 12], Vlasov et al. [26] and Toenjes et al. [23] that phase oscillators in starnetworks exhibit coexisting coherent and incoherent states caused by shear-induced degreedependence of the oscillator frequency.We obtain results on the persistence of mixed chimeras states obtained by weaklycoupling two identical stars. In the uncoupled limit ( ε = M = M + × M − , where M + denotes the invariant manifold associated withsynchronous motion of the first star and M − the invariant manifold associated to theasynchronous attractor of the other star. We show that this invariant manifold is normallyhyperbolic (Theorem 4.1) which persists for su ffi ciently small intra-star coupling ( ε > M ε that is close to M . Thechimera states persist for all time only if the corresponding trajectories of the system remainon M ε for all time, which may not be the case if M ε has a boundary ∂ M ε through whichtrajectories can escape. We may observe one of the following, see also Figure 2 for a sketch: • Stable chimera: the trajectory cannot escape through the boundary ∂ M ε . • Metastable chimera: the trajectory leaves the M ε in finite time by traversing ∂ M ε . ‡ Here O ( · ) stands for Landau’s symbols for order functions. So, t = O (cid:0) δ ( ε ) (cid:1) is: there exists ε ≥ K ≥ ≤ t ≤ K | δ ( ε ) | for 0 ≤ ε ≤ ε . ONTENTS Drifting
Trapping
Figure 2: The uncoupled stars have an invariant manifold product denoted by M = M s × M I with an invariant boundary ∂ M determined by U ⊂ T N . Once both stars are coupled undersmall coupling strength ε , the invariant manifold M is perturbed to a nearby invariant manifold M ε with its boundary ∂ M ε . Simulations suggest two long-time behaviours: (i) trapping: theexistence of a chimera state holds for large finite time; (ii) drifting: chimera state ends up overfinite time due to the dynamics leaves out the invariant manifold boundary ∂ M ε .Our main results rely on estimates from normally hyperbolic invariant manifold and averagingtheory to obtain lower bounds on the drift on M ε . Before discussing the proofs of the mainresults in more detail, in the next section we first present the model in detail and presentnumerical observations that have motivated the theory developed in this paper.
2. The coupled stars network
In this section we describe the model and present numerical observations that motivated themain results obtained in this paper.
Consider two coupled star networks labelled “ + ” and “ − ”, with hub states ϕ ± ∈ S and leafstates ϕ ± i ∈ S , i = , . . . , N , where S (cid:39) R / π Z , with equations of motion˙ ϕ + = ω + λ N (cid:88) j = H ( ϕ + j , ϕ + ) , ˙ ϕ + i = ω + λ H ( ϕ + , ϕ + i ) + (cid:15) h ( ϕ + i − ϕ − i ) , (1)and likewise for the “ − ” star. Here, ω is the natural frequency of the oscillators, λ representsthe corresponding coupling strength and h is a di ff use pairwise coupling function between thetwo stars. We note that considering a general smooth h for pairwise coupling will not a ff ectthe results. The function H , representing the coupling between the hub and leaves within each ONTENTS H ( ϕ j , ϕ i ) = c + sin( ϕ j − ϕ i + δ ) , (2)where c (cid:44) δ ∈ (0 , π/
4) is the phase frustration. This coupling is relatedto the phase reduction description for a Hopf-Andronov bifurcation [14]. Shear causes thefrequency of the oscillators to depend on node degree.As a measure of synchronization, we consider the order parameters z + = N N (cid:88) j = e i ( ϕ + j − ϕ + ) z − = N N (cid:88) j = e i ( ϕ − j − ϕ − ) , (3)where r + : = | z + | and r − : = | z − | are the absolute value. We perform the variable substitution ϕ j (cid:55)→ ˜ ϕ j (cid:66) ϕ j + ω t for all j = , . . . , N and a timerescaling t (cid:55)→ τ (cid:66) λ ct . This rescales Eq. (1) and turns the degree into the hub node’s naturalfrequency, and gives ˙ ϕ + = N + c N (cid:88) j = sin( ϕ + j − ϕ + + δ ) , ˙ ϕ + i = + c sin( ϕ + − ϕ + i + δ ) + (cid:15)λ c h ( ϕ + i − ϕ − i ) , (4)where we abuse notation by re-using variables ϕ and denoting di ff erentiation with respect to τ by · again. Let us denote σ = / c , ε = (cid:15)σ/λ and allow for a slight generalization, where β > k i = i = , . . . , N . When β = N we recover theprevious Eq. (4). Here β works as parameter for the frequency gap. So, the network dynamicswe analyze is given as˙ ϕ + = β + β σ N N (cid:88) j = sin( ϕ + j − ϕ + + δ )˙ ϕ + i = + σ sin( ϕ + − ϕ + i + δ ) + ε h ( ϕ + i − ϕ − i ) , i = , . . . , N , (5)and likewise for the “ − ” star. We begin by reviewing the synchronization diagram for a single star and then we proceed tonumerically explore chimera states for two coupled stars.
ONTENTS . . . . . . . σ . . . . . . r σ b σ f forwardbackward Figure 3: Synchronization diagram for a single star with N =
200 leaves. The parametersare β = δ = . ∆ σ = .
02. The region between the two vertical dashed lines is where both coherent andincoherent states coexist. The green solid line corresponds to the separatrix.
We first obtain the bifurcation diagram of the order parameteras a function of σ . Figure 3 shows the synchronization diagram for a star. Numerical procedure.
We integrate Equations (5) of one star using an implicit Runge–Kutta method of order five (we use a solver described in [2]). Starting at σ = , π ) and evaluate the order parameter r in the stationary regime.Then we increase adiabatically the coupling by ∆ σ = .
02 and, using the outcome of thelast run as the initial condition, calculate the new value of the stationary order parameter r at σ + ∆ σ , repeating these steps until a maximal value σ is reached. This curve is called theforward continuation. We can see that the order parameters increases slowly and smoothlyuntil σ = σ f , where it jumps discontinuously to r = ∆ σ from the synchronous state. We use the outcome of the last run as theinitial condition for the calculation of r we add to the initial conditions a small random uniformnumber drawn from the interval (0 , . Main findings.
In Section 3.2.2, we prove the existence of such synchronization diagramvia the theory of invariant manifolds. Here we summarize the main points. In Figure 3, theregion between the dashed lines corresponds to the values of the coupling critical couplings σ b = β − + β cos(2 δ ) and σ f = β − (cid:112) + β cos(2 δ ) . The backward critical point σ b is calculated in Section 3.3, where we prove that it is aconsequence of requiring the existence and stability conditions for the synchronous state. The ONTENTS σ f is calculated in Section 3.4 analyzing the di ff erential equation of theM¨obius group α parameter. This allows us to obtain the solid lines as well. The black solidline corresponds to the − sign and the green solid line to the + in the following equation( β − ± (cid:112) ( β − − σ (1 + β cos(2 δ )) σ (1 + β cos(2 δ )) . Moreover, in Section 3.4.4 we prove through the averaging principle that for small couplingstrength σ <
We consider two stars coupled as inFigure 1. Each star has N leaves, with the coupling strength between its leaves and its hubgiven by σ . The coupling between them is mediated by connecting, with strength ε , the k thleave of star + to the k th leaf of star − .The mechanism for generating a chimera state is a consequence of the coexistencebetween synchronous and incoherent states. From identical stars, we select a star to be atthe synchronous state and another to be at the incoherent. Both stars are coupled under aninter-coupling strength ε and function h ( ϕ ± j , ϕ ∓ i ) = sin( ϕ ± j − ϕ ∓ i + δ ) , i , j ∈ { , . . . N } . Starting from this prepared initial conditions, and from the arguments in Section 4, we showthat the coupled system presents similar macroscopic behavior as long as the perturbation caused by the coupling is small enough.Examples of such chimera-like state are depicted in Figure 1. In Figure 1 where everyleaf of the top star is connected to one of the bottom star, with strength ε = . breathing chimera [1]. Star motifs are building blocks for complexnetworks with scale-free networks [4]. Although our results for star graphs do not apply tothe scale-free network, we can use them as heuristics to build chimera-like states in scale-freenetworks. Consider the following equations of motion on complex networks,˙ ϕ i = k i + σ N (cid:88) j = A i j sin (cid:16) ϕ j − ϕ i + δ (cid:17) , i = , . . . , N , (6)where σ is the coupling strength, k i is the degree of the vertex i and A i j is the adjacencymatrix. The synchronization diagram for a Barab´asi–Albert [4] network with 1000 vertices,degree distribution exponent γ ≈ . (cid:104) k (cid:105) =
6, performed with the samemethod as was used in Figure 3 is shown in Figure 4.We coupled two Barab´asi–Albert networks with the same degree distribution exponent γ and mean degree (cid:104) k (cid:105) . The k th vertex of one network is coupled to (and only to) to the k th ONTENTS . . . . . . . σ . . . . . . r forwardbackward Figure 4: Synchronization diagram for a single Barab´asi–Albert network with N = γ ≈ . (cid:104) k (cid:105) = δ = .
03, and at each step the coupling strength was increased(decreased) by ∆ σ = .
02. There is a coherent and incoherent states coexistence region.
200 250 300 t . . . r − , r + σ = 1 . ε = 0 .
50 75 100 tσ = 1 . ε = 0 . Figure 5: Time series for the order parameters r + and r − for two coupled Barab´asi–Albertnetworks with N = γ ≈ . (cid:104) k (cid:105) =
6. Weindicate the intra-coupling strength σ and inter-coupling strength ε .vertex of the other network. The coupling function used is the Kuramoto–Sakaguchi coupling.Again σ denotes the intra-coupling strength and the ε inter-coupling. Figure 5 shows chimera-like states in the two coupled Barab´asi–Albert networks. We chose initial conditions insidethe hysteresis loop in Figure 4. The chimera breaking is observed in the right column.
3. Mathematical analysis
In the remainder of this paper, the main results, as informally discussed above in Section 1and summarized in Theorems 4.1-4.4, are proved. To this end we sketch the strategy to prove
ONTENTS
First Step:
Constructing the Synchronous Manifold of a single star . There exists a manifoldwhere all phases are locked. In Section 3.3 we prove that for certain intervals of theparameters the manifold is normally attracting, see Proposition 3.5.
Second Step:
Asynchronous manifold for a single star.
We construct it as follows:
Reduction in terms of M¨obius actions.
We rewrite the model in terms of relativephases. In these coordinates, dynamics are described in a low-dimensional manner bythe action of a three-parameter M¨obius group, see Section 3.4.1. The Watanabe–Strogatzapproach provides a di ff erential equation for the parameters of the M¨obius group α ∈ C and ψ ∈ S keeping invariant a set of constants of motion. Invariant Manifold for asynchronous dynamics.
To construct the manifold for theasynchronous dynamics we perform the following steps:(i) Construct open sets of initial conditions such that the order parameter z is C k close to α , see Section 3.4.2 in Lemma 3.7. The di ff erential equation for α depends on z and by a small perturbation we obtain an equation solely in terms of α , see Section 3.4.3in Theorem 3.8. By averaging theory, in Section 3.4.4, we prove that the perturbeddynamics of α enters in finite time a neighborhood of zero, see Theorem 3.9. Thisdetermines that the bifurcation diagram of a single star network is in fact discontinuousfor su ffi ciently small coupling strength.(ii) By persistence of hyperbolic manifolds we obtain the invariant manifoldcorresponding to the asynchronous dynamics in the full equations, see Proposition 3.13. Third Step:
Constructing Chimeras in coupled stars . For the coupled stars, the constantsof motion may drift to the boundary of the asynchronous manifold. When this happenswe lose control over the dynamics because the order parameter may no longer be in aneighborhood of α . So first, we prove the existence of chimera states in the couple starssystem, see Theorem 4.1. And finally, we show three distinct results about the time scalethese chimera states exists depending on how generic we couple both stars. Short metastability.
We estimate an upper bound for how long it takes the perturbeddynamics of the constants of motion drift out the invariant manifold.
Long metastability.
We assume the coupling function between both stars hassinusoidal form as an extra assumption. This implies the coupling term can be splitup into two terms and restricted to the invariant manifold. This restriction introduces ane ff ective coupling term of order ε to the constants of motion dynamics. Asymptotically stable metastability.
We show that coupling a single leaf maygenerate a drift of the conserved quantities, but this will happen only for a single constant.That is, the drift does not cascade to other leaves. Hence, if only a small fraction of leavesare coupled this drift is immaterial and the order parameter will stay in a C k neighborhoodof α and the chimera state is stable. ONTENTS In preparation for the subsequent analysis, we briefly recall some established results from thetheory normally hyperbolic invariant manifolds and averaging.
Our first ingredient is the theoryof normally hyperbolic invariant manifolds [8, 10, 7]. Such manifolds can be viewed asgeneralizations of hyperbolic fixed points and we use their property that they persist undersmall perturbations of the dynamics. We shall only consider normally attracting invariantmanifolds (NAIMs), i.e. those which only have stable normal directions, not unstable ones.Let ˙ x = f ( x ) with x ∈ R n define a dynamical system; we use R n for simplicity ofpresentation, but it can be replaced by a smooth manifold. A submanifold M ⊂ R n is saidto be a normally attracting invariant manifold for f if it is invariant under the flow Φ t of f and furthermore the linearized flow along M contracts normal directions, and does so morestrongly than it contracts directions tangent to M . Our precise definition below correspondsto that of “eventual absolute normal hyperbolicity” according to [10]. Definition 3.1 (Normally attracting invariant manifold) . Let ˙ x = f ( x ) with x ∈ R n and f ∈ C describe a dynamical system with flow Φ : R × R n → R n . Then a submanifold M ⊂ R n is calleda normally attracting invariant manifold of f if all of the following conditions hold true:(i) M is invariant, i.e. ∀ t ∈ R : Φ t ( M ) = M;(ii) there exists a continuous splitting T M R n = T M ⊕ N (7) of the tangent bundle T R n over M into the tangent and a (stable) normal bundle, withcontinuous projections π M , π N and this splitting is invariant under the tangent flowT Φ t = T Φ tM ⊕ T Φ tN ;(iii) there exist real numbers a < b ≤ and C > such that the following exponential growthconditions hold on the subbundles: ∀ t ≤ , ( m , v ) ∈ T M : (cid:107) T Φ tM ( m ) v (cid:107) ≤ C e b t (cid:107) v (cid:107) , ∀ t ≥ , ( m , v ) ∈ N : (cid:107) T Φ tN ( m ) v (cid:107) ≤ C e a t (cid:107) v (cid:107) . Note that if we have two dynamical systems f , f which each have a NAIM M , M respectively and max( a , a ) < min( b , b ), then M × M is a NAIM for the product system.In our case, the dynamics on the invariant manifolds will be (close to) neutral, so b , ≈ f . That is, we have the following theorem, see [8, Thm. 1] or [10, Thm. 4.1]. Theorem 3.2 (persistence of NAIMs) . Let M ⊂ R n be a compact normally attracting invariantmanifold for f . Then there exists an ε > such that for any vector field ˜ f with (cid:107) ˜ f − f (cid:107) C ≤ ε ,there exists a unique manifold ˜ M that is di ff eomorphic and O (cid:0) (cid:107) ˜ f − f (cid:107) C (cid:1) -close to M andinvariant under ˜ f . Furthermore, ˜ M is a NAIM for ˜ f .ONTENTS M to M can in general be expressed as follows: ˜ M can be described bythe graph of a section of the normal bundle of M , and this section is C -small. In our case M can explicitly be given by the graph of a function, hence ˜ M will be given by the graph of a C -small perturbation of this function.In our case the NAIM M will have a boundary. Theorem 3.2 still holds true as long asthe vector field f is pointing strictly outward at the boundary ∂ M (also called ‘overflowinginvariant’), see [8]. In our case M will not be overflowing invariant. This can be overcomeby artificially modifying f near the boundary ∂ M . A slightly simpler approach uses [7,Thm. 4.8]: we modify f such that it preserves a manifold S that transversely intersects M with S ∩ M = ∂ M . For example, S could be the local normal bundle of M restricted to ∂ M . The averaging principle will play role in our analysis. Considerthe original system given by˙ x = ε F ( x , η ) , ˙ η = Ω ( x , η ) + ε G ( x , η ) , ( x (0) , η (0)) = ( x , η ) , x ∈ D , η ∈ S . (8)and the averaged system given by˙ y = ε F ( y ) , ˙ ζ = Ω ( y , ζ ) , ( y (0) , ζ (0)) = ( x , ζ ) ∈ U × S , (9)where F ( y ) = π (cid:82) π F ( y ,ϑ ) Ω ( y ,ϑ ) d ϑ π (cid:82) π d ϑ Ω ( y ,ϑ ) is the first order averaged vector field. If the averaged system has an asymptotically stablefixed point, solutions of the original system are approximated by averaged ones. That is,we have the following adapted version of the theorem, see [6, 5, Thm. 7.2] or [22, 7,Lemma. 7.10.1-2, Thm.7.10.4]. Theorem 3.3. [Averaging theorem with attraction ]
Let D ⊂ R n be an open set and considerthe equation (8) . Suppose that the averaged system Equation (9) has an asymptotically stablefixed point y = , F ∈ C ( D ) and has a domain of attraction D ⊂ D. For any compactK ⊂ D there exists a ε and c > such that for all x ∈ K and for each ε < ε (cid:107) x ( t ) − y ( t ) (cid:107) ≤ c ε, ≤ t < ∞ . Bifurcation diagram of an isolated star network
The synchronous state is characterized by r =
1. The phase locking manifold M φ = { φ = . . . = φ N = φ ∈ S , φ − ϕ = ∆ φ ∈ S } ⊂ T N + satisfies this condition on the order parameter. The stability of M φ depends on the frequencymismatch, the phase frustration δ and intra-coupling strength σ [25]. The followinghypothesis is necessary to prove existence of the coherent manifold in Proposition 3.5. ONTENTS Hypothesis 3.4.
Assume that σ > σ b : = β − + β cos(2 δ ) and arg ( v ) − sin − (cid:16) − βσ (cid:107) v (cid:107) (cid:17) ∈ (cid:16) − π , π (cid:17) , (10) where v = ( β sin(2 δ ) , + β cos(2 δ )) ∈ R , and hence (cid:107) v (cid:107) = (cid:112) + β + β cos(2 δ ) . Proposition 3.5.
Let σ > , β > and δ ∈ (0 , π ) . The synchronized stateM C = { φ = . . . = φ N = φ C ∈ S } (11) is a normally attracting invariant manifold for the dynamics (16) if and only if Hypothesis 3.4holds and when φ C is given by φ C = δ − π + arctan (cid:32) + β cos(2 δ ) β sin(2 δ ) (cid:33) + arccos (cid:32) β − σ (cid:107) v (cid:107) (cid:33) . (12) Proof.
First we search for the value of φ = φ C such that M C is an invariant manifold.Evaluating (16) on M C we see that this amounts to˙ φ = − β − σ sin( φ + δ ) − βσ sin( φ − δ ) = . This can be rewritten as (cid:104) v , (cid:0) cos( φ C − δ ) , sin( φ C − δ ) (cid:1) (cid:105) = − β − σ (13)with v = ( β sin(2 δ ) , + β cos(2 δ )). Hence we have solutions if (cid:107) v (cid:107) = (cid:112) + β + β cos(2 δ ) ≥ β − σ . (14)To evaluate stability of M C , let ( ˙ φ , . . . , ˙ φ N , ˙ ϕ H ) = F ( φ , . . . , φ N , ϕ ) denote thesystem (16). Then we have the total derivative DF | M C = − σ cos( φ C + δ ) ∅ . . . ∅ − σ β N cos( φ C − δ ) · · · ... . . . ... ... · · · − · · · − . Thus, we find that (0 , . . . , ,
1) is an eigenvector with eigenvalue 0, that is, the direction along M C is neutral. In the directions transversal to M C we have the eigenvector v = (cid:16) , . . . , , − β cos( φ C + δ )cos( φ C − δ ) + β cos( φ C + δ ) (cid:17) with eigenvalue λ = − σ (cid:2) cos( φ C − δ ) + β cos( φ C + δ ) (cid:3) and N − v i = (0 . . . , , − , . . . ,
0) with eigenvalue λ = − σ cos( φ C − δ ). So M C is normally attractingif σ > φ C − δ ) > φ C − δ ) + β cos( φ C + δ ) > . ONTENTS (cid:104) w , (cid:0) cos( φ C − δ ) , sin( φ C − δ ) (cid:1) (cid:105) > , (15)with w = (1 + β cos(2 δ ) , − β sin(2 δ )) equal to v rotated clockwise over 90 degrees. Interpretingthe conditions (13) and (15) geometrically with φ C ∈ S , see Figure 6, we see that we haveat most two synchronized states. Since v ⊥ w , precisely one satisfies the second stabilitycondition of (15) if and only if (14) holds as a strict inequality. Then that state M C is actuallystable if cos( φ C − δ ) > φ C yields (12).The backward critical coupling σ b is obtained solving equation cos( φ C − δ ) = φ C expression obtained. This yields the first condition at (10). (cid:3) v w δ − ( β − ωσ (cid:107) v (cid:107) Figure 6: A geometric representation of the existence and stability of the synchronized state.The dots indicate the fixed points for φ C ∈ S and the red arc is the unstable region. Thevectors v and w and the two dashed lines are perpendicular. Again let us consider a single star network ( ε =
0) in system (5). We change to a new set ofcoordinates defined by φ k = ϕ k − ϕ . The φ k measure the phase di ff erence between the leavesand the hub, and the resulting equations are˙ φ i = − β − σ sin( φ i − δ ) − β σ N N (cid:88) j = sin( φ j + δ ) , ˙ ϕ = β + β σ N N (cid:88) j = sin( φ j + δ ) . (16)In the original coordinates the system had a global S -symmetry by which we have e ff ectivelyreduced the system, since the equations for the relative angles φ i are now decoupled from ϕ .To find a stable asynchronous state, we shall change the equations for a single star fromrelative phase form (16) to new variables introduced by Watanabe–Strogatz [27]. ONTENTS Watanabe and Strogatz [27], showed that a large classof systems of N identical coupled oscillators can be described by a coupled set of only threeODEs. They achieved this through a time dependent change of variables, that explicitly showsthat the system possesses N − N identical coupled phase oscillators with equations of motion that can be put in the form˙ φ j = f e i φ j + g + ¯ f e − i φ j . (17)Here f is a smooth complex-valued function of the phases φ , φ , . . . , φ N ∈ S , ¯ f is its complexconjugate, and g is in its turn a real-valued function of the phases φ , φ , . . . , φ N ∈ S . Thephases φ j evolve in time according to the action of the M¨obius § group, e i φ j ( t ) = M α ( t ) ,ψ ( t ) ( e i θ j ) , (18)where θ j are a set of constant (time independent) angles, and we use the group parametrization M α,ψ and its action on w ∈ C given by M α,ψ ( w ) = α + e i ψ w + ¯ α e i ψ w , (19)with ψ ∈ S and α ∈ D = { z ∈ C | | z | < } .The time evolution of the phases in (18) satisfies the equations (17) as long as α ( t ) and ψ ( t ) are solutions of the di ff erential equations˙ α = i (cid:16) f α + g α + ¯ f (cid:17) , ˙ ψ = f α + g + ¯ f ¯ α. (20)The functions f and g still implicitly depend on φ , φ , . . . , φ N ∈ S , so we have to use (18) toexpress the φ j ’s in terms of α, ψ , and the θ j ’s and obtain a closed system of equations for α, ψ .For suitable models, we can express f and g in terms of the order parameter z , i.e. f = f ( z , ¯ z )and g = g ( z , ¯ z ), where z is given by z ( α, ψ, θ ) : = N N (cid:88) j = e i φ j = N N (cid:88) j = M α,ψ (cid:0) e i θ j (cid:1) . (21)We now reformulate the Watanabe–Strogatz change of variables in a more formalgeometric setting. Let ˙ φ = X ( φ ) with X ∈ X ( T N ) denote the vector field of (17) and let( ˙ α, ˙ ψ, ˙ θ ) = ˆ X ( α, ψ, θ ) with ˆ X ∈ X ( Q ) denote the vector field associated to the Watanabe–Strogatz coordinates ( α, ψ, θ ) ∈ Q : = D × S × T N , i.e. we have ˙ θ = α, ˙ ψ given by (20).Let Φ t and ˆ Φ t denote the associated flows. Then we have the commuting diagram Q π (cid:15) (cid:15) ˆ Φ t (cid:47) (cid:47) Q π (cid:15) (cid:15) T N Φ t (cid:47) (cid:47) T N (22) § The complete M¨obius group consists of all fractional linear transformations M ( w ) = aw + bcw + d of the complexplane such that ad − bc (cid:44) ONTENTS π : Q → T N is the submersion φ = π ( α, ψ, θ ) defined by (18). Since π ( α, ψ, · ) : T N → T N is given by the diagonal action of M α,ψ , it is a di ff eomorphism and hence π indeed asubmersion. Note that although we loosely speak of π as a coordinate transformation, strictlyspeaking it is not, since it is not injective. The diagram (22) implies that T π ◦ ˆ X = X ◦ π andhence X = T π ◦ ˆ X ◦ π − for any right-inverse π − . On the other hand, a vector field Y ∈ X ( T N )can be lifted to a vector field ˆ Y = R ◦ X ◦ π ∈ X ( Q ) satisfying (22), given a choice of aright-inverse R of T π . We now choose R ( α, ψ, θ ) = T θ (cid:2) π ( α, ψ, · ) − (cid:3) : T φ T N → T θ T N ⊂ T ( α,ψ,θ ) Q (23)with φ = π ( α, ψ, θ ). Let M (cid:48) α,ψ denote the derivative of the action on S , then we have R ( α, ψ, θ ) = diag (cid:0) M (cid:48) α,ψ ( θ ) − , . . . , M (cid:48) α,ψ ( θ N ) − (cid:1) . (24)The upshot of this is that if we perturb the original vector field X to X + ε Y , for example whenintroducing a coupling between the stars, then we can lift Y to ˆ Y ∈ X ( Q ) such that ˆ Y is givenby ˙ α = , ˙ ψ = , ˙ θ j = M (cid:48) α,ψ ( θ j ) − Y j ( π ( α, ψ, θ )) . (25)In the case that only some of the components Y j are nonzero, it follows that only the associated θ j ’s are not conserved anymore by ˆ Y , while the other θ j ’s are still conserved.Note that formalizing the Watanabe–Strogatz method itself this way, leads to a di ff erentlifted vector field; indeed, with our choice of lift we would obtain ˙ α = ψ = X ∈ X ( T N ) is given by˙ φ j = X j ( φ ) = g ( φ ) + f ( φ ) e i φ j + f ( φ ) e − i φ j , (26)then the Watanabe–Strogatz approach tells us that ˆ X ∈ X ( Q ) given by (20) and ˙ θ = X . If X depends on extra parameters then this directly translates to the same parameterdependence for ˆ X . The upshot is that we can apply di ff erent lifts to parts of a vector field X ,thereby obtaining the most convenient choice of lifted vector field ˆ X . We write the sine functionsin (16) in exponential form. This shows that we can recast the equations in the Watanabe–Strogatz form (17) with f and g given by f ( z , ¯ z ) = − σ e − i δ i , g ( z , ¯ z ) = − β − σβ (cid:32) ze i δ − ¯ ze − i δ i (cid:33) . (27)For a more detailed analysis of a stable incoherent state in the system, we shall use that z ≈ α when the θ ’s are close to being uniformly distributed on S (called ‘splay states’), andwhen we are close to the thermodynamic limit, i.e. when the number of leaves N is su ffi cientlylarge. Let us make this more precise. ONTENTS Definition 3.6.
We say that θ ∈ T N is uniformly distributed if there exists a ϕ ∈ S such that θ j = π jN + ϕ up to permutation. We denote the set of all such θ by Θ . Lemma 3.7.
Fix < r < and consider the disc D r = { α ∈ C | | α | < r } . For each k ≥ and ε > there exists for all N su ffi ciently large an open neighborhood U of Θ ⊂ T N suchthat the map D r → C , α (cid:55)→ z ( α, ψ, θ ) − α , with z ( α, ψ, θ ) given by (21) , is smaller than ε inC k -norm, uniformly for all ( ψ, θ ) ∈ S × U.Proof.
Consider an element θ ∈ Θ and wlog. assume that it is of the form as in Definition 3.6without permutation of the indices j and that these run from 0 to N − + x ) − = (cid:80) l ≥ ( − x ) l with x = ¯ α e i ψ e i θ j ,which converges for any | α | <
1. We find z ( α, ψ, θ ) = N N − (cid:88) j = α + e i ψ e i θ j + ¯ α e i ψ e i θ j = N N − (cid:88) j = α (cid:16) + α − e i ( π jN + ϕ + ψ ) (cid:17) (cid:88) l ≥ (cid:16) − ¯ α e i (cid:0) π jN + ϕ + ψ (cid:1)(cid:17) l = α N N − (cid:88) j = (cid:16) + α − e i ( π jN + ϕ + ψ ) (cid:17)(cid:32) N − (cid:88) l = ( − ¯ α ) l e i (cid:0) π jlN + l ( ϕ + ψ ) (cid:1) + R N − , j ( α ) (cid:33) = α N N − (cid:88) l = N − (cid:88) j = (cid:104) ( − ¯ α ) l e i (cid:0) π jlN + l ( ϕ + ψ ) (cid:1) + ( − ¯ α ) l α e i (cid:0) π j ( l + N + ( l + ϕ + ψ ) (cid:1)(cid:105) + α N N − (cid:88) j = (cid:16) + α − e i ( π jN + ϕ + ψ ) (cid:17) R N − , j ( α ) = α N N − (cid:88) l = δ l , N − (cid:88) j = ( − ¯ α ) l e i (cid:0) π jlN + l ( ϕ + ψ ) (cid:1) + α N N − (cid:88) j = (cid:16) + α − e i ( π jN + ϕ + ψ ) (cid:17) R N − , j ( α ) = α + N N − (cid:88) j = (cid:16) α + e i ( π jN + ϕ + ψ ) (cid:17) R N − , j ( α ) . To obtain the stated result we have to bound the second term in C k -norm. We have R N − , j ( α ) = (cid:88) l ≥ N − (cid:16) − ¯ α e i (cid:0) π jN + ϕ + ψ (cid:1)(cid:17) l = ( − ¯ α ) N − (cid:88) l ≥ ( − ¯ α ) l (cid:16) e i (cid:0) π jN + ϕ + ψ (cid:1)(cid:17) l + N − and note that the sum still defines a function with radius of convergence one, so this functionand its derivatives up to order k are uniformly bounded for α ∈ D r by some number B k > (cid:107) R N − , j (cid:107) C k ≤ kr N − B k . Since the function F j ( α ) = α + e i ( π jN + ϕ + ψ ) multiplying R N − ( α ) ONTENTS C -norm bounded by 2 + r and all higher derivatives zero, we find (cid:13)(cid:13)(cid:13)(cid:13) N N − (cid:88) j = F j · R N − , j (cid:13)(cid:13)(cid:13)(cid:13) C k ≤ max j ∈ [0 , N − k (cid:88) m = (cid:16) (cid:107) F j (cid:107) C (cid:107) R N − , j (cid:107) C k + m (cid:107) F j (cid:107) C (cid:107) R N − , j (cid:107) C k − (cid:17) ≤ (cid:16) k + k ( k + (cid:17) max j ∈ [0 , N − (cid:107) R N − , j (cid:107) C k ≤ k ( k + + r ) r N − B k . Thus, for fixed r < k ≥ ε >
0, we see that there exists an N > N ≥ N (cid:107) α (cid:55)→ z ( α, ψ, θ ) − α (cid:107) C k < ε D r × S × Θ . Next, since this function is C ∞ (and actually analytic w.r.t. α, ¯ α, ψ, θ ) it follows that there existsan open neighborhood U ⊃ Θ such that (cid:107) α (cid:55)→ z ( α, ψ, θ ) − α (cid:107) C k < ε on D r × S × U , which completes the proof. (cid:3) The Watanabe–Strogatz equationswith functions f , g given by (27) depend on α, ψ (and θ as constant of motion) through thefunction z ( α, ψ, θ ) in (21). This makes it di ffi cult to find fixed points and analyse their stability.Our approach is to first make the assumption that z = α , thereby ‘closing the equation’ to asimple form. Secondly, after we find a stable fixed point, we use Lemma 3.7 to prove that itpersists after reinserting the original function z ( α, ψ, θ ).Let us therefore start analyzing the equations (20) with f and g given by (27) but z replaced by α . In this case the system reduces to a skew-product flow; the equation for α doesnot depend on ψ anymore and becomes˙ α = − σ (cid:16) e − i δ + β e i δ (cid:17) α + i (1 − β ) α + σ (cid:16) β | α | e − i δ + e i δ (cid:17) . (28)with 0 < δ < π/
4. The next proposition characterises the asynchronous branch.
Theorem 3.8.
Let < σ < σ f : = β − (cid:112) + β cos(2 δ ) . (29) Then equation (28) has an exponentially attracting fixed point α I given by | α I | = β − − (cid:112) ( β − − σ (1 + β cos(2 δ )) σ (1 + β cos(2 δ )) and arg( α I ) = − π + δ. (30) Moreover, α I is the unique branch of fixed points satisfying lim σ → α I = .ONTENTS Proof.
Representing equation (28) in polar coordinates, α = re i η , we obtain˙ r = σ − r ) cos( η − δ ) , (31a)˙ η = − β − σβ r sin( η + δ ) − σ + r r sin( η − δ ) . (31b)Because of the product structure of equation (31a), we readily see that fixed points correspondto either r = η − δ ) =
0. Since we are looking for the branch α I of fixed points thatconverge to zero as σ →
0, we discard the solution r =
1. Hence we substitute cos( η − δ ) = η − δ ) = ± ± ( β − + σβ r cos(2 δ ) + σ + r r = . Solving for r , we note that sin( η − δ ) = η = − π/ + δ and r ± = ( β − ± (cid:112) ( β − − σ (1 + β cos(2 δ )) σ (1 + β cos(2 δ )) . (32)These fixed points exist for σ < σ f = ( β − (cid:112) + β cos(2 δ ) , where the square root is well-defined since δ ∈ (0 , π ). The branch r − uniquely satisfies thecondition that it converges to zero as σ → J ( r − , − π + δ ) = (cid:32) σ (1 − r − ) a − βσ r − sin(2 δ ) (cid:33) with a = σ (cid:34) + β cos(2 δ ) − r − (cid:35) . (33)Considering (32) and interpreting ( β −
1) and σ (cid:112) + β cos(2 δ ) as the sides of a triangle, weapply a triangle inequality estimate and conclude that1 + β cos(2 δ ) < r − . Therefore, it follows that for δ ∈ (0 , π ) and σ > a < J ) < J ) > , hence both eigenvalues of J have negative real part, which completes the proof. (cid:3) We are interested in the large β regime. In thiscase, σ bc → β → ∞ . The Jacobian eigenvalues at the incoherence state α I are complexand their imaginary part is proportional to β whereas their real parts converge to a constantdepending on δ and σ . As a consequence, near the fixed point solutions spiral towards thefixed point with high frequency. We can also construct a Lyapunov function near α I but itsdomain is rather small. This shows that the nonlinear picture for | α | ≈ α we average out the fast oscillatory behaviour andanalyze the averaged system. ONTENTS Theorem 3.9.
Let z be the order parameter of the system (21) corresponding to initialconditions in an open neighborhood U of Θ ⊂ T N . Given k ∈ N , δ > , < δ < π/ , < σ < and ε > small enough, there is N ∈ N , β > such that for all N > N , β > β ,and z (0) ∈ D δ there is T > such that in τ : = β t | z ( τ ) − α I | ≤ ε, for all τ > T . The proof of this theorem is presented below. Roughly speaking, it means that if wefollow the synchronized manifold M C decreasing slowly the intra-coupling strength once itlooses stability the order parameter will drop and stay near zero. Thus, only the asynchronousbranch is stable for 0 < σ < β . We will obtain this result as a consequence ofthree Lemmas. To this end we first need some preliminary setting.The planar system (31) in polar coordinates can be rewritten in the form˙ r = F ( r , η ) and ˙ η = G ( r , η )where F ( r , η ) : = σ − r ) cos( η − δ ) , (34a) G ( r , η ) : = − β [1 + σ r sin( η + δ )] + (cid:32) − σ + r r sin( η − δ ) (cid:33) . (34b)We can parametrize time using the parameter β and introduce an averaged system on the newtime τ : = β t . First, we need some auxiliary results. Lemma 3.10.
Fix δ ∈ (0 , and consider R : [ δ − , − δ ] → R asR ( x ) = P ( x ) Q ( x ) , where P ( x ) = − π (cid:90) π sin ϑ + x sin ϑ d ϑ and Q ( x ) = π (cid:90) π d ϑ + x sin( ϑ ) . Then R is a smooth function satisfyingR (0) = , R ( x ) > for x ∈ (0 , − δ ] , and dR (0) dx = . Proof.
First we calculate P ( x ). We Taylor expand the function 1 / (1 + x sin ϑ ) to obtain P ( x ) = − π (cid:90) π sin ϑ ∞ (cid:88) n = ( − n ( x sin ϑ ) n d ϑ notice that (cid:82) π sin k + ϑ d ϑ = k + ϑ is an odd function in [0 , π ] and thus − π (cid:90) π sin( ϑ )1 + x sin( ϑ ) d ϑ = ∞ (cid:88) k = a k x k − ONTENTS a k = π (cid:90) π sin k ϑ d ϑ > , k ≥ P is a monotone function. Moreover, a = , a = , and a = . Repeating the same procedure Q ( x ) = π (cid:90) π ∞ (cid:88) n = ( − n ( x sin ϑ ) n d ϑ = + ∞ (cid:88) k = a k x k = + x P ( x ) . Therefore, R ( x ) is positive, and the claim follows. (cid:3) This Lemma will play a role in the averaging approximation for the dynamics of α Lemma 3.11.
Fix < δ < and consider D δ : = { z ∈ C : | z | < − δ } . For any < σ < , < δ < π/ , α ∈ D δ there exist constants c > , and β > such that for each β > β thesolutions α ( t ) of Equation (34) in τ : = β t with initial condition α ∈ D δ satisfy (cid:107) r ( τ ) − ρ ( τ ) (cid:107) < c β for τ ≥ , (35) where ρ ( τ ) is the trajectory of the averaged system ρ (cid:48) = β F ( ρ ) (36) with initial condition ρ ( η ) = r , where F : [ δ − , − δ ] → R is smooth and satisfyingF (0) = , F ( ρ ) < for ρ ∈ (0 , − δ ] and dF (0) d ρ = − σ sin 2 δ. Proof.
Performing a change of time scale t (cid:55)→ τ = β t , we obtain: r (cid:48) = ε F ( r , η ) and η (cid:48) = − − σ r sin( η + δ ) + ε G ( r , η ) , (37)where (cid:48) denotes di ff erentiation with respect to τ and ε = /β . The transformation induces aslow-fast structure where the small parameter is ε . Notice that r = r andthus solutions starting with r < r > r = η (cid:48) is bounded away from zero and we can use η as our new time. To apply the averagingprinciple, we define the averaged vector field as F ( r ) : = π (cid:82) π H ( r , ϑ, d ϑ π (cid:82) π d ϑ − − σ r sin ( ϑ + δ ) , (38) ONTENTS H ( r , η, ε ) = σ (1 − r ) cos( η − δ ) − [1 + σ r sin( η + δ )] + ε (cid:16) − σ + r r sin( η − δ ) (cid:17) . (39)The denominator of Equation (38) resembles the Q ( x ) function of Lemma 3.10. So, replacingEquation (39) into Equation (38), we calculate the numerator: (cid:90) π cos( η − δ )1 + σ r sin( η + δ ) d η = cos(2 δ ) (cid:90) π + δδ cos( ϑ )1 + σ r sin( ϑ ) d ϑ + sin(2 δ ) (cid:90) π + δδ sin( ϑ )1 + σ r sin( ϑ ) d ϑ, where we used the change of variables ϑ = η + δ and the trigonometric relation for the cosine.Since both integrals are along a full cycle of ϑ , the first integrals is zero. Then − π (cid:90) π cos( η − δ )1 + σ r sin( η + δ ) d η = sin(2 δ ) P ( σ r ) , From this observation we obtain the vector field for ρ . Thus we obtain F ( ρ ) = − sin 2 δ σ (1 − ρ ) R ( σρ ) . From Lemma 3.10 it follows that since R ( σρ ) > ρ > F ( ρ ) < F < ρ < δ su ffi ciently small. We showthat this convergence is exponential as dF (0) d ρ <
0. The result follows by the principle oflinearized stability: given δ > ρ < δ converge to theorigin exponentially fast. Next applying averaging Theorem 3.3 our claim follows. (cid:3) Lemma 3.12.
Fix < δ < . Then, for any < σ < , < δ < π/ , α ∈ D δ , there exists β > such that for each β > β the solutions α ( t ) of Equation (34) with initial condition α ∈ D δ converge to α I .Proof. Given ε > c and β suchthat for every β > β we have (cid:107) α ( τ ) (cid:107) = (cid:107) r ( τ ) (cid:107) ≤ ε , since 0 is a hyperbolic fixed point for ρ .Equation (34) has a stable fixed point α I satisfying | α I | → β → ∞ as can be observedin Equation (32). Moreover, when β → ∞ the real part of the Jacobian eigenvalues of α I tendsto a constant independent of β . Thus, by the principle of linearized stability there is δ I > β such that any initial condition in the open ball B ( α I , δ I ) converges to α I .Now we choose ε < δ I . This is guaranteed as δ I is independent of β . Thus, when β is largeenough, solutions α starting in D δ enter B ( α I , δ I ) and thus must converge to α I . (cid:3) Proof of Theorem 3.9.
By Lemma 3.7 there is N such that for all N > N | z ( τ ) − α ( τ ) | ≤ ε/ . Next, by the averaging principle there exists c and β such that for all β > β we find constants M ( β ) and µ ( β ) such that | α ( τ ) − α I | ≤ | r ( τ ) − ρ ( τ ) | + | ρ ( τ ) | + | ρ ( τ ) − r − |≤ c β + M | ρ | e − µτ + M | ρ − r − | e − µτ . ONTENTS T we can make it below ε/
2. Applying the triangleinequality for t > T we obtain | z ( t ) − α I | ≤ | z ( t ) − α ( t ) | + | α ( t ) − α I | and the claim follows. (cid:3) The final step is to construct the invariant manifold relative to the asynchronousdynamics. To this end we use a persistence argument over the equation of α and ψ . Proposition 3.13.
Let σ strictly satisfy inequality (29) . Given ε > there exist N suchthat for each N ≥ N there exists a normally attracting invariant manifold with invariantboundary, M I , such that the order parameter z on M I satisfies | z − α I | ≤ ε, where α I is given by (30) .Proof. To construct the incoherent invariant manifold, we change to coordinates α, ψ using theWatanabe–Strogatz approach 3.4.1 while viewing the θ i ’s as fixed parameters. First considerthe modified system (20) with α as argument of the functions f , g . Then the equation for α decouples and we find a stable fixed point α I ∈ C , with | α I | < ψ is neutral, it follows that M WS = { α = α I , ψ ∈ S } ⊂ D × S isa NAIM. Lemma 3.7 shows that given a k ≥ B r ( α I ) ⊂ D , there exists an N andfor all N ≥ N a neighborhood U ⊂ T N of the uniformly distributed states Θ , such that thesubstitution of α by z ( α, ψ, θ ) is a C k -small perturbation of Eq. (20).We choose k = z also has a NAIM˜ M WS ( θ ) = { α = α I + g I ( ψ, θ ) , ψ ∈ S } ⊂ D × S (40) C -close to M WS . Thus, since U is precompact, we have uniformly for all θ ∈ U that (cid:107) g I (cid:107) C ≤ ε when N is su ffi ciently large, and thus on ˜ M WS ( θ ) also | z − α I | ≤ ε holds. Now we consider θ ∈ T N as dynamical variables again; that is, we lift the system according to (22) to coordinates( α, ψ, θ ). Since the θ i ’s have no dynamics, also M I = (cid:91) θ ∈ U ˜ M WS ( θ ) × { θ } ⊂ D × S × T N (41)is a NAIM, but with boundary (cid:107) . Since θ has no associated dynamics, the boundary is invariant. (cid:3) (cid:107) Here we mean the topological boundary, that is, ∂ M I = M I \ M I . ONTENTS
4. Statement and Proof of the Main Theorem
In this section, we state and prove the main result of this paper 4.2. First we must recapitulatethe results we obtained in Section 3.Consider two star networks with a large number of leaves N . The first star we shallassume to be in a synchronized, coherent state (labeled + and with coordinates labeled with + ) while the second network we assume to be in an incoherent state (labeled − and withcoordinates labeled with − ). Both the coherent and incoherent states correspond to normallyattracting invariant submanifolds of each star network: • The coherent invariant manifold is given in φ + i , ϕ + coordinates with i = , . . . , N by M + = { φ + = . . . = φ + N = φ C , ϕ + ∈ S } where φ C is fixed and given by (12). • The incoherent invariant manifold is given in ϕ − , α, ψ, θ ∈ Q − : = S × D × S × T N coordinates by M − = { α = α I + g I ( ψ, θ ) , ψ ∈ S , θ ∈ U ⊂ T N , ϕ − ∈ S } with α I given by (30) and g I a C -small function for N su ffi ciently large, followingProposition 3.13.Note that the product system has M = M + × M − as normally attracting invariant manifold.Moreover, the order parameters z + M : T N + → C and z − M : Q − → C are understood as therestriction on M and this is extended when we couple the two stars with intra-coupling ε guaranteeing the existence of chimera states as presented in Theorem 4.1. Theorem 4.1 (Normally Attracting Invariant Manifold for weakly coupled stars) . Considerthe coupled star network system (5) . For any h ∈ C and η > , there exist ε > , N > , β > and an open set I ⊂ R such that for any N > N , β > β and < ε ≤ ε , σ ∈ I thesystem has a normally attracting invariant manifold with boundary M ε that is O ( ε ) -close toM = M + × M − ⊂ T N + × Q − such thatr + M ε ( ϕ + , φ + , . . . , φ + N ) > − η and r − M ε ( α, ψ, θ, ϕ − ) < η. The choice of coupling function h between each star a ff ects the order of time thesechimera states exist. So, we split the results in Theorem 4.2 and Theorem 4.4. Theorem 4.2 (At least metastable Chimera States) . (i) [General Coupling] If h ∈ C in Theorem 4.1 is arbitrary, there exists an open setU ⊂ T N + of initial conditions θ ∈ U on M ε such that the coupled dynamics staysin M ε for a time of order /ε , uniformly for θ ∈ U. In other words, chimera states surelyexist for a time of order /ε ONTENTS (ii) [Kuramoto–Sakaguchi Coupling] In particular, if the coupling h ∈ C ∞ is of the formh ( φ , φ ) = c sin( φ − φ + ϑ ) . Then there exists an open set U ⊂ T N + of initial conditions θ ∈ U on M ε such that thecoupled dynamics stays in M ε for a time of order /ε , uniformly for θ ∈ U. If extra information on the coupling is given or the coupling structure is not fullysymmetric, then the chimeras can be asymptotically stable. This motivates us to highlightand discuss how general is our results:
Remark 4.3.
Let A ∈ R N + × N + be the adjacency matrix corresponding to the intercouplingbetween both stars. In a general flavor the networks dynamics we analyse is given as ˙ ϕ + = β + β σ N N (cid:88) j = sin( ϕ + j − ϕ + + δ ) + ε N (cid:88) j = A j h ( ϕ + − ϕ − j ) , ˙ ϕ + i = + σ sin( ϕ + − ϕ + i + δ ) + ε N (cid:88) j = A i j h ( ϕ + i − ϕ − j ) , i = , . . . , N , (42) and likewise for the “ − ” star. Both previous theorems, Theorem 4.1 and Theorem 4.2-(i), areextended to this general coupling structure between stars without any extra condition, sincethe persistence argument goes through. In particular, for the case of small number of connections between both stars, we stateTheorem 4.4.
Theorem 4.4 (Stable Chimera States) . Consider the coupled star network system (42) . Forany h ∈ C η > there exists ε > , N > , β > and an open set I ⊂ R such that forany β > β , N > N and ≤ k / N (cid:28) , where only k of the leaves are coupled, that is, A ii = for i ∈ I and A i j = otherwise, with | I | = k. Then for any < ε ≤ ε , σ ∈ I and initialconditions on M ε the system stays in a product of coherent and incoherent states for all time,even though it may leave M ε .4.1. NAIM for Coupled Stars When we couple the two stars we introduce a coupling of size ε between the two systems,with ε ≤ ε su ffi ciently small, then this perturbation of the product system will again have anormally attracting invariant manifold M ε close to M . This manifold can be described as M ε = (cid:110) φ + i = φ C + ˜ g Ci ( ϕ ± , ψ, θ ) , α = α I + ( g I + ˜ g I )( ϕ ± , ψ, θ ) , ϕ ± , ψ ∈ S , θ ∈ U ⊂ T N (cid:111) , (43)with functions (cid:107) ˜ g C (cid:107) C , (cid:107) ˜ g I (cid:107) C ∈ O ( ε ) describing the perturbation of M ε away from M . Notethat since the combined system is invariant under a global phase shift, the functions ˜ g C , ˜ g I only depend on ϕ ± through their phase di ff erence ϕ + − ϕ − . Finally, since M ε is a NAIM anda normally hyperbolic invariant manifold, it follows that M ε has an invariantly foliated stablemanifold W s ( M ε ), see [10, Thm. 4.1]. This implies that the orbits of points in a leaf W s ( m )shadow and exponentially in time converge to the orbit of m ∈ M ε . ONTENTS M ε has a boundary, both the persistence result and existence of the invariantstable foliation hold true on M with an arbitrarily small neighborhood of its boundaryremoved. That is, we modify the perturbed dynamics in a vertical neighborhood over ∂ U ,such that the modified coupled dynamics leaves a vertical section S = T N × D × T N + × ∂ U over ∂ M invariant. Hence persistence of M to M ε is well-defined, where M ε again has an invariantboundary with ∂ M ε ⊂ S , see Figure 7. This follows from Theorem 4.8 in [7, Sec. 4.2]. Since M ε is invariant, it again has an invariant stable foliation. As long as solutions of the coupledsystem on M ε stay away from ∂ M ε , they satisfy the unmodified dynamics and the conclusionsdrawn above remain true. This proves Theorem 4.1. M ε T N ∂ U ∂ U ∂ U δδ T N × D T N + Figure 7: Flow along the invariant manifold M ε . For point i we have to study how the coupling may perturb the trivial dynamics ˙ θ = U , where normal hyperbolicitysubsequentially may be lost. In view of (43), we can express the invariant manifold M ε for thecoupled system as a graph, where the coordinates φ + i and α depend on the other coordinates ϕ + , ϕ − , ψ, θ ∈ T N + × U , see Figure 7. Now the boundary of M ε sits over T N + × ∂ U . Solutioncurves that leave M ε do so by θ leaving U through its boundary ∂ U .In the uncoupled setup there was no dynamics along θ ∈ U ⊂ T N , that is, ˙ θ =
0, thus,with coupling we have ˙ θ ∈ O ( ε ). We choose open neighborhoods U , U such that Θ ⊂ U ⊂ U ⊂ U ⊂ T N with B ( U , δ ) ⊂ U and B ( U , δ ) ⊂ U for some fixed δ >
0, see also Figure 7. We may also assume that the modifications we madefor Theorem 4.1 are contained outside U . Solution curves with initial conditions in M ε and θ ∈ U will take O (cid:0) δε (cid:1) time to cross the gap U \ U . As the vector field is unmodified in thispart of phase space, these conclusions are true for the unmodified system. This proves point i . We continue to prove point ii . Consider Equation (5) and assume that the coupling issinusoidal, that is, h = c sin + c cos. Thus the equations for φ − become˙ φ − i = (1 − β ) − σ sin( φ − i − δ ) − β σ N N (cid:88) j = sin( φ − j + δ ) + ε h ( φ + i − φ − i + ϕ + − ϕ − ) . (44) ONTENTS h ( φ + i − φ − i + ϕ + − ϕ − ). Wedenote the associated vector field by X ∈ X ( T N ), but note that X also depends on φ + , ϕ ± asexternal variables. We can write this as h ( φ + i − φ − i + ϕ + − ϕ − ) = h ( − φ − i ) h ( φ + i + ϕ + − ϕ − ) , where h , h again denote sinusoidal functions. We cannot apply the Watanabe–Strogatzlift (26) to X , since h ( φ + i + ϕ + − ϕ − ) is not identical for all i . However, when we restrictto the invariant manifold M ε , then φ + i = φ C + ˜ g Ci ( ϕ + − ϕ − , ψ, θ ) with ˜ g Ci ∈ O ( ε ). Hence we find h ( φ + i − φ − i + ϕ + − ϕ − ) = h ( − φ − i ) h ( φ C + ϕ + − ϕ − ) + O ( ε ) . Now we can lift the vector field X (cid:48) associated to h ( φ − i ) h ( φ C + ϕ + − ϕ − ) using the Watanabe–Strogatz lift (26). This implies that ˆ X (cid:48) leaves θ constant, and since the remaining part ofthe coupling vector field is of order ε , it follows that ˙ θ ∈ O ( ε ). This proves claim ii ofTheorem 4.2. Finally, we prove Theorem 4.4. We now assume that the coupling matrices are of the form A ii = i ∈ I ⊂ N and zero otherwise, where | I | = k (cid:28) N . The coupling function h ∈ C is arbitrary; since h is defined on a compact set, (cid:107) h (cid:107) C is bounded.When only a small fraction of the leaves of the two stars are coupled, we can improve theprevious argument to obtain that the chimera state is fully stable for all time. To this end, weapply the two di ff erent lifts of the vector field for the incoherent star. We apply the Watanabe–Strogatz lift (26) to the original uncoupled vector field and the lift (25) to the coupling terms.Without coupling term, the dynamics of the unsynchronized star are given in Watanabe–Strogatz lifted coordinates by (20), (21) and ˙ θ =
0. The extra coupling term is given in φ coordinates by ˙ φ ± i = ε h ( φ ∓ i − φ ± i ± Γ ) for i ∈ I , (45)where Γ = ϕ + − ϕ − is the di ff erence of the hub phases. Let us denote by ε Y ∈ X ( T N ) thevector field describing the evolution of the φ − i variables according to (45). We shall use (25)to find the lift ˆ Y ∈ X ( D × S × T N ) of Y that gives equivalent dynamics in Watanabe–Strogatzcoordinates. First we calculate M (cid:48) α,ψ ( θ ) = ∂φ∂θ ( α, ψ, θ ) with φ = π ( α, ψ, θ ). From e i φ = M α,ψ (cid:0) w (cid:1) with w = e i θ we obtain ie i φ ∂φ∂θ = ∂∂ w (cid:34) α + e i ψ w + ¯ α e i ψ w (cid:35) ie i θ ⇐⇒ α + e i ψ w + ¯ α e i ψ w ∂φ∂θ = (cid:34) e i ψ + ¯ α e i ψ w − α + e i ψ w (1 + ¯ α e i ψ w ) ¯ α e i ψ (cid:35) w ⇐⇒ ( α + e i ψ w ) ∂φ∂θ = (1 + ¯ α e i ψ w ) e i ψ w − α ¯ α e i ψ w − ¯ α e i ψ w + ¯ α e i ψ w ⇐⇒ ∂φ∂θ ( α, ψ, θ ) = (1 − α ¯ α ) e i ( ψ + θ ) (1 + ¯ α e i ( ψ + θ ) )( α + e i ( ψ + θ ) ) = − | α | | α + e i ( ψ + θ ) | . ONTENTS Y is given by˙ θ i = ε (cid:32) ∂φ∂θ (cid:33) − · h (cid:0) φ − i − φ + i − Γ (cid:1) = ε | α + e i ( ψ + θ i ) | − | α | h (cid:0) π i ( α, ψ, θ i ) + ϕ − − ϕ + i (cid:1) for i ∈ I , (46)˙ θ i = i (cid:60) I , ˙ α = ψ =
0. Although the particular form of (46) should be useful formore detailed calculations using averaging theory, we did not pursue this. We shall simplyuse the fact that the θ i with i (cid:60) I are still conserved. We use a slightly modified version ofLemma 3.7 to show that if k (cid:28) N , then the result of the Lemma still holds true with the valuesof θ i with i ∈ I arbitrary, hence there exists an invariant manifold with invariant boundary forthe incoherent star network.For simplicity, let us assume by relabelling that I = { N − k + , . . . , N } . Define the set¯ Θ : = (cid:8) ( ϑ, ξ ) ∈ T N − k × T k | ∃ θ ∈ Θ : ∀ ≤ i ≤ N − k : ϑ i = θ i (cid:9) ⊂ T N . (47)Note that ¯ Θ is a simple extension of Θ with the first N − k entries still uniformly distributedaccording to Definition 3.6 and the last k arbitrary. The coupled dynamics leaves ¯ Θ invariantsince ˙ ϑ =
0, still. We now have the following modified version of Lemma 3.7; here we fix a C -norm, but this is all we needed anyways. Lemma 4.5.
Fix < r < and consider the disc D r = { α ∈ C | | α | < r } . For each ε > thereexists for all N su ffi ciently large and k (cid:28) N an open neighborhood U = U (cid:48) × T k of ¯ Θ ⊂ T N such that the map D r → C , α (cid:55)→ z ( α, ψ, θ ) − α , with z ( α, ψ, θ ) given by (21) , is smaller than ε in C k -norm, uniformly for all ( ψ, θ ) ∈ S × U.Proof.
The proof is a straightforward extension of the proof of Lemma 3.7. Given θ (cid:48) = ( ϑ, ξ ) ∈ ¯ Θ , let θ ∈ Θ as in (47). For z ( α, ψ, θ (cid:48) ) we can estimate the di ff erence with z ( α, ψ, θ ) as z ( α, ψ, θ (cid:48) ) − z ( α, ψ, θ ) = N N (cid:88) j = N − k + (cid:32) α + e i ψ e i ξ j + ¯ α e i ψ e i ξ j − α + e i ψ e i θ j + ¯ α e i ψ e i θ j (cid:33) . Since | α | < r is bounded away from 1, all terms are uniformly bounded in C -norm, say by B . Hence (cid:107) z ( α, ψ, θ (cid:48) ) − z ( α, ψ, θ ) (cid:107) C ≤ kN B , and thus for k (cid:28) N this yields a su ffi ciently smallcontribution to the overall C -norm of α (cid:55)→ z ( α, ψ, θ (cid:48) ) − α relative to the estimate for z ( α, ψ, θ )in the proof of 3.7. (cid:3) To prove Theorem 4.4, we now use Lemma 4.5 instead of Lemma 3.7. Note that U = U (cid:48) × T k is invariant under the coupled dynamics since its boundary ∂ U = ∂ U (cid:48) × T k is invariant. Thus, we find that also the persistent manifold M ε has an invariant boundary. Connecting two stars using a general coupling mechanism invalidates the use of theWatanabe–Strogatz approach and the description of the dynamics through the action of theM¨obius group. An important exception to this rule is when only the hubs of the stars are
ONTENTS ff erence of the hubs acts as an external field to theequations describing the order parameters z .Consider two stars with their hubs coupled sinusoidally, with strength ε ,˙ ϕ ± k = + σ sin( ϕ ± − ϕ ± k + δ ) , ˙ ϕ ± = β + σβ N N (cid:88) j = sin( ϕ ± j − ϕ ± H + δ ) + ε sin( ϕ ∓ − ϕ ∓ + δ ) , (48)Again we change to coordinates φ ± k = ϕ ± k − ϕ ± that measure the relative phases of the leaveswith respect to the hubs. Secondly, we introduce Γ = ϕ + − ϕ − , the di ff erence between the hubphases. With these new variables, the equations become˙ φ ± k = − β − σ sin( φ ± k − δ ) − βσ N N (cid:88) j = sin( φ ± j + δ ) + ε sin( ∓ Γ + δ ) , ˙ Γ = σβ N N (cid:88) j = (cid:0) sin( φ + k + δ ) − sin( φ − k + δ ) (cid:1) − ε cos( δ ) sin( Γ ) . For each star the equations of motion are in the appropriated form and we can apply the WSapproach, where the functions f ± and g ± are now given by f ± = − σ e − i δ i , g ± = − β − βσ (cid:32) z ± e i δ − ¯ z ± e − i δ i (cid:33) + ε sin( ∓ Γ + δ ) . The di ff erence here is that the functions g ± now depend on the phase di ff erence Γ of the hubs,that, in its turn, depends on the relative phases of the leaves of both populations. Nevertheless,the dependence is the same for all leaves and therefore follows the WS approach. Theequations for the WS variables α ± are thus given by˙ α ± = − σ e − i δ + β e i δ )( z ± ) + i (1 − β ) (cid:32) + ε sin( ∓ Γ + δ )1 − β (cid:33) z ± + σ β e − i δ | z ± | + σ e i δ , where the di ff erence now is that where before would appear only the natural frequency of theleaves ω =
1, now we have a term 1 + ε sin( ∓ Γ+ δ )1 − β that comes from the coupling. The equationfor Γ can be written as ˙ Γ = σβ (cid:61) (( z + − z − ) e − i δ ) − ε cos( δ ) sin( Γ )We can view the e ff ect of the coupling as inducing oscillations of amplitude ε/ (1 − β ) to thenatural frequency ω and therefore we expect the persistence of the fixed points. ONTENTS (cid:104) τ (cid:105) as function of inter-coupling strength ε . The (dots) correspond toaverage over fifty initial conditions chosen nearby M and (blue shaded) area is the standarddeviation. The (dashed) lines represent the estimated scaling ε − and ε − . Panel a) displaysthe case of coupling function does not depend on δ and panel b) include it. The parametersare β = σ = . η = .
5. Open question and Conclusions
We presented chimera states born from the coexistence of synchronous and asynchronousdynamics in the mean-field of a single star graph. We showed that the chimera states incoupled star graphs correspond to the existence of an invariant manifold with boundary. Weassociated their breaking with dynamics running into the boundary of the invariant manifoldand provided lower bounds for their survival. Predicting whether such chimeras break to acomplete synchronous dynamics of the network or to asynchronous remains an open question.The coupling strength between the stars provides the time scale the chimera exists – eitherof order 1 /ε (in the general inter-coupling function) or 1 /ε (for sinusoidal). We performedseveral numerical experiments to check these predictions. We will discuss this now. To thisend, we introduce a parameter η which splits the absolute value of the order parameter inTheorem (4.1) measured on the coherent and the incoherent star. The chimera lifetime τ isthe minimum time it takes to violate the splitting condition from Theorem (4.1), so: | z − − a I | > η or r + < − η. The lifetime depends on the initial conditions and parameters of the network. We selectinitial conditions starting in a neighborhood of M ¶ . And for each value of ε , we compute τ as shown in Figure 8. The behavior of (cid:104) τ (cid:105) as a function of ε does not entirely seems toplateau after two decades. In interesting question is whether this plateau really exist or isa numerical artefact. As we only provide an lower bound for the survival time, this issueremains open. To obtain description of such behavior we need to take a closer look at the ¶ We select initial conditions satisfying the restrictions for system be in M and we add a small random uniformvector with each coordinate drawn from the interval (0 , . ONTENTS
Acknowledgements
Jaap Eldering was supported by FAPESP grant 2015 / / \ R1 \ / References [1] D M Abrams, R Mirollo, S H Strogatz, and D A Wiley,
Solvable model for chimera states of coupledoscillators , Phys. Rev. Lett. (2008), 084103.[2] Gerrit Ansmann, E ffi ciently and easily integrating di ff erential equations with JiTCODE, JiTCDDE, andJiTCSDE , Chaos (2018), 043116.[3] P Ashwin and O Burylko, Weak chimeras in minimal networks of coupled phase oscillators , Chaos (2015), 013106.[4] A-L Barab´asi and R Albert, Emergence of scaling in random networks , Science (1999), 509–512.[5] C Bick and P Ashwin,
Chaotic weak chimeras and their persistence in coupled populations of phaseoscillators , Nonlinearity (2016), 1468.[6] C Chicone, Ordinary di ff erential equations with applications , 2nd ed., Texts in applied mathematics, v.34,Springer, New York, 2006.[7] J Eldering, Normally hyperbolic invariant manifolds — the noncompact case , Atlantis Series in DynamicalSystems, vol. 2, Atlantis Press, Paris, September 2013.[8] N Fenichel,
Persistence and smoothness of invariant manifolds for flows , Indiana Univ. Math. J. (1971 / Self-organized alternating chimera states in oscillatory media ,Scientific Reports (2015), no. 9883.[10] M W Hirsch, C C Pugh, and M Shub, Invariant manifolds , Lecture Notes in Mathematics, vol. 583,Springer-Verlag, Berlin, 1977.[11] F P Kemeth, S W Haugland, L Schmidt, I G Kevrekidis, and K Krischer,
A classification scheme forchimera states , Chaos (2016), 094815.[12] T-W Ko and G B Ermentrout, Bistability between synchrony and incoherence in limit-cycle oscillators withcoupling strength inhomogeneity , Phys. Rev. E (2008), 026210.[13] , Partially locked states in coupled oscillators due to inhomogeneous coupling , Phys. Rev. E (2008), 016203.[14] Y Kuramoto, Chemical Oscillations, Waves, and Turbulence , vol. 19, Springer Berlin Heidelberg, 1984.[15] Y Kuramoto and D Battogtokh,
Coexistence of coherence and incoherence in nonlocally coupled phaseoscillators. , Nonlinear Phenom. Complex Syst. (2002), no. 4, 380–385.[16] L Larger, B Penkovsky, and Y Maistrenko,
Laser chimeras as a paradigm for multistable patterns incomplex systems , Nature Communications (2015), no. 7752.[17] E A Martens, S Thutupalli, A Fourri`ere, and O Hallatschek, Chimera states in mechanical oscillatornetworks , Proc. Natl. Acad. Sci. USA (2013), 10563–10567.[18] S A Marvel, R E Mirollo, and S H Strogatz,
Identical phase oscillators with global sinusoidal couplingevolve by M¨obius group action , Chaos (2009), 043104, 11.[19] O E Omel’chenko, Coherence–incoherence patterns in a ring of non-locally coupled phase oscillators ,Nonlinearity (2013), 2469.[20] , The mathematics behind chimera states , Nonlinearity (2018), R121–R164.[21] M J Panaggio and D M Abrams, Chimera states: coexistence of coherence and incoherence in networks ofcoupled oscillators , Nonlinearity (2015), R67–R87. ONTENTS [22] J A Sanders, F Verhulst, and J Murdock, Averaging methods in nonlinear dynamical systems , AppliedMathematical Sciences, Springer New York, 2007.[23] R Toenjes, C Fiore, and T Pereira,
Network induced coherence resonance , to appear (2020).[24] J F Totz, J Rode, M R Tinsley, K Showalter, and H Engel,
Spiral wave chimera states in large populationsof coupled chemical oscillators , Nature Physics (2018), 282–285.[25] V Vlasov, A Pikovsky, and E E N Macau, Star-type oscillatory networks with generic kuramoto-typecoupling: A model for “japanese drums synchrony” , Chaos (2015), 123120.[26] V Vlasov, Y Zou, and T Pereira, Explosive synchronization is discontinuous , Phys. Rev. E (2015), no. 1,012904.[27] S Watanabe and S H Strogatz, Constants of motion for superconducting josephson arrays , Phys. D (1994), 197–253.[28] M Wolfrum and E Omel’chenko, Chimera states are chaotic transients , Physical Review E (2011),015201.[29] M Wolfrum, O E Omel’chenko, S Yanchuk, and Y L Maistrenko, Spectral properties of chimera states ,Chaos21