Bifurcations of Clusters and Collective Oscillations in Networks of Bistable Units

We investigate dynamics and bifurcations in a mathematical model that captures electrochemical experiments on arrays of microelectrodes. In isolation, each individual microelectrode is described by a one-dimensional unit with a bistable current-potential response. When an array of such electrodes is coupled by controlling the total electric current, the common electric potential of all electrodes oscillates in some interval of the current. These coupling-induced collective oscillations of bistable one-dimensional units are captured by the model. Moreover, any equilibrium is contained in a cluster subspace, where the electrodes take at most three distinct states. We systematically analyze the dynamics and bifurcations of the model equations: We consider the dynamics on cluster subspaces of successively increasing dimension and analyze the bifurcations occurring therein. Most importantly, the system exhibits an equivariant transcritical bifurcation of limit cycles. From this bifurcation, several limit cycles branch, one of which is stable for arbitrarily many bistable units.

The interaction between simple individual units can lead to the spontaneous occurrence of spatial-or/and temporal patterns. A lot of focus has been on the collective dynamics of coupled individual (microscopic) oscillators. However, in some important experimental contexts, the individual units are one-dimensional and bistable rather than oscillatory. Examples include phase transition cathodes in Li-ion batteries and the electrocatalytic oxidation of COon-Pt-nanoparticles or arrays of Pt-microelectrodes. For the microelectrode array, one can observe that fixing the total current through the electrodes leads to clustered dynamics of two or three clusters where all units within one cluster take on the same state. Moreover, such clustered states can exhibit sustained collective oscillations. These dynamic features are captured in a simple mathematical model 1 . Here we elucidate the dynamics and bifurcations that arise in the model equations as parameter are varied. We exploit the fact that in a stationary state, each of the N bistable units takes on only one of three different values and reduce the N-dimensional system to two or three degrees of freedom that describe the clustered dynamics. This allows us to perform a bifurcation analysis of the steady states and limit cycles, which arise in Hopf bifurcations. Furthermore, using knowledge of the stability of the clusters in the full system 1 , we can conclude that stable oscillations also exist in thermodynamically large systems with 10 10 individual bistable units, for example, in electrocatalytical systems composed of nanoparticular catalyst particles on some support. a) Electronic mail: munir.salman@tum.de b) Electronic mail: krischer@tum.de

I. INTRODUCTION
Networks of interacting bistable units are an important class of dynamical systems that are relevant in diverse physical contexts. While coupled units with oscillatory or chaotic dynamics have been widely considered (see, e.g., the textbook Ref. 2), multi-component systems in which each unit has a single degree of freedom have received much less attention. Zanette 3 considered globally coupled bistable elements which interacted linearly through their common mean. The global coupling was diffusive in the sense that it vanished for uniform states. For random initial condition and a sufficiently large coupling strength, the global interaction led to a coherent motion of the entire ensemble towards one of the two stable states. The impact of additive noise to the bistable units in related network configurations was analyzed in Refs. 4-7, while Ref. 8-10 considered the collective response of bistable units that were connected through different coupling topologies. A recent study also examined the transition from a local to a global coupling topology, focusing on the mathematical perspective 11 .
Kouvaris et al. 12 investigated the effect of a global feedback on the dynamics of networks of diffusively coupled bistable units. The global feedback was chosen such that it altered the excitation threshold, i.e., the position of the saddle point, of the local one-dimensional reaction dynamics. In this case, localized stationary activation patterns form, the size of which can be adjusted by varying the feedback strength. The impact of a global time-delayed feedback and additive noise on an ensemble of bistable units was examined in Refs. 4 and 13, while in Ref. 14 the impact of heterogeneities on the synchronization state of globally coupled bistable electronic circuits is discussed.
Here, we consider bistable units that interact through a macroscopic observable which is forced to take a fixed value as a global constraint. This observable can be the sum of the intrinsic state variables of the individual units, or it can be the sum of a function of the individual state variables. The global constraint is enforced by allowing a parameter that controls the state of an individual unit on the equilibrium branch to adapt; thus the global constraint can be seen as a bifurcation parameter of the coupled system. Indeed, such coupling between bistable units is realized in many physical systems: These include Li-ion batteries with phase-transition cathodes 15,16 and bistable electrochemical reactions, such as the CO oxidation on an array of Pt electrodes [17][18][19] . In case of Liion batteries, the cathodes consist of billions of nano-particles that can be considered as bistable units. Each nanoparticle can be in a Li-rich or a Li-poor state depending on the chemical potential or equivalently the voltage. Thus, the battery can be seen as interacting bistable units coupled through a global constraint: When charging or discharging the battery slowly, a constraint is set to the time evolution of the total charge, while the voltage adjusts accordingly. Similarly, the CO-oxidation on Pt is a prototypical electrocatalytic reaction which exhibits bistable reaction rates, and thus a bistable current-voltage characteristic. In technological applications, electrodes consist of billions of catalytically active nanoparticles on a nonreactive support that interact globally when a set current is passed through the electrode. A setup that enables the measurement of the state of an individual bistable component is an array of Pt-electrodes. Such measurements reveal key dynamical properties: The bistable units activate sequentially upon a slow parameter ramp, they form clusters with most electrodes in the two stable states and at most one electrode on the unstable state of the individual electrode, and they may exhibit collective oscillations of the entire ensemble [17][18][19] . We recently derived a general necessary condition when a system of bistable units subject to a global constraint can become unstable in a Hopf bifurcation 1 and validated this condition in a simplified CO-electrooxidation model.
In this paper, we analyze this CO-electrooxidation model and elucidate its bifurcation properties. The behavior of the coupled units is particularly interesting when the value of the constraint is in the transitional regime between the two outer 'active' and 'passive' states of an uncoupled individual unit. As a globally and identically coupled system of identical units, the bifurcations are constrained by the symmetry properties of the system: The dynamical equations are S Nsymmetric (equivariant), where S N denotes the group of permutations of N symbols that acts by permuting component indices; see Ref. 20 for a general introduction to dynamical systems with symmetry. The presence of symmetries gives rise to dynamically invariant subspaces, that correspond, for example, to cluster configurations where the state of some components coincide. Symmetric systems can exhibit bifurcation behavior that is nongeneric for a general dynamical system as their bifurcation behavior is constrained by the symmetry 20,21 . We will consider symmetry breaking bifurcations where the configuration of all units being synchronized (the one-cluster configuration) loses stability. Indeed, since the state space of individual unit is one-dimensional for the system we consider, collective oscillations cannot occur if all units are synchronized. Thus the emergence of collective oscillations necessarily requires symmetry-breaking bifurcations away from full synchrony.
In our bifurcation analysis we concentrate specifically on the structure imposed by the invariant subspaces due to symmetry. Here the one-, two-, and three cluster subspaces are crucially important as they contain all equilibria of the system. We give a detailed bifurcation analysis of the steady states in these subspaces and highlight bifurcations that give rise to and stabilize limit cycles as they induce collective oscillations observed in experiments. Specifically, the paper is organized as follows. In Section II we summarize the model equations. In Section III, we consider the one-dimensional, fully synchronized dynamics where the states of all units take the same value and form one cluster. In Section IV, we discuss symmetry breaking bifurcations away from full synchrony to two-cluster equilibria and bifurcations of the two-cluster equilibria; these can give rise to collective oscillations (that are transversely unstable for most cluster sizes). In Section V, we outline bifurcations of three cluster equilibria and indicate that the transversely unstable collective oscillations within the two-cluster subspace can be stabilized in a transcritical bifurcation. We conclude with some remarks in Section VI.

A. Model equation
We consider the following simplified model for COelectrooxidation on Pt-microelectrodes. Its derivation from an established, more detailed model version can be found in the supplement of Ref. 1. The intrinsic state x k of unit k evolves according tȯ where y imposes a global constraint and a and b are positive constants; throughout this paper, we fix a = 0.05, b = 0.01 as in Ref. 1. In terms of the CO dynamics, x k is the COcovered fraction of the surface and u the electrode potential. The macroscopic observable y tot is the total CO oxidation current, and y denotes the mean current per electrode (or nanoparticle). In our setting, y is the bifurcation parameter, while u adjusts such that the system can attain the preset mean value. Note that this type of coupling corresponds to a control that fixes the value of a function of the mean state of the ensemble and allows the 'natural control parameter'-in our case u-to change. This situation is different from the situation considered in Ref. 3, where the main bifurcation parameter of the bistable unit is also the bifurcation parameter of the globally coupled ensemble, and thus set to a constant value. First consider a single uncoupled unit, N = 1, with state x 1 = x; the corresponding equilibrium dynamicsẋ = 0 are shown in Fig. 1. When u is treated as a parameter, x is bistable in u with equilibrium branch as in Fig. 1a; there are two saddle-node (SN) bifurcations that bound the region of bistability. For ease of notation, we subdivide the branch of equilibria into three segments specified by the points where these SN bifurcations occur and the branch turns over: The active segment ("a"), where the current is high, the middle ("m") segment, and the passive segment ("p"), where the current is low (see Fig. 1). We denote the value of u corresponding to the transitions between the passive and middle segment with u p→m , and the transition between the middle and active segment with u m→a . The transitions between segments play an explicit role in the bifurcation analysis as we will discuss below. For the equilibrium values of x, there is a single-valued, but non-monotonic relationship between u and y, depicted in Fig. 1b. By contrast the value of x at equilibrium depends monotonically on y (Fig. 1c). The dynamical equations (1) can be rearranged to the explicit formẋ with the electrode potential which is determined by the intrinsic states x k and the bifurcation parameter y. If we define the two functions the dynamical equations take the forṁ Configurations where the states of the units form three clusters are particularly interesting because all equilibria of Eq. (5) necessarily lie in a three cluster subspace. To see this, consider the fixed points x * k of Eq. (1) with u = u * which satisfy Since this cubic polynomial has at most three distinct roots, the coefficients of any equilibrium of (5) take at most three distinct values and thus form a cluster configuration.

B. Cluster subspaces and full synchrony
As a globally coupled network of N identical bistable units, the system equation (5) are symmetric (equivariant) with respect to permutations of units. In other words, if the group S N of permutations of N symbols acts on the bistable units by permuting indices, the equations of motion remain unchanged. This implies that cluster configurations, where the state of all units in each cluster is identical, are also invariant under the dynamics.
Suppose that the units form M clusters such that the kth cluster contains N k units with identical state ξ k . Evidently, we have to have N = ∑ M j=1 N j . Define the relative cluster size n k = N k N and write C (n 1 ,...,n M ) for the cluster subspace where the first N 1 oscillators form cluster 1, the next N 2 oscillators cluster 2, and so on. The effective dynamics on C (n 1 ,...,n M )and by symmetry any other cluster configuration with these cluster sizes-are M-dimensional and the state ξ k of cluster k ∈ {1, . . . , M} evolves according tȯ y.
While n k only takes finitely many values for finite networks, in the limit of infinitely many units, N → ∞, we can see n k as a continuous parameter.

III. FULLY SYNCHRONIZED DYNAMICS
The simplest cluster configuration C (1) is a single cluster, that is, all units are synchronized with x k = ξ 1 = ξ . With P, Q as above, the synchronized dynamics of (5) are given bẏ Eq. (8) has only one single equilibrium point for any y ∈ (0, 1), which is identical to the curve of an individual element depicted in Fig. 1. Restricted to the synchronized subspace, the system is monostable for all parameter values, i.e., for each y there is exactly one single equilibrium. (If u is given instead of y, however, there can be three equilibria.) Linearizing Eq. (8) at the equilibrium ξ * yieldṡ and thus ξ * is always stable within C (1) . The stability of ξ * as an equilibrium of the full system (5) depends on the N − 1 eigenvalues transverse to C (1) . Due to the symmetry, all transversal eigenvalues are the same; they are all negative on the "a" and "p" segments, and positive on the "m" segment. This implies that the equilibrium is stable on the "a" and "p" segments and unstable on the "m" segment. The change of stability of the equilibrium ξ * happens in S N -equivariant transcritical bifurcations at u p→m and u m→a where all eigenvalues transverse to C (1) pass through zero simultaneously. Note that at these transcritical bifurcations an uncoupled individual bistable unit undergoes a SN bifurcation if u is considered as a bifurcation parameter (cf. Eq. (1a)). That is, the bifurcation happens where all synchronized units jointly transition from the passive to the middle segment of the branch or from the middle to the active segment of the branch.

IV. TWO-CLUSTER DYNAMICS AND BIFURCATIONS
Now, we consider two-cluster configurations C (n 1 ,n 2 ) with relative sizes n 1 and n 2 , respectively. The two-dimensional dynamics on C (n 1 ,n 2 ) are given bẏ The intersection of all possible two-cluster subspaces contain the fully synchronized subspace {ξ 1 = ξ 2 }. Note that for n 1 and n 2 fixed, the system has a parameter symmetry which corresponds to exchanging the clusters. Since n 1 +n 2 = 1 the parameter n 1 fully determines the system and the parameter symmetry can also be written as Thus, without loss of generality, we may assume n 1 ≥ n 2 . If the clusters are of equal size, n 1 = n 2 , the parameter symmetry yields a system symmetry (ξ 1 , ξ 2 ) → (ξ 2 , ξ 1 ) on the cluster subspace.  Fig. 1. This branch interacts with the two-cluster equilibria (blue line) in transcritical bifurcations (yellow dot). The labels "a", "m", "p" indicate the active, middle, and passive segments of the branch as detailed in the main text with the first letter representing ξ 1 , the second ξ 2 . A Hopf bifurcation (green dot) gives rise to oscillatory dynamics. The shaded parameter range is considered in detail in Fig. 3.
The branch of fully symmetric equilibria in C (1) interacts with branches of two-cluster solutions in C (n 1 ,n 2 ) in transcritical bifurcations. Note that we say that two (distinct) branches of equilibria or limit cycles interact if they are involved in the same bifurcation; for example, two branches of equilibria interact in a simple transcritical bifurcation. All branches and their bifurcations described throughout the paper were calculated numerically using AUTO-07P 22 . In the following, we will first focus on the specific relative cluster sizes n 1 = 0.8, n 2 = 0.2 to illustrate the bifurcation behavior. Fig. 2 shows the equilibria branches in the (u, y)-plane as the bifurcation parameter y is varied. The black line shows the fully symmetric equilibrium ξ = = (ξ * , ξ * )-with ξ * as in Section III -within C (1) ⊂ C (n 1 ,n 2 ) ; this branch is the same as the one shown in Fig. 1b. The locations where the fully symmetric equilibrium solution undergoes transcritical bifurcations are indicated by yellow dots. At these bifurcation points, the symmetric branch interacts with the branch of two-cluster equilibria (dark blue line). It is useful to subdivide the branch of two-cluster equilibria into segments: We write "aa" if both clusters are on the active segment, "am" if the first cluster is on the active and the second one on the middle segment, etc. This leads to the labelling shown in Fig. 2 and the boundaries of the segments are given by the local maxima and minima corresponding to u ∈ {u p→m , u m→a }. This labelling yields some intuition of what happens at the transcritical bifurcations: At the first transcritical bifurcations the two-cluster equilibrium transitions from "pm" to "mp" and the clusters switch roles as the majority cluster (n 1 = 0.8) switches from "p" to "m" and the minority cluster (n 2 = 0.2) from "m" to "p". Similarly, at the second transcritical bifurcation there is a transition from "ma" to "am". Following the two-cluster equilibrium branch further, it folds over twice in SN bifurcations (green dots) and undergoes a Hopf bifurcation at y ≈ 0.8 (purple dot). This gives rise to a branch of limit cycles (light blue) that is stable within C (n 1 ,n 2 ) and ends in a homoclinic bifurcation close to the transcritical bifurcation at y ≈ 0.9.
There is actually a series of bifurcations close to the second transcritical bifurcation as the parameter y is varied as shown in Fig. 3. When the Hopf bifurcation occurs at y ≈ 0.8, there are three equilibria in C (n 1 ,n 2 ) , the symmetric equilibrium ξ = in C (1) ⊂ C (n 1 ,n 2 ) as well as nonsymmetric equilib- the former gives rise to the stable limit cycle. At the transcritical bifurcation ξ − interacts with ξ = and satisfies ξ − 1 < ξ − 2 after the bifurcation point. This equilibrium then collides with the limit cycle in a homoclinic bifurcation at y ≈ 0.91. Finally, the two equilibria ξ − and ξ + meet in the fold where the nonsymmetric branch in Fig. 2 folds over itself. Note that the transcritical bifurcation must happen before the homoclinic bifurcation to allow for this bifurcation scenario: It cannot occur on the "ma" segment before the transcritical bifurcation but has to occur on the "am" segment (but before ξ − and ξ + merge in the SN bifurcation).
Considering transverse stability beyond C (n 1 ,n 2 ) we realize that the branch of two-cluster equilibria undergoes transcritical bifurcations whenever u ∈ {u p→m , u m→a }. Hence, there are six in total, the two indicated by yellow dots in Fig. 2, where all bifurcating branches lie within C (n 1 ,n 2 ) , and four additional ones where the two-cluster equilibria interact with three-cluster equilibria; we will discuss these in more detail in the next section.
To elucidate the dynamics and bifurcations for other relative cluster sizes, we continued the bifurcations in two parameters, y and n 1 . The resulting two-parameter bifurcation lines are shown in Fig. 4. Note that the symmetry of the bifurcation diagram reflects the parameter symmetry (12). The transcritical bifurcation lines are horizontal since the bifurcation condition does not depend on the relative cluster sizes: Consider equilibria on C (n 1 ,n 2 ) , that is, solutions of Eq. (6) with two distinct real roots ξ * 1 = ξ * 2 (a nonsymmetric equilibrium on C (n 1 ,n 2 ) ). Now suppose that these roots coincide at a bifurcation for u = u * , defined in (3), which means that (6) has a triple root ξ * 1 = ξ * 2 = ξ * (an equilibrium on C (1) ⊂ C (n 1 ,n 2 ) ). Then, the associated value of y is independent of n 1 . Furthermore, for any given relative cluster sizes (n 1 , n 2 ), there is a branch of equilibria in C (n 1 ,n 2 ) . These branches arise in saddle-node bifurcations and interact with the one-cluster equilibrium in the transcritical bifurcation. The more asymmetric the clusters are, the larger is the distance of the saddle-node bifurcation to the transcritical bifurcation in the bifurcation parameter y. The Hopf and homoclinic bifurcations of the two-cluster states described above extend to a large range of cluster sizes from the most asymmetric distribution n 1 = 1 (or, equivalently, n 1 = 0) to n 1 ≈ 0.4 (respectively n 1 ≈ 0.6). They meet in a codimension-two Bogdanov-Takens bifucation point at (y, n 1 ) ≈ (0.907, 0.34) ((y, n 1 ) ≈ (0.907, 0.66), respectively). For increasing (respectively decreasing) cluster size n 1 , the Hopf bifurcation curve approaches the saddle node bifurcation curve again. In this region of parameter space, the limit cycle solution grows extremely quickly as in a Canard phenomenon 23 before it is destroyed in a branch of homoclinic bifurcations. In the full system (5) all these bifurcations happen simultaneously in a large number of invariant subspaces. On the one hand, for fixed (n 1 , n 2 ) the bifurcations happen simultaneously in the symmetric copies of C (n 1 ,n 2 ) . On the other hand, recall that (n 1 , n 2 ) parameterize invariant subspaces of different cluster sizes. The transcritical bifurcations for all the different (n 1 , n 2 ) happen at the same value of y given by (13). Thus, the symmetric equilibrium undergoes a transcritical bifurcation simultaneously in all these invariant subspaces as one would expect in our symmetric setting. In the case of symmetric clusters, n 1 = n 2 = 0.5, these degenerate to a symmetric pitchfork bifurcation, in line with the reflection symmetry (ξ 1 , ξ 2 ) → (ξ 2 , ξ 1 ) on C (n 1 ,n 2 ) ; cf. Fig. 4.

A. Equivariant bifurcation of equilibria
Since any equilibrium is constrained to a three-cluster configuration due to condition (6), any branch of steady states is contained in a cluster subspace C (n 1 ,n 2 ,n 3 ) . Fig. 5 shows the equilibrium branches in C (n 1 ,n 2 ,n 3 ) for relative cluster sizes n 1 = 0.8, n 2 = 0.16, and n 3 = 0.04. We have FIG. 3. Phase portraits elucidate the dynamics in C (n 1 ,n 2 ) due to the successive transcritical bifurcation (yellow dot), homoclinic bifurcation (light green dot), and saddle node bifurcation (green dot) for n 1 = 0.8, n 2 = 0.2, for a = 0.05, b = 0.01. Gray solid lines indicate the nullclines; solutions and bifurcation points are colored as in Fig. 2. As the parameter y is increased, the limit cycle (light blue circle) that arises from the Hopf bifurcation of ξ + (upper blue dot in the first panel) grows quickly. At the same time, the two-cluster equilibrium ξ − (lower blue dot in the first panel) crosses through the invariant subspace C (1) -the diagonal in the (ξ 1 , ξ 2 ) plane-in the transcritical bifurcation before it interacts with the limit cycle in a homoclinic bifurcation. Finally, the equilibria ξ + and ξ − merge in a saddle node bifurcation. C (1) ⊂ C (n 1 ,n 2 ,n 3 ) , so the branch of symmetric equilibria in Fig. 1 appears (black line). Note that n 2 + n 3 = 0.2 and thus C (0.8,0.2) ⊂ C (0.8,0.16,0.04) for {ξ 2 = ξ 3 } and the equilibrium branches shown in Fig. 2 Fig. 5 (blue line). In addition, C (0.8,0.16,0.04) also contains the two-cluster subspaces {ξ 1 = ξ 3 } and {ξ 1 = ξ 2 }; these correspond to C (0.84,0. 16) and C (0.96,0.04) respectively. The steady state bifurcations can be read in Fig. 4 and the corresponding branches of equilibria are depicted in Fig. 5 as purple and teal lines, respectively. Finally, there is one branch of three-cluster equilibria. Taken together, for u ∈ [u p→m , u m→a ] there are 3 3 = 27 equilibria because each x j can take one of three values whileẋ j = 0,u = 0).

reappear in
These equilibria now (ex)change their stability properties as they bifurcate. In addition to the transcritical bifurcations of the fully synchronized equilibrium, the two-cluster equilibria also undergo transcritical bifurcations in which threecluster equilibria are involved. More specifically, each of the three two-cluster subspaces contained in C (0.8,0.16,0.04) , that is, ξ 1 = ξ 2 , ξ 2 = ξ 3 , ξ 1 = ξ 3 , contains a transcritical bifurcation for each of the two double root of the equilibrium equation (6). This yields six transcritical bifurcations in addition to transcritical bifurcations of the fully synchronized equilibrium that already appeared in Fig. 2. These bifurcations link to the physical interpretation as in the previous section: At each transcritical bifurcation, the value of u is identical with one of the values at which the individual unit undergoes a saddle node bifurcation when u is considered as a bifurcation parameter. Indeed, it can be shown that when the transcritical bifurcation involves an equilibrium with the larger of the two uvalues, i.e., u = u p→m the components of the three-cluster state which are on the middle and passive branch segment switch the segment while the component on the active segment remains unchanged. For example, if ξ x j , j ∈ {1, 2, 3} denotes the jth component of a three cluster that is on the segment x ∈ {p, m, a}, then an equilibrium ξ * = (ξ a 1 , ξ m 2 , ξ We can evaluate the values at which the transcritical bifurcations occur explicitly as in Sect. IV; in general, they depend on the cluster sizes n 1 , n 2 , n 3 . Consider equilibria on C (n 1 ,n 2 ,n 3 ) , that is, solutions of Eq. (6) with three distinct real roots ξ * 1 = ξ * 2 = ξ * 3 = ξ * 1 (a nonsymmetric equilibrium on C (n 1 ,n 2 ,n 3 ) ). Now suppose two of these roots coincide at a bifurcation for u = u * , defined in (3) so that Eq. (6) has one real double root ξ d and one real single root ξ s . Without loss of generality we assume that ξ * 1 = ξ s , ξ * 2 = ξ * 3 = ξ d on FIG. 4. A codimension-2 Bogdanov-Takens bifurcation organizes the bifurcation behavior in C (n 1 ,n 2 ) as both y and the relative cluster size n 1 is varied; the remaining parameters are n 2 = 1 − n 1 , a = 0.05, and b = 0.01. The symmetry of the bifurcations is due to the parameter symmetry (12). The upper panel shows a magnification close to the codinmension-2 Bogdanov-Takens point where saddle-node, Hopf, and homoclinic bifurcation branches meet. For n 1 = n 2 = 1 2 , the saddle node and transcritical bifurcations merge in a pitchfork symmetry breaking bifurcation. C (n 1 ,n 2 +n 3 ) ⊂ C (n 1 ,n 2 ,n 3 ) . Then, the associated value of y is which is linear in n 1 . The cluster at ξ s is what distinguishes this from the previous lower-dimensional case, where this "bystander" was not present. By symmetry, the condition for the transcritical bifurcation in the other two-cluster subspaces depends linearly on n 2 , n 3 , respectively. Thus the transcritical bifurcations are given by lines as the cluster sizes are varied, shown in Fig. 6. In this figure, we see that eight transcritical bifurcations exist for every n 1 . The value n 1 = 0.8 corresponds to Fig. 5. We see the two transcritical bifurcations in C (n 1 ,n 2 +n 3 ) following Eq. (15). Furthermore, n 2 and n 3 are chosen to be linearly dependent on n 1 (see caption Fig. 6), so the corresponding bifurcations lie on straight lines, too. For the full system with N identical units this means that for each ∈ {2, . . . , N} (i.e., the combined size of cluster 2 and 3, which fall together at the bifurcation) there are N ∑ −1  Fig. 1 and explained in the main text, the labeling with letters "a", "m", "p" indicate the active, middle, and passive segments of the equilibrium branch in cluster 1, 2, and 3, respectively. The black line shows the fully symmetric equilibrium (as in Figs. 1 and 2), blue lines of different color the three different two-cluster equilibrium branches (the dark blue with ξ 2 = ξ 3 is the same one shown in Fig. 2), and the gray line corresponds to a branch of asymmetric equilibria. Two-and three-cluster equilibria interact in symmetric transcritical bifurcations (yellow dots) and other bifurcation points are indicated as above. The shaded parameter range is considered in detail in Fig. 7. is only one such bifurcation for each combination of ∈ {2, . . . , N} (i.e., combined size of cluster 2 and 3) and k ∈ {1, . . . , − 1} (i.e., size of cluster 2). Furthermore, many of these transcritical bifurcations lie on the same points in parameter space and in phase space such that in total there are actually just N − 1 multi-branch bifurcations-points at every u * , i.e. one for each possible size of cluster 1 (the cluster at ξ s ), i.e., one for each ∈ {2, . . . , N}. The number of branches is 2 + 2 ∑ −1 k=1 1 2 k . These bifurcations are invariant under intracluster index-permutations and are very similar to the equivariant bifurcations we saw in lower dimensions. Here, however, they can occur at various different values of y. In fact, for N → ∞ they occur at each of the two u * -values in the entire y-intervals from the transcritical bifurcations almost to the respective fully symmetric state.
As seen above for the two-cluster subspaces, also threecluster equilibria might undergo a Hopf bifurcation. In Fig. 6 the locations of the Hopf-bifurcations of the three-cluster and the different two-cluster equilibria are depicted in the (y, n 1 )plane. The labels of the different curves indicate the segments of the three components of the equilibria from which the limit cycles bifurcated. All of them contain at least one cluster on the middle segment. As shown in Ref. 1, in the full system, any steady state with more than one element on the middle branch is unstable. Hence, even if the Hopf-bifurcations are FIG. 6. Numerical continuation shows two-parameter bifurcation as both the bifurcation parameter y and the relative cluster sizes are varied. Here transcritical bifurcation branches (yellow lines; cf. Eq. (15)) as well as Hopf bifurcation branches (purple lines) occurring in three-cluster subspaces are shown for the relative cluster sizes n 2 = 0.8 · (1 − n 1 ), n 3 = 0.2 · (1 − n 1 ) and fixed parameter a = 0.05, b = 0.01. The letters "a", "m", "p" indicate how the clusters distribute on the three segments (see Fig. 1) of the equilibria that undergo the Hopf bifurcation.
supercritical and the limit cycle branch off equilibria that are stable in the cluster subspaces, they are unstable in the full system. The limit cycles which emerge in a Hopf bifurcation will therefore not be observable in systems containing a large number of individual elements. However, as we will demonstrate in the next subsection, the three-cluster limit cycles might interact with a two-cluster limit cycle in a transcritical bifurcation, thereby stabilizing the limit cycle in the two-cluster subspace.

B. Equivariant bifurcation of limit cycles
Consider the phase portraits in the three-cluster subspace C (0.8,0.16,0.04) depicted in Fig. 7. At y = 0.78625, there are two two fixed points and two limit cycles which were both born at slightly smaller values of y in supercritical Hopf bifurcations. One limit cycle and one fixed point are threecluster limit sets, while the other two limit sets lie in the invariant C (n 1 ,n 2 +n 3 ) subspace, i.e., ξ 2 = ξ 3 . The two-cluster limit cycle is stable within the two-cluster subspace but unstable in the direction perpendicular to it, while the threecluster limit cycle is stable within the invariant three-cluster subspace. As y is increased to y = 0.79, we see that the limit cycles have increased in size and have moved closer together. At y ≈ 0.79625 there is a transcritical bifurcation of the limit cycles, i.e., the they coincide on the two-cluster subspace C (n 1 ,n 2 +n 3 ) . Then, at y = 0.8, the three-cluster limit cycle is on the other side of the C (n 1 ,n 2 +n 3 ) plane. In the transcrit- FIG. 7. Variation of the parameter y indicates a transcritical bifurcation of limit cycles as illustrated in phase portraits; the other parameters are n 1 = 0.8, n 2 = 0.16, n 3 = 0.04, a = 0.05, b = 0.01. (a): Two limit cycles have just branched from a Hopf bifurcation at slighly smaller y. One of these two limit cycles lies on a plane, which is the two-cluster subspace C (n 1 ,n 2 +n 3 ) ⊂ C (n 1 ,n 2 ,n 3 ) . The other lies above the plane in the three cluster subspace C (n 1 ,n 2 ,n 3 ) . (b): The two limit cycles grow as y is further increased. (c): The two limit cycles interact in a tanscritical bifurcation of limit cycles. (d): The three-cluster limit cycle is now below the two-cluster plane.
ical bifurcation, the limit cycles have exchanged their stabilities perpendicular to the two-cluster subspace. Hence, the two-cluster limit cycle is now stable within the three-cluster subspace.
Analogous to the equivariant transcritical bifurcation of equilibria, in the full system the transcritical bifurcation of limit cycles is an equivariant transcritical bifurcation. As above, the cluster sizes are not true bifurcation parameters but rather enumerate different invariant subspaces in the same system. Thus, the transcritical bifurcation takes places simultaneously in all the different subspaces which are related to each other by index permutation, and in the corresponding subspaces which correspond to different cluster sizes n 2 and n 3 at fixed n 2 + n 3 = 1 − n 1 . In other words, in the full system with N identical components and sufficiently large N at y = 0.79625 multiple limit cycles from different threecluster subspaces meet in the same two cluster limit cycle, each changing the stability of a different Floquet exponent of the two-cluster limit cycle. Since all transversal eigenvalues (i.e., intra-cluster eigenvalues) of the Jacobian matrix are degenerate due to symmetry 1 , it actually affects all of them in the same way. As a consequence, in the full system, the equivariant transcritical bifurcation of two-cluster limit cycles stabilizes the latter and macroscopic two-cluster oscillations should be observable in arbitrarily large systems, i.e. even in the thermodynamic limit. We verified this conclusion numerically for a system with N = 20 units and the other parameters as in Fig. 7.

VI. DISCUSSION AND CONCLUSIONS
In the CO-electrooxidation model we analyzed in this paper, the fully synchronized one-cluster equilibrium loses stability in an S N -equivariant transcritical bifurcation where N − 1 real eigenvalues change sign simultaneously. Such bifurcations can arise for example as a perturbation of a pitchfork bifurcation 21 where the invariant subspace persists. Indeed, the pitchfork bifurcation for identical cluster sizes (n 1 = n 2 ) perturbs to a pair of saddle node and transcritical bifurcation for nonidentical cluster sizes. While the symmetric transcritical bifurcation happens simultaneously for all cluster sizes, the parameter values of the saddle node bifurcation are distinct for different relative cluster sizes as seen in Fig. 4: The most asymmetric two-cluster state is generated first, then increasingly symmetric ones, and finally the symmetric ( 1 2 , 1 2 )cluster in a pitchfork bifurcation. The symmetry breaking bifurcations from two-to three-cluster equilibria occur in a similar fashion, just that there is an additional cluster as an 'bystander'. More generally, if the phase-space geometry allows, S N -equivariant transcritical bifurcation can also have a global flavor 24 .
Equations that describe globally-coupled identical onedimensional units also directly relate to more general contexts, such as coupled oscillators. If a general S N -equivariant system undergoes a S N -equivariant transcritical bifurcation, one can describe the dynamics on the corresponding center manifold through a suitable normal form. This approach has been used to describe the dynamics of identical Stuart-Landau oscillators at a symmetry breaking bifurcation [25][26][27] . Truncating to appropriate order, one obtains a cubic S N -equivariant system of coupled one-dimensional variables x k 20,26-30 that describe the dynamics on the center manifold. In the limit of two large clusters and one cluster of vanishing size, the dynamics is characterized by a web of heteroclinic orbits between the two-cluster equilibria 27 . These cubic vector fields show the same bifurcation behavior close to the symmetry breaking bifurcations of the system analyzed here, including the transcritical and saddle node bifurcations, as well as the secondary bifurcations where the two-cluster solutions lose stability; cf. Refs. 25 and 26. A system of coupled onedimensional units with a cubic vector field, however, is unable to capture some of the more intricate secondary bifurcation behavior, such as Hopf and homoclinic bifurcations, described here. These require the inclusion of higher order terms, such as the coupling of (1) of our CO-electrooxidation model that is beyond cubic order, or a higher-order approximation on the center manifold.
The dynamical features described in this paper arise in a mathematical model that captures the experimental control of a physical system. From this perspective, the bifurcations discussed are common to many real-world systems. Yet, we made one simplification which will never be strictly fulfilled in a physical system, namely that all bistable units are identical. While the permutational symmetry is broken if the units are slightly heterogeneous, hyperbolic equilibria and limit cycles-which are the main focus of our analysis-will persist. The exact behavior of the branches at the (symmetric) bifurcation points will of course differ as these bifurcations break up into a series of generic bifurcation points; how exactly depends on the heterogeneity. However, the qualitative behavior that is relevant for real-world systems is preserved: In all ensemble steady states, each of the bistable units will be on one of its three branch segments. We can still define a cluster as a group of units that are on the same segment, and a unit transits between the segments at its turning points of u, i.e., u p→m or u m→a . However, for the individual units, the turning points do not necessarily occur at identical values of u as the branches deform. We verified that Hopf and homoclinic bifurcations still occur in a model of heterogeneous bistable units. Consequently, our bifurcation analysis facilitates a qualitative understanding of experimental observations also for nearlysymmetric systems and the conditions for the occurrence of oscillations: Consider, e.g., Li-ion batteries. Here, a global constraint, corresponding to the constant current in our model, is imposed on the total charge. The hysteresis in the cell voltage U over the total charge which is measured when charging and discharging Li-ion batteries with phase transition cathodes even at infinitely slow rates [14] could be seen as the 'shadow' of the transcritical bifurcations at the turning points of the charge vs voltage curve of the individual nano-sized storage particles in the symmetric system.
Another physical system that shows similar symmetrybreaking bifurcations, including a Hopf bifurcation, are certain electronic circuits, which, in addition, can be considered as superlattice models 31 . While in that case the governing equations of the individual elements were twodimensional, our studied revealed that stable, macroscopic collective oscillations may also exist in ensembles of coupled one-dimensional bistable units. This might lead to a reinterpretation of the origin of oscillations observed in manyparticle systems, where macroscopic oscillations have been interpreted as the result of synchronization of oscillating units.

DATA AVAILABILITY STATEMENT
The data that support the findings of this study are available from the corresponding author upon reasonable request.