2-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center Manifold
Felix P. Kemeth, Bernold Fiedler, Sindre W. Haugland, Katharina Krischer
22-C
LUSTER F IXED -P OINT A NALYSIS OF M EAN -C OUPLED S TUART -L ANDAU O SCILLATORS IN THE C ENTER M ANIFOLD
A P
REPRINT
Felix P. Kemeth
Department of Chemical and Biomolecular EngineeringWhiting School of Engineering, Johns Hopkins University3400 North Charles Street, Baltimore, MD 21218, USA [email protected]
Bernold Fiedler
Institut für MathematikFreie Universität Berlin14195 Berlin, Germany
Sindre W. Haugland
Physik-Department, Nonequilibrium Chemical Physics,Technische Universität München,85748 Garching, Germany
Katharina Krischer
Physik-Department, Nonequilibrium Chemical Physics,Technische Universität München,85748 Garching, Germany A BSTRACT
We reduce the dynamics of an ensemble of mean-coupled Stuart-Landau oscillators close to thesynchronized solution. In particular, we map the system onto the center manifold of the Benjamin-Feirinstability, the bifurcation destabilizing the synchronized oscillation. Using symmetry arguments, wedescribe the structure of the dynamics on this center manifold up to cubic order, and derive expressionsfor its parameters. This allows us to investigate phenomena described by the Stuart-Landau ensemble,such as clustering and cluster singularities, in the lower-dimensional center manifold, providingfurther insights into the symmetry-broken dynamics of coupled oscillators. We show that clustersingularities in the Stuart-Landau ensemble correspond to vanishing quadratic terms in the centermanifold dynamics. In addition, they act as organizing centers for the saddle-node bifurcationscreating unbalanced cluster states as well for the transverse bifurcations altering the cluster stability.Furthermore, we show that bistability of different solutions with the same cluster-size distributioncan only occur when either cluster contains at least / of the oscillators, independent of the systemparameters. K eywords Globally coupled oscillators · Center manifold reduction · S N -equivariant systems Long-range interactions play a crucial role in various dynamical phenomena observed in nature. In a swarm of flashingfireflies, they may act as a synchronizing force, causing the swarm to flash in unison. Analogously, in an audienceclapping, the acoustic sound of the clapping can be recognized by each individual, leading to clapping in unison. Inthese cases, long-range interactions lead to the synchronization of individual units [1].On the other hand, long-range interactions may also lead to a split up of the individuals into two or more groups, alsocalled dynamical clustering. In electrochemistry, a stirred electrolyte or a common resistance may induce long-rangecoupling, leading to spatial clustering on the electrode [2, 3, 4, 5, 6, 7, 8]. In biology, this may explain the formation ofdifferent genotypes in an otherwise homogeneous environment [9, 10]. a r X i v : . [ n li n . C D ] O c t -Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center ManifoldThe individual units which experience this long-range or global coupling may be oscillatory, as in the case of flashingfireflies or in a clapping audience, or, as in the case of sympatric speciation, stationary genotypes. Here, we focus onthe former case of oscillatory units with long-range interactions.Clustering in oscillatory systems with long-range interactions has been subject to theoretical investigation for manyyears [11, 12, 13, 14, 15]. See also Ref. [16] for a recent review on globally coupled oscillators. In particular whenthe long-range interactions are weak compared to the intrinsic dynamics of the oscillator, it suffices to describe thephase evolution of each unit, and the analysis greatly simplifies [17, 18, 19]. If, however, the influence of the couplingis strong, as in the case considered here, such a reduction is no longer feasible and the amplitude dynamics must beconsidered. Our work aims to add to the theoretical understanding of clustering in this case of strong coupling.From the view-point of symmetry, if the coupling between N identical oscillators is global (i.e. all-to-all), then thegoverning equations are equivariant under the symmetric group S N . This means that the evolution equations f commutewith elements σ from the symmetry group, f ( σx ) = σf ( x ) ∀ σ ∈ S N . (1)In addition, this implies that the system has a trivial solution which is invariant under S N , that is, in which all oscillatorsare synchronized. Cluster states composed of two clusters, also called 2-cluster states, can then be viewed as stateswith the reduced symmetry S N × S N , with N and N being the number of oscillators in each cluster. Using theequivariant branching lemma, it can then be shown that these 2-cluster states bifurcate off the trivial solution [9, 20].The bifurcation at which the synchronized motion becomes unstable and the 2-cluster branches (also called primarybranches) emerge is commonly referred to as the Benjamin-Feir instability [21, 22].The intrinsic dimensionality of each oscillatory unit may range from d = 2 for FitzHugh-Nagumo [23] and Van der Poloscillators [24], via d = 3 for the Oregonator [25] to d = 4 for the original Hodgkin-Huxley model [26], and even higherfor more detailed physical models [27]. A system composed of N of these oscillators thus lives in a d · N -dimensionalphase space, making its full investigation unfeasible even for small d and N . One can, however, circumvent this problemof increasingly large dimensions by restricting the dynamics to the center manifold of certain bifurcations. In particular,it is known that the center space of the Benjamin-Feir instability is N − dimensional [20, 28], and thus a reduction tothe center manifold at this bifurcation allows for reducing the dimension of the problem to N − and thus by a factorof ≈ d . As we show below, such a reduction lets us reveal invariant sets and bifurcation curves analytically – a difficulttask in the original d · N -dimensional space.In this work, we focus on a particular example of a globally coupled system, in which the network is composed ofoscillating units called Stuart-Landau oscillators, each represented by a complex variable W k ∈ C . As opposed tophase oscillators, each Stuart-Landau oscillator has two degrees of freedom, i. e. an amplitude and a phase. With alinear global coupling, the dynamics are then given by ˙ W k = W k − (1 + iγ ) | W k | W k + ( β r + iβ i ) ( (cid:104) W (cid:105) − W k ) (2)with the complex coupling constant β r + iβ i and the real parameter γ , also called the shear [29]. (cid:104)·(cid:105) indicates theensemble mean and ˙ W = dW/dt . Bold face W indicates a vector containing the ensemble values [ W , W , . . . , W N ] .For (cid:107) β r + iβ i (cid:107) = 0 the ensemble is decoupled, and each Stuart-Landau oscillator oscillates with unit amplitudeand angular velocity − γ . For (cid:107) β r + iβ i (cid:107) (cid:54) = 0 , however, a plethora of different dynamical states can be observed.These states include fully synchronized oscillations, in which all oscillators maintain an amplitude equal to one andhave a mutual phase difference of zero [30], cluster states, in which the ensemble splits up into two or more sets ofsynchrony [31, 32, 14], and a variety of quasi-periodic and chaotic dynamics [33, 12].2-cluster states can be born and destroyed at saddle-node bifurcations if the number of oscillators in each cluster aredifferent, that is, when they are unbalanced [13]. Balanced solutions with N = N emerge from the synchronizedsolution at the Benjamin-Feir instability. For N = 16 oscillators and for γ = 2 , the saddle-node bifurcations fordifferent unbalanced cluster distributions N (cid:54) = N and the Benjamin-Feir instability are depicted in Fig. 1, as afunction of the coupling parameters β r and β i . Here, all the 2-cluster solutions exist locally in parameter space belowtheir respective saddle-node bifurcation curve, that is for smaller β r values. Up to the Benjamin-Feir instability theycoexist with the stable synchronized solution. Descending from large β r values, notice that the most-unbalanced clusterstate with N : N = 1 : 15 is created first. The more balanced cluster states are born subsequently, depending on theirdistribution, until eventually the balanced cluster state N : N = 8 : 8 is born at the Benjamin-Feir instability. At β r = − (1 − √ γ ) / , β i = ( − γ − √ / , there exists a codimension-two point where the saddle-node bifurcationsof all cluster distributions coincide. This point is called a cluster singularity [32]. Note that the qualitative picture inFig. 1 does not change when increasing the total number of oscillators N . For large numbers N → ∞ of infinitelymany oscillators we expect a bow-tie-shaped band of saddle-node bifurcation curves, ranging from the saddle-nodebifurcation of the most unbalanced cluster state to the Benjamin-Feir instability. As argued in Ref. [32], the cluster2-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center Manifoldsingularity can thus be viewed as an organizing center. By projecting the dynamics close to the Benjamin-Feir instabilityonto its center manifold, we aim to obtain further insights into the properties of this organizing center, and to elucidatethe clustering behavior near it. − . − . − . − . − . − . − . − . − . β i1 . . . . . . . β r Cluster Singularity γ = 2Benjamin-Feir7 : 96 : 105 : 114 : 123 : 132 : 141 : 15 Figure 1: The Benjamin-Feir instability involving the cluster (dark blue) and the different saddle-node curvescreating the unbalanced cluster solutions, N (cid:54) = N , in the β i , β r plane with γ = 2 and N = 16 . Each curve belongs toa particular cluster distribution N : N , and is obtained with numerical continuation using auto-07p [34, 35]. Note theposition of the cluster singularity at β r = − (1 − √ γ ) / ≈ . , β i = ( − γ − √ / ≈ − . as indicated.The remainder of this article is organized as follows: In Sec. 2, we pass to a corotating frame and introduce theaverage amplitude R , the deviations from the average amplitude r k , the deviations from the mean phase ϕ k . Usingthis corotating system, we discuss how one can describe the dynamics in the center manifold, see Sec. 3. In Sec. 4we derive the parameters for the dynamics of x k . Detailed calculations are provided in Appendix C, for convenience.Based on the parameters in the center manifold, we study the bifurcations of 2-cluster states and the role of the clustersingularity in the center manifold, in Sec. 5. We conclude with a detailed discussion of our results and an outlook offuture work. For a detailed mathematical analysis of the dynamics of 2-cluster states in the center manifold, see thecompanion paper Ref. [36]. Notice that Eq. (2) is invariant under a rotation in the complex plane W k → W k exp iφ . This invariance can beeliminated by choosing variables in a corotating frame, thus effectively reducing the dimensions of the system from N to N − . 3-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center ManifoldIn particular, we express the complex variables W k in log-polar coordinates W k = exp( R k + i Φ k ) . Then Eq. (2) turnsinto ˙ R = 1 − e R (cid:104) e r (cid:105) + Re (cid:0) ( β r + iβ i ) (cid:0) (cid:104) e z (cid:105)(cid:104) e − z (cid:105) − (cid:1)(cid:1) (3a) ˙ r k = − e R (cid:103) e r k + Re (cid:16) ( β r + iβ i ) (cid:16) (cid:104) e z (cid:105) (cid:103) e − z k (cid:17)(cid:17) (3b) ˙ ϕ k = − γe R (cid:103) e r k + Im (cid:16) ( β r + iβ i ) (cid:16) (cid:104) e z (cid:105) (cid:103) e − z k (cid:17)(cid:17) (3c)with k = 1 , . . . , N − , the abbreviations shown in Tab. 1 and the new coordinates summarized in Tab. 2 (see Appendix Afor a derivation). Hereby, (cid:101) · symbolizes the deviation from the ensemble mean (cid:104)·(cid:105) , and R and Φ are the ensemble meanTable 1: Abbreviations (cid:104) x m (cid:105) = 1 /N (cid:80) Nj =1 x mj (cid:102) x mk = x mk − (cid:104) x m (cid:105)(cid:104) e x (cid:105) = 1 /N (cid:80) Nj =1 e x j (cid:102) e x k = e x k − (cid:104) e x (cid:105) Table 2: Coordinate transformations R = (cid:104) R (cid:105) r k = (cid:102) R k ⇒ (cid:104) r (cid:105) = 0Φ = (cid:104) Φ (cid:105) ϕ k = (cid:102) Φ k ⇒ (cid:104) ϕ (cid:105) = 0 z k = r k + iϕ k ⇒ (cid:104) z (cid:105) = 0 logarithmic amplitude and phase, respectively. The logarithmic amplitude and phase deviation of each oscillator fromtheir averages are r k and ϕ k . Notice that through this construction, the averages of these deviations vanish. Furthermore,bold face of a variable, e.g. x , symbolizes the set of the respective ensemble variables { x , x , . . . , x N } .To simplify notation, r k + iϕ k is abbreviated by the complex variable z k . The transformation into Eqs. (3a) to (3c)has the advantage that the resulting equations are independent of the mean phase Φ . A change of Φ corresponds toa uniform phase shift of the whole ensemble in the complex plane, which in turn means that periodic orbits in theStuart-Landau ensemble, Eq. (2), correspond to stationary solutions in the transformed system, Eqs. (3a) to (3c). Thus,we can ignore the mean phase Φ in our subsequent analysis.Synchronized oscillations correspond to R = 0 , r k = 0 , Φ = − γt and ϕ k = 0 . The stability of this equilibrium can beinvestigated using the eigenspectrum of the Jacobian evaluated at this point. Due to the S N -symmetry of the solutionand the S N -equivariance of the governing equations, the Jacobian becomes block-diagonal, and thus has a degenerateeigenvalue spectrum [21, 15], see Appendix B:• There is one singleton eigenvalue λ = − < , corresponding to an eigendirection affecting all oscillatorsidentically. That is, this direction (cid:126)v shifts the amplitude of the synchronized motion but does not alter itssymmetry.• There is the eigenvalue λ + = − − β r + (cid:112) − β i − β i γ =: − − β r + d which becomes zero at theBenjamin-Feir instability and is of geometric multiplicity N − . The corresponding directions correspond to2-cluster states, with each direction corresponding to one cluster distribution N : N . Up to conjugacy, wearrange here the units such that the first N oscillators correspond to the same cluster. All 2-clusters with thesame distribution but different assignments of the oscillators then belong to the same conjugacy class.• Finally, there is the eigenvalue λ − = − − β r − d which is negative close to the synchronized solution, whichhas a geometric multiplicity of N − and whose eigendirections also have S N × S N -symmetry.Hereby d = (cid:112) − β i − β i γ abbreviates the root of the discriminant where we assume − β i − β i γ > , i.e. real λ ± .Notice that the Benjamin-Feir instability λ + = 0 , alias β r = d − , i.e. the dark blue curve in Fig. 1, is of codimensionone. 4-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center Manifold In the following, we calculate an expansion to third order of the ( N − -dimensional center manifold which correspondsto the Benjamin-Feir instability at λ + = 0 = − − β r + d . In order to do so, it is useful to introduce the coordinates x k = − r k + d +1 γ (cid:48) ϕ k d (4) y k = r k + d − γ (cid:48) ϕ k d (5)such that r k = (1 − d ) x k + (1 + d ) y k (6) ϕ k = γ (cid:48) x k + γ (cid:48) y k . (7)Here we use the notations γ (cid:48) = 2 γ + β i and d as defined above. See Appendix B for a derivation. The variables x k describe the dynamics in the ( N − -dimensional center manifold tangent to y k = 0 ∀ k , while y k together with R describe the dynamics in the stable manifold tangent to x k = 0 ∀ k .Note that the center-manifold must be S N -invariant. In addition, the global restrictions (cid:104) r (cid:105) = (cid:104) ϕ (cid:105) = 0 and thus (cid:104) x (cid:105) = (cid:104) y (cid:105) = 0 must hold. Therefore, the general form of the center manifold up to quadratic order must follow y k = y k ( x ) = a (cid:102) x k + O (cid:0) x k (cid:1) (8) R = R ( x ) = b (cid:104) x (cid:105) + O (cid:0) x k (cid:1) (9)with the coefficients a = a ( β i , γ ) and b = b ( β i , γ ) . Here, we use tangency of our coordinates R and y k , that is, dd x k R (cid:12)(cid:12)(cid:12) x =0 = 0 and dd x k y k (cid:12)(cid:12)(cid:12) x =0 = 0 . Since the Benjamin-Feir instability β r = d − is of codimension one, thethree-dimensional parameter space ( β r , β i , γ ) becomes two-dimensional. The parameters in the center manifold thusonly depend on β i and γ . By S N -equivariance, the reduced dynamics ˙ x k in the center manifold, up to cubic order, mustbe of the form ˙ x k = λ + x k + A (cid:102) x k + B (cid:102) x k + C (cid:104) x (cid:105) x k + O (cid:0) x k (cid:1) , (10)see also Refs. [10, 20], with the parameters A = A ( β i , γ ) and B = B ( β i , γ ) and C = C ( β i , γ ) . a , b , A , B and C In this section, we discuss the approach to calculate the coefficients a , b , A , B and C for the dynamics in the centermanifold. See Appendix C for complete details.First, we determine b . In particular we observe that ˙ R = (cid:18) dd x k R (cid:19) ˙ x = 2 b (cid:104) x ˙ x (cid:105) + O (cid:0) x k (cid:1) = 2 bλ + (cid:104) x (cid:105) + O (cid:0) x k (cid:1) holds. Since λ + = 0 at the bifurcation, ˙ R up to second order in x k must vanish. Therefore, expressing z k = r k + iϕ k and r k , ϕ k in terms of x k in Eq. (3a), we can compute b by comparing the coefficients of the (cid:104) x (cid:105) : the terms in front of (cid:104) x (cid:105) must thereby vanish. This allows us to estimate b = b ( β i , γ ) as b = 1 − d (cid:0) γ (cid:48) + d + 4 d − (cid:1) (11)with γ (cid:48) and d as defined above.Analogously, we can calculate a using Eqs. 3a and 3b up to second order in x k and employing ˙ y k = (cid:18) dd x k y k (cid:19) ˙ x k = O (cid:0) x k (cid:1) . This means we can use d ˙ y k = ˙ r k + ( d − /γ (cid:48) ˙ ϕ k , substitute the z k with x k in Eqs. 3a and 3b and keep terms up to O (cid:0) x k (cid:1) . Comparing the coefficients in front of (cid:102) x k then results in a = (1 − d ) (cid:16) γ (cid:48) + (1 − d ) (cid:17) (cid:0) (cid:0) d − (cid:1) + γ (cid:48) (cid:1) d γ (cid:48) . (12)5-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center ManifoldFinally, we can calculate A , B and C using d ˙ x k = − ˙ r k + d + 1 γ (cid:48) ˙ ϕ k = λ + x k + A (cid:102) x k + B (cid:102) x k + C (cid:104) x (cid:105) x k . Taking Eqs. 3b and 3c and the coefficients a and b obtained above, we can evaluate this equality up to cubic order,yielding the coefficients A = ( d − (cid:16) γ (cid:48) + (1 + d ) (cid:17) (cid:16) γ (cid:48) − d − (cid:17) γ (cid:48) d (13) B = − ( d − (cid:16) γ (cid:48) + ( d − (cid:17) (cid:16) γ (cid:48) + ( d + 1) (cid:17) (cid:0) γ (cid:48) − γ (cid:48) d + 3 (cid:0) d − (cid:1)(cid:1) (cid:0) γ (cid:48) + 2 γ (cid:48) d + 3 (cid:0) d − (cid:1)(cid:1) γ (cid:48) d (14) C = ( d − d γ (cid:48) (cid:18) γ (cid:48) − γ (cid:48) (cid:0) d − d + 1 (cid:1) − γ (cid:48) (cid:0) d + d − d + 22 d + 1 (cid:1) − γ (cid:48) (cid:0) d + 5 d − d − d + 2 d + 11 d − (cid:1) + 9 (cid:0) d − (cid:1) (cid:19) . (15)Together with λ + , the expressions for A , B and C fully specify the dynamics in the center manifold based on theoriginal parameters γ , β r and β i . By rescaling time and x k in Eq. (10), the number of independent parameters can bereduced to two, see Ref. [36]. For simplicity, we use the unscaled equation as in Eq. (10) here. As shown in Fig. 1 for N = 16 oscillators, we observe a range of saddle-node bifurcations creating the different2-cluster states. The expressions for λ + , A , B and C above determine the corresponding parameter values in thecenter manifold. The respective λ + and A values for the numerical curves shown in Fig. 1 are depicted in Fig. 2 asdashed curves. Notice that the Benjamin-Feir curve corresponds to the line λ + = 0 . Furthermore, we can derivethe saddle-node curves creating unbalanced 2-cluster states in the center manifold analytically, see Appendix D. Inparticular, λ sn = A (1 − α ) B (1 − α + α ) + Cα ) (16)for unbalanced cluster solutions, with α = N /N . The respective analytical curves for N = 16 are shown as solidcurves in Fig. 2. Notice the close correspondence between the mapped bifurcation curves from the full system and thebifurcation curves determined in the center manifold. For more balanced solutions, the saddle-node curves obtainedfrom the Stuart-Landau ensemble depart more strongly from the saddle-node curves calculated analytically in the centermanifold. We expect this to be due to the cubic truncation of the flow in the center manifold, thus limiting its accuracyaway from the Benjamin-Feir curve.Note that to obtain the curves in Fig. 2, we fix γ = 2 and vary β i , β r . We then use the expressions for A ( β i , β r ) , B ( β i , β r ) and C ( β i , β r ) to get the parameters in the center manifold. Thus the parameters A , B and C lie on a two-dimensionalmanifold. For the curves shown in Fig. 2, we furthermore use Eq. (16), yielding one-dimensional curves. The curves are,however, not exactly parabolas, since B and C vary in addition to A , which is not shown in Fig. 2. For all subsequentfigures, we use the values of C = − and B = − / (2 √ − at the cluster singularity for γ = 2 , which can beobtained analytically. See Ref. [36] p. 36 for a derivation.Furthermore, from Fig. 2 we observe that A = 0 , in addition to λ + = 0 , at the cluster singularity. This means that thiscodimension-two point is distinguished by vanishing quadratic dynamics in the center manifold, cf. Eq. 10. In addition,it serves as an organizing center for the saddle-node bifurcations of the unbalanced cluster states: At the saddle-nodebifurcation, we have in the center manifold for a cluster state x ∗ = − A (1 − α )2 ( B (1 − α + α ) + Cα ) , with B < and C < for the range of β i , β r considered here (not shown), see Appendix D. This means that fornegative A values, the saddle-node curves occur at positive x , for positive A values at negative x , and for A = 0 ,at the cluster singularity, all saddle-node bifurcations occur at the synchronized solution x k = 0 . This behavior canindeed be observed in the Stuart-Landau ensemble, see Fig. 6 of Ref. [32].6-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center Manifold − . − . − . . . . . A − . − . − . − . . . λ + γ = 2 Benjamin-Feir λ sn (7 : 9) λ sn (6 : 10) λ sn (5 : 11) λ sn (4 : 12) λ sn (3 : 13) λ sn (2 : 14) λ sn (1 : 15)7 : 96 : 105 : 114 : 123 : 132 : 141 : 15 Figure 2: The Benjamin-Feir instability (blue, λ + = 0 ) and the different saddle-node curves creating the unbalancedcluster solutions in the A , λ + plane. The dashed curves belong to particular cluster distribution N : N obtainedby projecting the curves from the Stuart-Landau ensemble shown in Fig. 1 using the expressions for λ + ( β r , β i , γ ) = − − β r + d , and A ( β r , β i , γ ) , cf. Eq. (13). The solid curves λ sn indicate the saddle-node bifurcations of the unbalancedcluster states obtained analytically in the center manifold, see Eq. (16). Note that close to the cluster singularity, whereanalytical expansions work best, numerical continuation fails due to the concentration of solutions in phase space.The unbalanced cluster states do, in general, not emerge as stable states from the saddle-node bifurcations. Rather, oneof the two branches created at the saddle-node bifurcation is subsequently stabilized through transverse bifurcationsinvolving 3-cluster solutions with symmetry S N × S N × S N , also called secondary branches [10]. For a moredetailed discussion on secondary branches, see also Refs. [37, 28]In order to explain this in more detail, we follow Ref. [37] Section 4. Note that each N : N S N × S N . From this, it follows that one can block-diagonalise the Jacobian atthe 2-cluster solutions S N × S N . In doing so, one can calculate the ( N − -degenerate eigenvalue µ describingthe intrinsic stability of cluster Ξ , that is its stability against transverse perturbations. Note, however, that a clusterof size 1 cannot be broken up. Following Ref. [10] p. 23 and using isotopic decomposition, the eigenvalue µ can beexpressed as µ = J | Ξ − J | Ξ . Here, J ij | Ξ denotes ∂f i /∂x j , with the respective x i and x j in cluster Ξ and f i being the right hand side of Eq. (10).Without loss of generality, we assume in the following that Ξ is the cluster with the smaller number of oscillators, thatis, N ≤ N or α ≤ . Evaluating the Jacobian, one obtains that the eigenvalue µ changes sign at λ + , = (1 − α ) B − αC ( α − B A . (17)7-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center ManifoldAnalogously, the transverse stability of cluster Ξ is described by µ = J | Ξ − J | Ξ which changes sign at λ + , = ( α − B − C (4 α − α + 1) B αA . (18)Hereby, µ describes the intrinsic stability of cluster Ξ . Furthermore notice that for the balanced cluster α = 1 andtherefore λ + , = λ + , . Since both clusters contain an equal number of units, their respective intrinsic stabilities changesimultaneously. − . − . − . − .
05 0 .
00 0 .
05 0 .
10 0 .
15 0 . A − . . . . . . λ + γ = 2 B = -2 / (2 √ − Benjamin-Feir λ + , | (8 : 8) λ sn (4 : 12) λ + , (4 : 12) λ + , (4 : 12) Figure 3: The bifurcation curves λ + , ( µ = 0 , dotted orange) and λ + , ( µ = 0 , dash-dotted orange) for the cluster state in the A , λ + plane and the parameters B = − / (2 √ − , C = − . The saddle-node curve creating the cluster is shown as solid orange curve. The Benjamin-Feir line is shown in blue, with the λ + , = λ + , curve forthe balanced cluster state depicted as dotted blue curve. The cluster is stable in the two regions between therespective λ + , = 0 and λ + , = 0 curve. The balanced cluster state is stable above the dotted blue curve.In Fig. 3, λ sn , λ + , and λ + , are shown as solid, dotted and dash-dotted orange curves, respectively, for the λ + = 0 , and the transverse bifurcation curve λ + , = λ + , , where the balanced cluster state is stabilized, is drawn as adotted blue curve. See Fig. 4 for the respective curves for a range of cluster distributions.Fig. 3 can be interpreted as follows: Coming from negative λ + values, the unbalanced cluster state is bornat λ sn (4 : 12) (solid orange). However, this 2-cluster state is unstable for the parameter values considered here: thecluster Ξ with units is intrinsically unstable with µ > and µ < . At the dotted orange curve, µ changes sign,rendering the cluster state stable. Subsequently, at the dash-dotted orange curve, µ changes sign, leaving thecluster Ξ with 12 units intrinsically unstable and thus the cluster unstable.8-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center ManifoldThe qualitatively same behavior can be observed for any cluster distribution, cf. Fig. 4, except for the most unbalancedstate ( ). There, cluster Ξ cannot be intrinsically unstable, since it contains only one unit. This means that thiscluster solution is born stable in its saddle-node bifurcation, and becomes unstable only at λ when µ = 0 . See alsothe bottom right plot in Fig. 4. In particular, λ sn = A / B for α = 0 , see Eq. (16), coincides with λ + , = A / B for α = 0 , cf. Eq. (17). Furthermore, it is worth noting that the stable patches in parameter space overlap for differentcluster distributions. This means that there is a multistability of different stable 2-cluster states.Notice that these results are in close correspondence with the behavior observed in the full Stuart-Landau ensemble,compare, for example, Fig. 3 with Figs. 4b and 5b in Ref. [32]. λ sn and λ + , are continuous functions of α . For N → ∞ , this means that there are continuous bands of bifurcationcurves: Going from λ sn ( α = 0) = A / B to λ sn ( α = 1) = 0 , there is band of saddle-node bifurcations creating theunbalanced cluster solutions. This band becomes infinitesimally thin at the cluster singularity A = 0 , giving it a bow-tielike shape. From λ + , ( α = 0) = A / B to λ + , ( α = 1) = ( − B − C ) A /B , the transverse bifurcations of thesmaller cluster stretch from the saddle-node curve of the most unbalanced cluster state to the transverse bifurcations ofthe balanced cluster state where λ + , is maximal, again yielding a bow-tie like shape in the A , λ + plane. Since λ + , has a pole at α = 1 / , the interpretation is a bit more involved. First, for the balanced cluster state α = 1 : λ + , ( α = 1) = ( − B − C ) A /B = λ + , ( α = 1) , and thus λ + , and λ + , coincide. For the most unbalanced solution α = 0 : λ + , ( α = 0) = 0 . This means the largercluster of the most unbalanced solution becomes unstable exactly when the balanced solution is born, that is, at theBenjamin-Feir instability λ + = 0 . For intermediate α values, however, the λ + , curve becomes steeper and infinitelysteep at α = 1 / , with the tip reaching to the cluster singularity. This can also be observed in Fig. 4, where the parabolabecomes thinner when going from the to the cluster states, and subsequently broadens again until the cluster. Altogether, the λ + , curves fill out the half plane λ + ≥ except the line A = 0 .These three bow-tie like regions of λ sn , λ + , and λ + , become infinitesimally thin and thus singular only at the clustersingularity λ + = 0 , A = 0 .The bifurcation scenario can be better visualized by plotting the λ sn , λ + , and λ + , as a function of the cluster size N /N , see Fig. 5. It depicts the λ + values of the saddle-node bifurcations creating the 2-cluster states ( λ sn , blue) andthe two transverse bifurcations (Eqs. (17) and (18)) altering the stability of the 2-clusters, with λ + , in green and λ + , in orange.When increasing λ + coming from negative values, all cluster states with N /N (cid:54) = 1 / are born in the saddle-nodebifurcation λ sn . Note that in fact two solutions for each N /N are created this way. In Fig. 5, one can observe that forthe most unbalanced state N /N → , the transverse bifurcation stabilizing the smaller cluster λ + , occurs immediatelyafter the saddle-node bifurcation creating that cluster. This bifurcation alters the stability of one of the two solutionsborn in the saddle-node bifurcation, and in particular renders the smaller of the two clusters in that solution stable totransverse perturbations. For the parameter regime considered here ( A = − . , B = − / (2 √ − and C = − ),this solution is in fact stabilized at this bifurcation, that is for λ + > λ + , .For N /N < / , the respective 2-cluster solution remains stable until λ + , , where the larger cluster becomes unstable,thus rendering the whole solution unstable. This can, for example, be observed for the cluster-size distribution,see Fig. 6(top). There, the variable of one cluster, x , is plotted as a function of the bifurcation parameter λ + . The bluedot on the left marks the saddle-node bifurcation wherein the two solutions are created. Initially, both solutionsare unstable. At λ + , (orange dot), one of them is stabilized and at λ + , (green dot), it is subsequently destabilized.For N /N > / , the scenario is different. There the solution that got stabilized at λ + , remains stable for all λ + > λ + , . The bifurcation λ + , instead occurs at the second cluster solution created at the saddle-node bifurcation.This is illustrated more clearly in Fig. 6(bottom) for the cluster solution. One of the two solutions becomes stableat λ + , , marked by an orange dot and as discussed above. Since N /N = 7 / > / , this solution remains stable forall λ + > λ + , . The second solution (upper curve in the bottom part of Fig. 6) first passes the synchronized solution atthe Benjamin-Feir bifurcation λ + = 0 and finally becomes stabilized at λ + , , marked by a green dot. λ + , diverges atthe pole N /N = 1 / , separating the two scenarios shown in Fig. 6. There the bifurcation switches from the solutionwith negative x (which, for λ + → ∞ , diverges to −∞ ) to the solution with positive x (which, for λ + → ∞ , divergesto + ∞ ).Notice how for the cluster distribution N /N = 7 / the two 2-cluster solutions are bistable for λ + > λ + , . Thatis, there exist two stable 2-cluster solutions with different x but the same cluster size ratio that are both stable.This, in fact, has also been observed in the Stuart-Landau ensemble, see for example Fig. 6 in Ref. [32]. Note that thesingularity of λ + , at N /N = 1 / ( α = 1 / ) is independent of the parameters A , B and C , see Eq. (18). This meansthat bistable solutions created as described above can in general only exist for N /N > / .9-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center Manifold . . . λ + λ sn (8 : 8) λ + , (8 : 8) λ + , (8 : 8) λ sn (7 : 9) λ + , (7 : 9) λ + , (7 : 9) . . . λ + λ sn (6 : 10) λ + , (6 : 10) λ + , (6 : 10) λ sn (5 : 11) λ + , (5 : 11) λ + , (5 : 11) A . . . λ + λ sn (4 : 12) λ + , (4 : 12) λ + , (4 : 12) A λ sn (3 : 13) λ + , (3 : 13) λ + , (3 : 13) − . − . . . . A . . . λ + λ sn (2 : 14) λ + , (2 : 14) λ + , (2 : 14) − . − . . . . A λ sn (1 : 15) λ + , (1 : 15) λ + , (1 : 15) B = -2 / (2 √ − Figure 4: The theoretical bifurcation curves λ + , ( µ = 0 , dotted) and λ + , ( µ = 0 , dash-dotted) for the differentcluster size distributions in the A , λ + plane and the parameters B = − / (2 √ − , C = − . The saddle-node curvescreating the unbalanced cluster solutions are represented as solid curves, which correspond to the shaded curves inFig. 2 with the same color coding. The Benjamin-Feir line is shown in blue. The unbalanced cluster states are stableabove the respective dotted curve and below the dash-dotted curve, except for the cluster, which is stable alreadyat the saddle-node bifurcation. For the cluster, the dotted and solid curves do not coincide but lie very close inparameter space. In this paper, we showed how one can map a system of globally coupled Stuart-Landau oscillators onto the ( N − -dimensional center manifold at the Benjamin-Feir instability. Thereby, we observed that the bifurcation curves at which2-cluster solutions are born closely resemble their counterparts in the original oscillatory system. This allowed us toinvestigate a codimension-two point called cluster singularity, from which all these bifurcation curves emanate. In thecenter manifold, we saw that this point corresponds to a vanishing coefficient A = 0 in front of the quadratic term ofthe equations of motion. Due to the reduced dynamics in this manifold, we were able to obtain stability boundariesfor 2-cluster states analytically. This allows for the more detailed investigation of the bow-tie-shaped cascade oftransverse bifurcations that govern the stability of these 2-cluster states, highlighting the role of the cluster singularityas an organizing center. The observed behavior is hereby independent of the oscillatory nature of each Stuart-Landauoscillator, but a result of the S N -equivariance of the full system. These findings may thus facilitate our understandingof this codimension-two point, and of clustering in general, even beyond oscillatory ensembles.Through this reduction to the center manifold, we could calculate the bifurcation curves creating the cluster solutions( λ sn ) and altering their stability ( λ + , and λ + , ) analytically. This allowed us to investigate when stable 2-clustersolutions exist more systematically, and in particular revealed when different solutions with the same cluster-sizedistribution are bistable (cf. Fig. 6). The relative cluster size N /N = 1 / seems to be a general lower limit for such a10-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center Manifold N /N − . . . . . . . λ + A = -0 .
2, B = -2 / (2 √ − Benjamin-Feir λ sn λ + , λ + , Figure 5: The bifurcation curves λ sn (blue), λ + , (orange) and λ + , (green) for the different cluster state distributionsin the λ + , N /N plane with A = − . and the parameters B = − / (2 √ − and C = − . The Benjamin-Feirinstability is indicated by the black solid line. The dashed magenta line indicates the location where λ + , diverges.The positions of the and cluster states are marked by the dotted vertical gray lines, see also Fig. 6 for therespective solution curves.bistable behavior. The bifurcation scenario of how states with different cluster size ratios N /N are created is therebydifferent from the Eckhaus instability [38] in reaction-diffusion systems. There, solutions of different wavelengths arecreated through supercritial pitchfork bifurcations at the trivial solution and subsequently stabilized through a sequenceof subcritical pitchfork bifurcations involving mixed-mode states. In our case, the different 2-cluster states are createdin saddle-node bifurcations and stabilized at λ + , at a single equivariant bifurcation point involving 3-cluster states.However, the detailed interaction between 2- and 3-cluster states still remains an open topic for future research.Note that the cubic truncation of the flow in the center manifold has a gradient structure [36]. This means that we canassign an abstract potential to each of the cluster distributions for a particular set of parameters λ + , A, B and C . Isthere a particular cluster distribution with a minimal potential value? What is its role in the dynamics between thesecluster distributions? The companion paper [36] addresses some of these dynamical questions.Here, we fixed the parameter γ = 2 in the full Stuart-Landau system, and varied the coupling parameters β r , β i . Thisrestricts our analysis to a small region in parameter space. It is important to mention that for different parameter regimes,a qualitatively different behavior close to the cluster singularity might be observed [36].As discussed in Sec. 2, the Stuart-Landau ensemble permits the transformation into a corotating frame. This turnslimit-cycle dynamics into fixed-point dynamics and thus greatly facilitates the reduction onto the center manifold.For more general oscillatory ensembles, such as systems composed of van der Pol or Hogdkin-Huxley type units, thetransformation to a corotating frame may be more cumbersome or not even possible. If the coupling between such unitsis of a global nature, we expect, however, that the nesting of bifurcation curves creating different cluster distributions,cf. Fig. 2, can also be observed in these systems. 11-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center Manifold − . − . − . . . . x synchronized solution4 : 12 λ sn λ + , λ + , .
00 0 .
01 0 .
02 0 .
03 0 .
04 0 .
05 0 . λ + − . − . . . . x synchronized solution7 : 9 λ sn λ + , λ + , Figure 6: The variable x of the cluster solution (top) and the cluster solution (bottom) as a function of thebifurcation parameter λ + with A = − . , and the parameters B = − / (2 √ − and C = − . Solid curves indicatethat the solution is stable for the respective range of parameters, dashed curves represent unstable solutions. The pointsmark the λ sn (blue), the λ + , (orange) and the λ + , (green) bifurcations. The synchronized solution x i = 0 ∀ i isindicated by the black horizontal line. See also Fig. 5 for the locations of the and cluster in the λ + , N /N planeThis directly links to the fact that we focused on oscillatory dynamics in this article. An exciting further question isthe possibility of equivalent dynamics, such as clustering and cluster singularities, in systems composed of bistable orexcitable units. Acknowledgement
FPK thanks BF for the hospitality and the exciting discussions at the Freie Universität Berlin. BF gratefully acknowl-edges the deep inspiration by, and hospitality of, his coauthors at München who initiated this work. This work has alsobeen supported by the Deutsche Forschungsgemeinschaft, SFB910, project A4 “Spatio-Temporal Patterns: Control,Delays, and Design”, and by KR1189/18 “Chimera States and Beyond”.
A Variable transformation
Using log-polar coordinates W k = exp ( R k + i Φ k ) , Eq. (2) turns into (cid:16) ˙ R k + i ˙Φ k (cid:17) e R k + i Φ k = e R k + i Φ k − (1 + iγ ) e R k e R k + i Φ k + ( β r + iβ i ) (cid:0) (cid:104) e R + i Φ (cid:105) − e R k + i Φ k (cid:1) . W k this becomes ˙ R k + i ˙Φ k = 1 − (1 + iγ ) e R k + ( β r + iβ i ) (cid:0) (cid:104) e R + i Φ (cid:105) e − R k − i Φ k − (cid:1) . We average over k and separate real and imaginary parts. The mean amplitude R and the mean phase Φ then satisfy ˙ R = 1 − (cid:104) e R (cid:105) + Re (cid:0) ( β r + iβ i ) (cid:0) (cid:104) e R + i Φ (cid:105)(cid:104) e − R − i Φ (cid:105) − (cid:1)(cid:1) ˙Φ = − γ (cid:104) e R (cid:105) + Im (cid:0) ( β r + iβ i ) (cid:0) (cid:104) e R + i Φ (cid:105)(cid:104) e − R − i Φ (cid:105) − (cid:1)(cid:1) . Substituting the variables listed in Tab. 2, one obtains (cid:104) exp (2 R ) (cid:105) = (cid:104) exp (2 r + 2 R ) (cid:105) = exp (2 R ) (cid:104) exp (2 r ) (cid:105) , and (cid:104) exp ( R + i Φ ) (cid:105) = (cid:104) exp ( r + R + i ϕ + i Φ) (cid:105) = exp ( R + i Φ) (cid:104) exp z (cid:105) . Therefore ˙ R = 1 − e R (cid:104) e r (cid:105) + Re (cid:0) ( β r + iβ i ) (cid:0) (cid:104) e z (cid:105)(cid:104) e − z (cid:105) − (cid:1)(cid:1) ˙Φ = − γe R (cid:104) e r (cid:105) + Im (cid:0) ( β r + iβ i ) (cid:0) (cid:104) e z (cid:105)(cid:104) e − z (cid:105) − (cid:1)(cid:1) . For the deviations r k = R k − R and ϕ k = Φ k − Φ one may write ˙ r k = ˙ R k − ˙ R = 1 − e R e r k + Re (cid:0) ( β r + iβ i ) (cid:0) (cid:104) e z (cid:105) e − z k − (cid:1)(cid:1) − ˙ R = − e R (cid:103) e r k + Re (cid:16) ( β r + iβ i ) (cid:16) (cid:104) e z (cid:105) (cid:103) e − z k (cid:17)(cid:17) ˙ ϕ k = ˙Φ k − ˙Φ= − γe R e r k + Im (cid:0) ( β r + iβ i ) (cid:0) (cid:104) e z (cid:105) e − z k − (cid:1)(cid:1) − ˙Φ= − γe R (cid:103) e r k + Im (cid:16) ( β r + iβ i ) (cid:16) (cid:104) e z (cid:105) (cid:103) e − z k (cid:17)(cid:17) The equations for ˙ R , ˙ r k and ˙ ϕ k then constitute the corotating system Eqs. (3a) to (3c). B Linearization
Linearizing the dynamics of the transformed system, Eqs. (3a) to (3c), at the equilibrium R = 0 , r k = ϕ k = 0 , z k = 0 ,and using the fact that (cid:104) r (cid:105) = 0 , (cid:104) z (cid:105) = 0 , see Tab. 2, one gets ˙ R ˙ r k ˙ ϕ k = (cid:32) − R − r k − Re ( kz k ) − γr k − Im ( kz k ) (cid:33) = (cid:32) − R − (2 + β r ) r k + β i ϕ k − (2 γ + β i ) r k − β r ϕ k (cid:33) = (cid:32) − − − β r β i − γ − β i − β r (cid:33) · (cid:32) Rr k ϕ k (cid:33) = J · (cid:32) Rr k ϕ k (cid:33) . The Jacobian thus has the eigenvalues• Eigenvalue λ = − with eigenvector (cid:126)v = (cid:16) ,(cid:126) ,(cid:126) (cid:17) .and two eigenvalues of geometric multiplicity N − given by the eigendecompositioneig − − β r β i − γ − β i − β r , which gives• the eigenvalue λ + = − − β r + (cid:112) − β i − β i γ = − − β r + d • and the eigenvalue λ − = − − β r − (cid:112) − β i − β i γ = − − β r − d .13-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center ManifoldHere, we assume − β i − β i γ > , that is real λ ± . For an analysis of the case − β i − β i γ < , see Ref. [39]. Theeigenvectors corresponding to these two eigenvalues can be obtained using − − β r β i − γ − β i − β r − λ ± ( N − × ( N − (cid:126)v ± = (cid:126) . For λ + , one thus obtains (cid:18) − − d β i − γ − β i − d (cid:19) (cid:126)v + = (cid:18) − − d β i − γ − β i − d (cid:19) (cid:18) r k ϕ k (cid:19) = (cid:18) ( − − d ) r k + β i ϕ k ( − γ − β i ) r k + (1 − d ) ϕ k (cid:19) = 0 . Choosing ϕ k = (1 + d ) /β i r k , (A.1)we get (cid:18) ( − − d ) r k + (1 + d ) r k ( − γ − β i ) r k + (cid:0) − d (cid:1) /β i r k (cid:19) = (cid:18) − (1 + d ) r k + (1 + d ) r k − (2 γ + β i ) r k + (2 γ + β i ) r k (cid:19) = (cid:126) , thus solving the equality above. The constraint Eq. (A.1), together with (cid:104) r (cid:105) = (cid:104) ϕ (cid:105) = 0 , defines an ( N − -dimensionalsubspace of R N − .For λ − , one thus obtains (cid:18) − d β i − γ − β i d (cid:19) (cid:126)v + = (cid:18) − d β i − γ − β i d (cid:19) (cid:18) r k ϕ k (cid:19) = (cid:18) ( − d ) r k β i ϕ k ( − γ − β i ) r k (1 + d ) ϕ k (cid:19) . Choosing ϕ k = (1 − d ) /β i r k , (A.2)solves the conditions above. In particular, (cid:18) ( − d ) r k + (1 − d ) r k ( − γ − β i ) r k + (cid:0) − d (cid:1) /β i r k (cid:19) = (cid:18) − (1 − d ) r k + (1 − d ) r k − (2 γ + β i ) r k + (2 γ + β i ) r k (cid:19) = (cid:126) . The constraint Eq. (A.2), together with (cid:104) r (cid:105) = (cid:104) ϕ (cid:105) = 0 define an ( N − -dimensional subspace of R N − .Now, one can define the eigencoordinates x k describing the dynamics in the space defined by the constraint Eq. (A.1),the center space of the bifurcation, and eigencoordinates y k , describing the dynamics in the space defined by theconstraint Eq. (A.2). These two sets of variables, together with R , can then be used to describe the full system. C Parameter Derivation
In this section of the appendix, we derive expressions for the parameters a , b , A , B and C as a function of the parameters γ , β r and β i from the Stuart-Landau ensemble. Hereby, we will use the condition that R and the y k are tangential, thatis, dd x k R (cid:12)(cid:12)(cid:12) x =0 = 0 and dd x k y k (cid:12)(cid:12)(cid:12) x =0 = 0 . C.1 a and b In order to calculate a and b , it is useful to write out the following expressions z k = r k + iϕ k = (1 − d ) x k + (1 + d ) y k + i ( γ (cid:48) x k + γ (cid:48) y k )= (1 − d + iγ (cid:48) ) x k + a (1 + d + iγ (cid:48) ) (cid:102) x k + O (cid:0) x k (cid:1) z k = ( r k + iϕ k ) = ((1 − d ) x k + (1 + d ) y k + i ( γ (cid:48) x k + γ (cid:48) y k )) = ((1 − d + iγ (cid:48) ) x k + (1 + d + iγ (cid:48) ) y k ) = (1 − d + iγ (cid:48) ) x k + 2 a (1 − d + iγ (cid:48) ) (1 + d + iγ (cid:48) ) x k (cid:102) x k + O (cid:0) x k (cid:1) z k = ( r k + iϕ k ) = ((1 − d + iγ (cid:48) ) x k + (1 + d + iγ (cid:48) ) y k ) = (1 − d + iγ (cid:48) ) x k + O (cid:0) x k (cid:1) y k and the notation γ (cid:48) = 2 γ + β i . Similarly, we expand the following parts and keep terms upto cubic order: e z k = 1 + z k + z k z k O (cid:0) x k (cid:1) e − z k = 1 − z k + z k − z k O (cid:0) x k (cid:1) (cid:104) e z (cid:105) = (cid:104) z + z z O (cid:0) x k (cid:1) (cid:105) = 1 + 12 (cid:104) z (cid:105) + 16 (cid:104) z (cid:105) + O (cid:0) x k (cid:1)(cid:103) e − z k = e − z k − (cid:104) e − z (cid:105) = 1 − z k + z k − z k − − (cid:104) z (cid:105) + 16 (cid:104) z (cid:105) + O (cid:0) x k (cid:1) = − z k + 12 (cid:101) z k − (cid:101) z k + O (cid:0) x k (cid:1) (cid:104) e z (cid:105)(cid:104) e − z (cid:105) = (cid:18) (cid:104) z (cid:105) + 16 (cid:104) z (cid:105) (cid:19) (cid:18) (cid:104) z (cid:105) − (cid:104) z (cid:105) (cid:19) + O (cid:0) x k (cid:1) = 1 + 12 (cid:104) z (cid:105) + 16 (cid:104) z (cid:105) + 12 (cid:104) z (cid:105) − (cid:104) z (cid:105) + O (cid:0) x k (cid:1) = 1 + (cid:104) z (cid:105) + O (cid:0) x k (cid:1) (cid:104) e z (cid:105) (cid:103) e − z k = (cid:18) (cid:104) z (cid:105) + 16 (cid:104) z (cid:105) (cid:19) (cid:18) − z k + 12 (cid:101) z k − (cid:101) z k (cid:19) + O (cid:0) x k (cid:1) = (cid:18) (cid:104) z (cid:105) (cid:19) (cid:18) − z k + 12 (cid:101) z k − (cid:101) z k (cid:19) + O (cid:0) x k (cid:1) = − z k + 12 (cid:101) z k − (cid:101) z k − z k (cid:104) z (cid:105) + O (cid:0) x k (cid:1) . With the expression for R , see Eq. (9), we can furthermore write e R = 1 + 2 R + O (cid:0) x k (cid:1) = 1 + 2 b (cid:104) x (cid:105) + O (cid:0) x k (cid:1) e r k = 1 + 2 r k + 2 r k + 43 r k + O (cid:0) x k (cid:1) (cid:104) e r (cid:105) = 1 + 2 (cid:104) r (cid:105) + 43 (cid:104) r (cid:105) + O (cid:0) x k (cid:1)(cid:103) e r k = e r k − (cid:104) e r (cid:105) = 2 r k + 2 (cid:101) r k + 43 (cid:101) r k + O (cid:0) x k (cid:1) e R (cid:104) e r (cid:105) = (cid:0) b (cid:104) x (cid:105) (cid:1) (cid:18) (cid:104) r (cid:105) + 43 (cid:104) r (cid:105) (cid:19) + O (cid:0) x k (cid:1) = 1 + 2 (cid:104) r (cid:105) + 2 b (cid:104) x (cid:105) + 43 (cid:104) r (cid:105) + O (cid:0) x k (cid:1) e R (cid:103) e r k = (cid:0) b (cid:104) x (cid:105) (cid:1) (cid:18) r k + 2 (cid:101) r k + 43 (cid:101) r k (cid:19) + O (cid:0) x k (cid:1) = 2 r k + 2 br k (cid:104) x (cid:105) + 2 (cid:101) r k + 43 (cid:101) r k + O (cid:0) x k (cid:1) . R up to second order in x k ˙ R = 1 − e R (cid:104) e r (cid:105) + Re (cid:0) ( β r + iβ i ) (cid:0) (cid:104) e z (cid:105)(cid:104) e − z (cid:105) − (cid:1)(cid:1) = 1 − (cid:0) (cid:104) r (cid:105) + 2 b (cid:104) x (cid:105) (cid:1) + Re (cid:0) ( β r + iβ i ) (cid:0) (cid:104) z (cid:105) − (cid:1)(cid:1) = − (cid:104) r (cid:105) − b (cid:104) x (cid:105) + Re (cid:0) ( β r + iβ i ) (cid:104) z (cid:105) (cid:1) = − − d ) (cid:104) x (cid:105) − b (cid:104) x (cid:105) + Re (cid:16) ( β r + iβ i ) (1 − d + iγ (cid:48) ) (cid:17) (cid:104) x (cid:105) = − − d ) (cid:104) x (cid:105) − b (cid:104) x (cid:105) + (cid:16) β r (cid:16) (1 − d ) − γ (cid:48) (cid:17) − β i ( γ (cid:48) (1 − d )) (cid:17) (cid:104) x (cid:105) = − β r (cid:104) x (cid:105) − b (cid:104) x (cid:105) + (cid:0) β r (cid:0) β r − γ (cid:48) (cid:1) − (cid:0) β r + 2 β r (cid:1) β r (cid:1) (cid:104) x (cid:105) = − (cid:0) β r − β r (cid:0) β r − γ (cid:48) (cid:1) + 2 (cid:0) β r + 2 β r (cid:1) β r − b (cid:1) (cid:104) x (cid:105) = − (cid:0) β r + β r + β r γ (cid:48) +2 b (cid:1) (cid:104) x (cid:105) Now, we use the tangential property of R . In particular, we can write ˙ R = (cid:18) dd x k R (cid:19) ˙ x = 2 b (cid:104) x ˙ x (cid:105) + O (cid:0) x k (cid:1) = 2 bλ + (cid:104) x (cid:105) + O (cid:0) x k (cid:1) . At λ + = 0 , ˙ R up to second order must vanish. This allows us to calculate b by comparing the terms in front of (cid:104) x (cid:105) in ˙ R , yielding ⇒ b = − β r (cid:0) γ (cid:48) + 6 β r + β r (cid:1) = 1 − d (cid:0) γ (cid:48) + d + 4 d − (cid:1) . We can derive the expression for a in a similar way. Here, we write out the dynamics of y k up to second order. Thisyields d ˙ y k = ˙ r k + d − γ (cid:48) ˙ ϕ k = − (cid:18) d − γγ (cid:48) (cid:19) e R (cid:103) e r k + Re (cid:18)(cid:18) − i d − γ (cid:48) (cid:19) ( β r + iβ i ) (cid:16) (cid:104) e z (cid:105) (cid:103) e − z k (cid:17)(cid:19) = − (cid:18) d − γγ (cid:48) (cid:19) (cid:16) r k + 2 (cid:101) r k (cid:17) + Re (cid:18)(cid:18) − i d − γ (cid:48) (cid:19) ( β r + iβ i ) (cid:18) − z k + 12 (cid:101) z k (cid:19)(cid:19) = − (cid:18) d − γγ (cid:48) (cid:19) (cid:16) r k + 2 (cid:101) r k (cid:17) + Re (cid:18)(cid:18) − i d − γ (cid:48) (cid:19) ( β r + iβ i ) (cid:18) − r k − iϕ k + 12 (cid:101) z k (cid:19)(cid:19) . The term of the coupling constant and its parameters in front can be summarized by (cid:18) − i β r γ (cid:48) (cid:19) ( β r + iβ i ) = β r + β i β r γ (cid:48) − i (cid:18) β r γ (cid:48) − β i (cid:19) = β r − β r + 2 β r γ (cid:48) − i (cid:18) β r γ (cid:48) + β r + 2 β r γ (cid:48) (cid:19) β r γγ (cid:48) = β r γ (cid:48) − β i γ (cid:48) = β r β r + 2 β r γ (cid:48) . This simplifies the expression for ˙ y k to d ˙ y k = − (cid:18) β r + β r + 2 β r γ (cid:48) (cid:19) (cid:16) r k + (cid:101) r k (cid:17) + Re (cid:18)(cid:18) β r − β r + 2 β r γ (cid:48) − i (cid:18) β r γ (cid:48) + β r + 2 β r γ (cid:48) (cid:19)(cid:19) (cid:18) − r k − iϕ k + 12 (cid:101) z k (cid:19)(cid:19) = − (cid:18) β r + β r + 2 β r γ (cid:48) (cid:19) (cid:16) r k + (cid:101) r k (cid:17) + Re (cid:18)(cid:18) β r − β r + 2 β r γ (cid:48) − i (cid:18) β r γ (cid:48) + β r + 2 β r γ (cid:48) (cid:19)(cid:19) (cid:18) − r k − iϕ k + 12 (cid:16) (1 − d ) − γ (cid:48) (cid:17) (cid:102) x k + i (1 − d ) γ (cid:48) (cid:102) x k (cid:19)(cid:19) = − (cid:18) β r + β r + 2 β r γ (cid:48) (cid:19) (cid:16) r k + (cid:101) r k (cid:17) + (cid:18) β r − β r + 2 β r γ (cid:48) (cid:19) (cid:18) − r k + 12 (cid:16) (1 − d ) − γ (cid:48) (cid:17) (cid:102) x k (cid:19) − (cid:18) β r γ (cid:48) + β r + 2 β r γ (cid:48) (cid:19) (cid:16) ϕ k − (1 − d ) γ (cid:48) (cid:102) x k (cid:17) = − β r + 1) r k − (cid:18) β r + β r + 2 β r γ (cid:48) (cid:19) (cid:101) r k + 12 (cid:18) β r − β r + 2 β r γ (cid:48) (cid:19) (cid:0) β r − γ (cid:48) (cid:1) (cid:102) x k − (cid:0) β r + β r (cid:1) ( x k + y k ) − (cid:0) β r + β r (cid:1) β r (cid:102) x k = − β r + 1) ( − β r x k + ( β r + 2) y k ) − (cid:18) β r + β r + 2 β r γ (cid:48) (cid:19) β r (cid:102) x k + 12 (cid:18) β r − β r + 2 β r γ (cid:48) (cid:19) (cid:0) β r − γ (cid:48) (cid:1) (cid:102) x k − (cid:0) β r + β r (cid:1) ( x k + y k ) − (cid:0) β r + β r (cid:1) β r (cid:102) x k = − β r + 1) y k − (cid:18) β r + β r + 2 β r γ (cid:48) (cid:19) β r (cid:102) x k + 12 (cid:18) β r − β r + 2 β r γ (cid:48) (cid:19) (cid:0) β r − γ (cid:48) (cid:1) (cid:102) x k − (cid:0) β r + β r (cid:1) β r (cid:102) x k = − β r + 1) y k − (cid:18) β r + β r + 2 β r γ (cid:48) (cid:19) β r (cid:102) x k + 12 (cid:18) β r − β r + 2 β r γ (cid:48) (cid:19) (cid:0) β r − γ (cid:48) (cid:1) (cid:102) x k = − β r + 1) y k − (cid:18) β r + 3 β r + 6 β r γ (cid:48) (cid:19) β r (cid:102) x k − (cid:0) β r γ (cid:48) − β r − β r (cid:1) (cid:102) x k = − β r + 1) a (cid:102) x k − (cid:18) β r + 3 β r + 6 β r γ (cid:48) (cid:19) β r (cid:102) x k − β r γ (cid:48) (cid:102) x k = − β r + 1) a (cid:102) x k − β r γ (cid:48) (cid:0) γ (cid:48) + 6 β r γ (cid:48) + 4 β r γ (cid:48) + 3 β r + 6 β r (cid:1) (cid:102) x k = − β r + 1) a (cid:102) x k − β r γ (cid:48) (cid:0) γ (cid:48) + β r (cid:1) (cid:0) β r ( β r + 2) + γ (cid:48) (cid:1) (cid:102) x k Similar to R , the y k are tangential to the center manifold. This translates into the fact that ˙ y k = (cid:18) dd x k y k (cid:19) ˙ x k vanishes up to second order in x k . Therefore, comparing the terms in front of the (cid:102) x k above yields a = − β r (cid:0) γ (cid:48) + β r (cid:1) (cid:0) β r ( β r + 2) + γ (cid:48) (cid:1) β r + 1) γ (cid:48) = (1 − d ) (cid:16) γ (cid:48) + (1 − d ) (cid:17) (cid:0) (cid:0) d − (cid:1) + γ (cid:48) (cid:1) d γ (cid:48) . C.2 A , B and C Finally, the coefficients A , B and C for the dynamics in the center manifold, cf. Eq. (10), can be obtained by expandingthe dynamics of x k , d ˙ x k = − (cid:18) − d + 1) γγ (cid:48) (cid:19) e R (cid:103) e r k + Re (cid:18)(cid:18) − − i d + 1 γ (cid:48) (cid:19) k (cid:16) (cid:104) e z (cid:105) (cid:103) e − z k (cid:17)(cid:19) , (C.2)17-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center Manifoldin powers of x k : The terms in front of (cid:102) x k , (cid:102) x k and x k (cid:104) x (cid:105) correspond to the coefficients A , B and C , respectively. Inorder to do so, we approximate several terms as follows: (cid:104) e z (cid:105) (cid:103) e − z k = − z k + 12 (cid:101) z k − (cid:101) z k − z k (cid:104) z (cid:105) + O (cid:0) x k (cid:1) z k = (1 − d + iγ (cid:48) ) x k + a (1 + d + iγ (cid:48) ) (cid:102) x k + O (cid:0) x k (cid:1) z k = (1 − d + iγ (cid:48) ) x k + 2 a (1 − d + iγ (cid:48) ) (1 + d + iγ (cid:48) ) x k (cid:102) x k + O (cid:0) x k (cid:1)(cid:101) z k = z k − (cid:104) z (cid:105) = (1 − d + iγ (cid:48) ) (cid:102) x k + 2 a (1 − d + iγ (cid:48) ) (1 + d + iγ (cid:48) ) (cid:16) x k (cid:102) x k − (cid:104) x (cid:102) x (cid:105) (cid:17) + O (cid:0) x k (cid:1) x k (cid:102) x k − (cid:104) x (cid:102) x (cid:105) = x k − x k (cid:104) x (cid:105) − (cid:104) x (cid:105) + (cid:104) x (cid:104) x (cid:105)(cid:105) = (cid:102) x k − x k (cid:104) x (cid:105) (cid:101) z k = (1 − d + iγ (cid:48) ) (cid:102) x k + 2 a (1 − d + iγ (cid:48) ) (1 + d + iγ (cid:48) ) (cid:16)(cid:102) x k − x k (cid:104) x (cid:105) (cid:17) + O (cid:0) x k (cid:1) z k = (1 − d + iγ (cid:48) ) x k + O (cid:0) x k (cid:1)(cid:101) z k = (1 − d + iγ (cid:48) ) (cid:102) x k + O (cid:0) x k (cid:1) z k (cid:104) z (cid:105) = (cid:16) (1 − d + iγ (cid:48) ) x k + a (1 + d + iγ (cid:48) ) (cid:102) x k (cid:17) · (cid:104) (1 − d + iγ (cid:48) ) x + 2 a (1 − d + iγ (cid:48) ) (1 + d + iγ (cid:48) ) x (cid:102) x (cid:105) + O (cid:0) x k (cid:1) = (cid:16) (1 − d + iγ (cid:48) ) x k + a (1 + d + iγ (cid:48) ) (cid:102) x k (cid:17) (cid:104) (1 − d + iγ (cid:48) ) x (cid:105) + O (cid:0) x k (cid:1) = (1 − d + iγ (cid:48) ) x k (cid:104) x (cid:105) + O (cid:0) x k (cid:1) . Using these terms, we can write (cid:104) e z (cid:105) (cid:103) e − z k = − z k + 12 (cid:101) z k − (cid:101) z k − z k (cid:104) z (cid:105) + O (cid:0) x k (cid:1) = − (1 − d + iγ (cid:48) ) x k − a (1 + d + iγ (cid:48) ) (cid:102) x k + 12 (1 − d + iγ (cid:48) ) (cid:102) x k + a (1 − d + iγ (cid:48) ) (1 + d + iγ (cid:48) ) (cid:16)(cid:102) x k − x k (cid:104) x (cid:105) (cid:17) −
16 (1 − d + iγ (cid:48) ) (cid:102) x k −
12 (1 − d + iγ (cid:48) ) x k (cid:104) x (cid:105) = − (1 − d + iγ (cid:48) ) x k + (cid:18)
12 (1 − d + iγ (cid:48) ) − a (1 + d + iγ (cid:48) ) (cid:19) (cid:102) x k + (cid:18) a (1 − d + iγ (cid:48) ) (1 + d + iγ (cid:48) ) −
16 (1 − d + iγ (cid:48) ) (cid:19) (cid:102) x k + (cid:18) −
12 (1 − d + iγ (cid:48) ) − a (1 − d + iγ (cid:48) ) (1 + d + iγ (cid:48) ) (cid:19) x k (cid:104) x (cid:105) e R (cid:103) e r k = (cid:0) b (cid:104) x (cid:105) (cid:1) (cid:18) r k + 2 (cid:101) r k + 43 (cid:101) r k (cid:19) + O (cid:0) x k (cid:1) = 2 r k + 4 br k (cid:104) x (cid:105) + 2 (cid:101) r k + 43 (cid:101) r k + O (cid:0) x k (cid:1) r k = (1 − d ) x k + (1 + d ) y k = (1 − d ) x k + (1 + d ) a (cid:102) x k r k = (1 − d ) x k + 2 a (1 − d ) (1 + d ) x k (cid:102) x k + O (cid:0) x k (cid:1) r k = (1 − d ) x k + O (cid:0) x k (cid:1)(cid:101) r k = r k − (cid:104) r (cid:105) = (1 − d ) (cid:102) x k + 2 a (1 − d ) (1 + d ) (cid:16)(cid:102) x k − x k (cid:104) x (cid:105) (cid:17)(cid:101) r k = (1 − d ) (cid:102) x k + O (cid:0) x k (cid:1) e R (cid:103) e r k = 2 r k + 4 br k (cid:104) x (cid:105) + 2 (cid:101) r k + 43 (cid:101) r k + O (cid:0) x k (cid:1) = 2 (1 − d ) x k + 2 a (1 + d ) (cid:102) x k + 4 b (1 − d ) x k (cid:104) x (cid:105) + 2 (1 − d ) (cid:102) x k + 4 a (1 − d ) (1 + d ) (cid:16)(cid:102) x k − x k (cid:104) x (cid:105) (cid:17) + 43 (1 − d ) (cid:102) x k + O (cid:0) x k (cid:1) = 2 (1 − d ) x k + (cid:16) a (1 + d ) + 2 (1 − d ) (cid:17) (cid:102) x k + (cid:18) a (1 − d ) (1 + d ) + 43 (1 − d ) (cid:19) (cid:102) x k + (4 b (1 − d ) − a (1 − d ) (1 + d )) x k (cid:104) x (cid:105) . We can now insert the different orders of x k from e R (cid:103) e r k and (cid:104) e z (cid:105) (cid:103) e − z k in Eq. (C.2) (here, we use sympy [40] to solvefor the coefficients), yielding d ˙ x k = ( d − (cid:16) γ (cid:48) + (1 + d ) (cid:17) (cid:16) γ (cid:48) − d − (cid:17) γ (cid:48) (cid:102) x k − ( d − (cid:16) γ (cid:48) + ( d − (cid:17) (cid:16) γ (cid:48) + ( d + 1) (cid:17) (cid:0) γ (cid:48) − γ (cid:48) d + 3 (cid:0) d − (cid:1)(cid:1) (cid:0) γ (cid:48) + 2 γ (cid:48) d + 3 (cid:0) d − (cid:1)(cid:1) γ (cid:48) d (cid:102) x k + (1 − d ) d (cid:18) γ (cid:48) − γ (cid:48) (cid:0) d − d + 1 (cid:1) − (cid:0) d + d − d + 22 d + 1 (cid:1) − γ (cid:48) (cid:0) d + 5 d − d − d + 2 d + 11 d − (cid:1) + 9 γ (cid:48) (cid:0) d − (cid:1) (cid:19) x k (cid:104) x (cid:105) . A = ( d − (cid:16) γ (cid:48) + (1 + d ) (cid:17) (cid:16) γ (cid:48) − d − (cid:17) γ (cid:48) dB = − ( d − (cid:16) γ (cid:48) + ( d − (cid:17) (cid:16) γ (cid:48) + ( d + 1) (cid:17) (cid:0) γ (cid:48) − γ (cid:48) d + 3 (cid:0) d − (cid:1)(cid:1) (cid:0) γ (cid:48) + 2 γ (cid:48) d + 3 (cid:0) d − (cid:1)(cid:1) γ (cid:48) d C = ( d − d γ (cid:48) (cid:18) γ (cid:48) − γ (cid:48) (cid:0) d − d + 1 (cid:1) − γ (cid:48) (cid:0) d + d − d + 22 d + 1 (cid:1) − γ (cid:48) (cid:0) d + 5 d − d − d + 2 d + 11 d − (cid:1) + 9 (cid:0) d − (cid:1) (cid:19) . D 2-Cluster states in the center manifold
For 2-cluster states, we can take N = N + N and write ˙ x k = λ + x k + A (cid:102) x k + B (cid:102) x k + C (cid:104) x (cid:105) x k + O (cid:0) x k (cid:1) = λ + x k + A (cid:18) x k − N (cid:0) N x + N x (cid:1)(cid:19) + B (cid:18) x k − N (cid:0) N x + N x (cid:1)(cid:19) + CN (cid:0) N x + N x (cid:1) x k with the constraint k ∈ { , } and N x + N x = 0 , that is, x = − ( N /N ) x . Note that ˙ x k must vanish at the2-cluster equilibria. The 2-cluster therefore satisfies λ + x + A (cid:18) x − N (cid:18) N x + N N x (cid:19)(cid:19) + B (cid:18) x − N (cid:18) N x − N N x (cid:19)(cid:19) + CN (cid:18) N x + N N x (cid:19) x = λ + x + A (cid:18) x − N N x (cid:19) + B (cid:18) x − N ( N − N ) N x (cid:19) + CN N x = λ + x + A N − N N x + B N − N ( N − N ) N x + CN N x , and writing α = N /N , λ + x + A (1 − α ) x + (cid:0) B (cid:0) − α + α (cid:1) + Cα (cid:1) x . This equation has the solutions x = 0 , x = 0 and x ± = 12 ( B (1 − α + α ) + Cα ) (cid:18) − A (1 − α ) ± (cid:113) A (1 − α ) − λ + ( B (1 − α + α ) + Cα ) (cid:19) x ± = − ( N /N ) x ± . The saddle-node curves creating the 2-cluster solutions are thus parametrized by the vanishing discriminant A (1 − α ) − λ + (cid:0) B (cid:0) − α + α (cid:1) + Cα (cid:1) ⇒ λ + = λ sn = A (1 − α ) B (1 − α + α ) + Cα ) for unbalanced cluster solutions, that is, α (cid:54) = 1 or N (cid:54) = N . Thus, at the saddle-node bifurcation x ± = x ∗ = − A (1 − α )2 ( B (1 − α + α ) + Cα ) . References [1] Steven H. Strogatz. From Kuramoto to Crawford: Exploring the onset of synchronization in populations ofcoupled oscillators.
Physica D: Nonlinear Phenomena , 143(1-4):1–20, 2000.[2] Vladimir García-Morales and Katharina Krischer. Normal-form approach to spatiotemporal pattern formation inglobally coupled electrochemical systems.
Physical Review E , 78(5):057201, 2008.20-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center Manifold[3] István Z. Kiss, Yumei Zhai, and John L. Hudson. Characteristics of cluster formation in a population of globallycoupled electrochemical oscillators: an experiment-based phase model approach.
Progress of Theoretical PhysicsSupplement , 161(nil):99–106, 2006.[4] Wen Wang, István Z. Kiss, and J. L. Hudson. Experiments on arrays of globally coupled chaotic electrochemicaloscillators: Synchronization and clustering.
Chaos: An Interdisciplinary Journal of Nonlinear Science , 10(1):248–256, 2000.[5] Hamilton Varela, Carsten Beta, Antoine Bonnefont, and Katharina Krischer. A hierarchy of global couplinginduced cluster patterns during the oscillatory H2-electrooxidation reaction on a Pt ring-electrode.
PhysicalChemistry Chemical Physics , 7(12):2429, 2005.[6] F. Plenge, H. Varela, and K. Krischer. Pattern formation in stiff oscillatory media with nonlocal coupling: anumerical study of the hydrogen oxidation reaction on Pt electrodes in the presence of poisons.
Physical ReviewE , 72(6):066211, 2005.[7] Konrad Schönleber, Carla Zensen, Andreas Heinrich, and Katharina Krischer. Pattern formation during theoscillatory photoelectrodissolution of n-type silicon: Turbulence, clusters and chimeras.
New Journal of Physics ,16(6):063024, 2014.[8] Minseok Kim, Matthias Bertram, Michael Pollmann, Alexander von Oertzen, Alexander S. Mikhailov, Harm Hin-rich Rotermund, and Gerhard Ertl. Controlling chemical turbulence by global delayed feedback: Pattern formationin catalytic CO oxidation on Pt(110).
Science , 292(5520):1357–1360, 2001.[9] Toby Elmhirst.
Symmetry and emergence in polymorphism and sympatric speciation . PhD thesis, University ofWarwick, Warwick, 2001.[10] Ian Stewart, Toby Elmhirst, and Jack Cohen.
Symmetry-Breaking as an Origin of Species , pages 3–54. Bifurcation,Symmetry and Patterns. Birkhäuser Basel, 2003.[11] Koji Okuda. Variety and generality of clustering in globally coupled oscillators.
Physica D: Nonlinear Phenomena ,63(3-4):424–436, 1993.[12] Naoko Nakagawa and Yoshiki Kuramoto. From collective oscillations to collective chaos in a globally coupledoscillator system.
Physica D: Nonlinear Phenomena , 75(1-3):74–80, 1994.[13] Murad Banaji. Clustering in globally coupled oscillators.
Dynamical Systems , 17(3):263–285, 2002.[14] Hiroaki Daido and Kenji Nakanishi. Aging and clustering in globally coupled oscillators.
Physical Review E ,75(5):056206, 2007.[15] Wai Lim Ku, Michelle Girvan, and Edward Ott. Dynamical transitions in large systems of mean field-coupledLandau-Stuart oscillators: Extensive chaos and cluster states.
Chaos: An Interdisciplinary Journal of NonlinearScience , 25(12):123122, 2015.[16] Arkady Pikovsky and Michael Rosenblum. Dynamics of globally coupled oscillators: Progress and perspectives.
Chaos: An Interdisciplinary Journal of Nonlinear Science , 25(9):097616, 2015.[17] Y. Kuramoto.
Chemical Oscillations, Waves and Turbulence , volume 19 of
Springer Series in Synergetics .Springer-Verlag Berlin Heidelberg, 1984.[18] Shinya Watanabe and Steven H. Strogatz. Integrability of a globally coupled oscillator array.
Physical ReviewLetters , 70(16):2391–2394, 1993.[19] Shinya Watanabe and Steven H. Strogatz. Constants of motion for superconducting Josephson arrays.
Physica D:Nonlinear Phenomena , 74(3-4):197–253, 1994.[20] Martin Golubitsky and Ian Stewart.
Linear Stability , pages 33–57. Birkhäuser Basel, Basel, 2002.[21] Vincent Hakim and Wouter-Jan Rappel. Dynamics of the globally coupled complex Ginzburg-Landau equation.
Physical Review A , 46(12):R7347–R7350, 1992.[22] T. Brooke Benjamin and J. E. Feir. The disintegration of wave trains on deep water part 1. Theory.
Journal ofFluid Mechanics , 27(3):417–430, 1967.[23] Richard FitzHugh. Mathematical models of threshold phenomena in the nerve membrane.
The Bulletin ofMathematical Biophysics , 17(4):257–278, 1955.[24] Balth. van der Pol. On relaxation-oscillations.
The London, Edinburgh, and Dublin Philosophical Magazine andJournal of Science , 2(11):978–992, 1926.[25] Richard J. Field and Richard M. Noyes. Oscillations in chemical systems. iv. limit cycle behavior in a model of areal chemical reaction.
The Journal of Chemical Physics , 60(5):1877–1884, 1974.21-Cluster Fixed-Point Analysis of Mean-Coupled Stuart-Landau Oscillators in the Center Manifold[26] A. L. Hodgkin and A. F. Huxley. A quantitative description of membrane current and its application to conductionand excitation in nerve.
The Journal of Physiology , 117(4):500–544, 1952.[27] Emilio Andreozzi, Ilaria Carannante, Giovanni D’Addio, Mario Cesarelli, and Pietro Balbi. Phenomenologicalmodels of nav1.5. a side by side, procedural, hands-on comparison between hodgkin-huxley and kinetic formalisms.
Scientific Reports , 9(1):17493, 2019.[28] Ana Paula S. Dias and Ana Rodrigues. Secondary bifurcations in systems with all-to-all coupling. part ii.
Dynamical Systems , 21(4):439–463, 2006.[29] D.G. Aronson, G.B. Ermentrout, and N. Kopell. Amplitude response of coupled oscillators.
Physica D: NonlinearPhenomena , 41(3):403–449, 1990.[30] Arkady Pikovsky, Michael Rosenblum, and Jürgen Kurths.
Synchronization , pages 222–235. Cambridge UniversityPress (CUP), 2001.[31] André Röhm, Kathy Lüdge, and Isabelle Schneider. Bistability in two simple symmetrically coupled oscillatorswith symmetry-broken amplitude- and phase-locking.
Chaos: An Interdisciplinary Journal of Nonlinear Science ,28(6):063114, 2018.[32] Felix P. Kemeth, Sindre W. Haugland, and Katharina Krischer. Cluster singularity: the unfolding of clusteringbehavior in globally coupled Stuart-Landau oscillators.
Chaos: An Interdisciplinary Journal of Nonlinear Science ,29(2):023107, 2019.[33] N. Nakagawa and Y. Kuramoto. Collective chaos in a population of globally coupled oscillators.
Progress ofTheoretical Physics , 89(2):313–323, 1993.[34] E. J. Doedel. AUTO: A program for the automatic bifurcation analysis of autonomous systems. In
Congressnumerantium , volume 30, page nil, 4 1981.[35] E. J. Doedel and X. J. Wang. AUTO-07P: Continuation and bifurcation software for ordinary differential equations.Technical report, Center for Research on Parallel Computing, California Institute of Technology, Pasadena CA91125, 2007.[36] Bernold Fiedler, Sindre W. Haugland, Felix Kemeth, and Katharina Krischer. Global 2-cluster dynamics underlarge symmetric groups, 2020. arXiv preprint arXiv:2008.06944.[37] Paula S Dias and Ian Stewart. Secondary bifurcations in systems with all-to-all coupling.
Proceedings of the RoyalSociety of London. Series A: Mathematical, Physical and Engineering Sciences , 459(2036):1969–1986, 2003.[38] Laurette S. Tuckerman and Dwight Barkley. Bifurcation analysis of the eckhaus instability.
Physica D: NonlinearPhenomena , 46(1):57–86, 1990.[39] Felix P. Kemeth.
Symmetry Breaking in Networks of Globally Coupled Oscillators: From Clustering to ChimeraStates . PhD thesis, Technische Universität München, Garching, 4 2019.[40] Aaron Meurer, Christopher P. Smith, Mateusz Paprocki, Ondˇrej ˇCertík, Sergey B. Kirpichev, Matthew Rocklin,AMiT Kumar, Sergiu Ivanov, Jason K. Moore, Sartaj Singh, Thilina Rathnayake, Sean Vig, Brian E. Granger,Richard P. Muller, Francesco Bonazzi, Harsh Gupta, Shivam Vats, Fredrik Johansson, Fabian Pedregosa, Matthew J.Curry, Andy R. Terrel, Štˇepán Rouˇcka, Ashutosh Saboo, Isuru Fernando, Sumith Kulal, Robert Cimrman, andAnthony Scopatz. Sympy: Symbolic computing in python.