Cluster synchronization on hypergraphs
CCluster synchronization on hypergraphs
Anastasiya Salova ∗ Department of Physics and Astronomy, University of California, Davis, CA 95616, USA andComplexity Sciences Center, University of California, Davis, CA 95616, USA
Raissa M. D’Souza
Complexity Sciences Center, University of California, Davis, CA 95616, USADepartment of Computer Science and Department of Mechanical and Aerospace Engineering,University of California, Davis, CA 95616, USA andSanta Fe Institute, Santa Fe, NM 87501, USA
We study cluster synchronization (a type of synchronization where different groups of oscillators inthe system follow distinct synchronized trajectories) on hypergraphs, where hyperedges correspondto higher order interactions between the nodes. Specifically, we focus on how to determine admissiblesynchronization patterns from the hypergraph structure by clustering its nodes based on the inputthey receive from the rest of the system, and how the hypergraph structure together with the patternof cluster synchronization can be used to simplify the stability analysis. We formulate our resultsin terms of external equitable partitions but show how symmetry considerations can also be used.In both cases, our analysis requires considering the partitions of hyperedges into edge clusters thatare induced by the node clusters. This formulation in terms of node and edge clusters providesa general way to organize the analysis of dynamical processes on hypergraphs. Our analysis hereenables the study of detailed patterns of synchronization on hypergraphs beyond full synchronizationand extends the analysis of cluster synchronization to beyond purely dyadic interactions.
I. INTRODUCTION
Patterns of synchronization in complex interdepen-dent dynamical systems are important to their func-tion. Often, such systems are modeled by networks withdyadic interactions between agents [1]. However, study-ing dyadic interactions is not always sufficient. Higherorder edges may be required to describe many systems,including chemical [2], biological [3], and coauthorshipinteractions [4, 5], and processes such as consensus dy-namics [6], making it necessary to go beyond the pairwiseinteraction analysis [7]. Structurally, we refer to suchsystems as hypergraphs (also known in the literature ashigher order networks or as simplicial complexes if extraconditions on hyperedges are satisfied [8]).Cluster synchronization is a type of synchronizationwhere different groups of oscillators in the system followdistinct synchronized trajectories. Studies of cluster syn-chronization have focused on networks with dyadic inter-actions, including a broad class of important behaviorssuch as remote synchronization and chimera states [9, 10]with wide areas of applicability from neuroscience andecological networks to opinion dynamics [11–15]. Ideasfrom graph and equivariant dynamical systems theorycan be applied to deduce admissible patterns of clustersynchronization on dyadic networks and simplify theirstability analysis [16, 17]. A more general question offinding the patterns of synchronization in coupled cellnetworks including the possibility of non-additive inter-actions is considered in Ref.[18], although the stability ofthese patterns is not analyzed. ∗ [email protected] Recently, many contributions to analyzing full syn-chronization on hypergraphs and simplicial complexeshave been made. For instance, several recent works studyfull synchronization in phase oscillator models, wherehigher order interactions extend the range of availablemodels and lead to new behaviors [19–24]. Performingstability analysis is crucial to determine which synchro-nization patterns can be observed in experiments or natu-ral systems. Some recent works extend the master stabil-ity function formalism, originally formulated for dyadicinteractions, to higher order structures to analyze fullsynchronization and its stability [25–27]. Finally, themethod for analyzing the stability of cluster synchroniza-tion states on hypergraphs for the special case of non-intertwined clusters is presented in Ref. [28]. However,the admissibility and stability of cluster synchronizationstates in arbitrary systems of oscillators coupled via hy-pergraphs has not yet been considered.In this manuscript, we formulate the conditions forcluster synchronization on hypergraph structures usingexternal equitable partitions. We formulate these con-ditions in terms of the node cluster partitions, but wealso need to consider how those node clusters partitionthe hyperedges into “edge clusters” for each order of theinteractions (e.g., dyadic edge clusters, triadic edge clus-ters, etc.). This allows us to simplify the stability calcu-lation by grouping the contributions to node dynamicsvia each of the distinct edge clusters. We extend theresults from Refs. [17, 29] (originally formulated on net-works with dyadic interactions) to a hypergraph settingand show how to block diagonalize the Jacobian to sim-plify the stability analysis beyond dyadic interactions andnon-intertwined clusters.The projection of a hypergraph onto a dyadic network a r X i v : . [ n li n . AO ] J a n with multiple edge types is a useful tool for studying com-plete synchronization on hypergraphs [25, 30]. However,we demonstrate that considering this reduced, projectedstructure is not generally sufficient to analyze the admis-sibility and stability of cluster synchronization patternson hypergraphs. This highlights the importance of for-mulating the analysis in terms of the hyperedges them-selves.The rest of the manuscript is organized as follows. Weintroduce relevant notation as well as the formal conceptof cluster synchronization in Section II. In Section III,we study cluster synchronization emerging from exter-nal equitable partitions. We show how stability analysiscan be simplified by finding a Jacobian block diagonaliza-tion in Section IV. The manuscript focuses on undirectedhypergraphs with Laplacian-like coupling. Appendix Aoutlines how similar analysis can be performed for othertypes of coupling functions. In Appendix B, we considerhow cluster synchronization arises from permutationalsymmetries of the incidence matrices (or, equivalently,adjacency tensors) corresponding to each hyperedge or-der. We show how group representation theory can beused to simplify the stability analysis in Appendix C. II. SETUP
First, we define the general form of the dynamics onhypergraphs that is being considered. A hypergraph isdefined by a set of N nodes and a set of hyperedges e j ∈E . In this work, we focus on undirected hyperedges. Let E i ⊂ E be the set of hyperedges that contain node i .Each hyperedge e j ∈ E i is characterized by a set of nodes e j = { i, j , ..., j m − } . The order of the hyperedge e j is m , which is the number of nodes including i that are partof it. Thus, m = 2 corresponds to dyadic edges, m = 3to triadic edges, etc.Using notation similar to Ref.[25], we can express theevolution of the state of each node in the system, x i ∈ R n ,as: ˙ x i = F i ( x i ) + (cid:88) e ∈E i G j ( x i , x e \ i ) . (1)Here, the function F i ( x i ) describes the evolution of un-coupled nodes, and the function G j ( x i , x e \ i ) is a couplingfunction corresponding to the influence of the hyperedge e on node i , where x i is the state of the node i itself, and x e \ i is the state of the rest of the edge. This setup is gen-eral, including the case when the interaction hypergraphis a simplicial complex and the additional requirementthat each subset of hyperedge nodes forms a hyperedgemust be satisfied.Often, some degree of homogeneity is present withinthe nodal dynamics, F i ( x i ), of different nodes i as wellas in the coupling dynamics, G j ( x i , x e \ i ). In that case,one can use the hypergraph structure to find nontrivialpartitions into sets of nodes that can fully synchronize. Inthe simplest case, all the self-dynamics are characterized by the same function F and the coupling dynamics ofa given order m are characterized by the same function G ( m ) . In that case, it is sufficient to consider adjacencystructures (e.g., adjacency tensors) with binary entries.For each edge order m , the adjacency structure canbe defined in terms of the collection of m incidence ma-trices I ( m ) , or, equivalently, adjacency tensors A ( m ) (asillustrated in Fig.3 of Ref.[7]). Let E ( m ) i be the set ofhyperedges of order m containing the node i . Then, thenonzero elements of the incidence matrix are [ I ( m ) ] i,e = 1if e ∈ E ( m ) i .On the other hand, the adjacency structure can berepresented as the adjacency tensor A ( m ) for each or-der m , where [ A ( m ) ] i,j ,...,j m − = 1 if the hyperedge e = i, j , ..., j m − , is present. These ways of representingthe adjacency structure are equivalent, and which one isused depends on which is most convenient [7]. Namely,[ I ( m ) ] i,e = [ A ( m ) ] i,j ,j ,...,j m − if e = { i, j , ..., j m − } .In this manuscript, we focus on noninvasive couplingfunctions. Specifically, we assume that the coupling func-tion is Laplacian for dyadic interactions and Laplacian-like for higher-order interactions, (and will use the term“Laplacian-like” in general). Laplacian coupling fordyadic interactions is of the form G ( x i , x j ) = H ( x j ) − H ( x i ). Laplacian-like coupling for edges of order m is ofthe form G ( m − (cid:80) l =1 x j l − ( m − x i ). Coupling functions ofthis form are natural, for instance, for higher order net-works of phase oscillators [23]. Additionally, we requirethe hypergraph to be undirected. Namely, if the m th or-der hyperedge A ( m ) i,j ,...,j m − = 1 is present, we require thatto be the case for all the index permutations (individualpermutations denoted by P ): A ( m ) i,j ,...,j m − = A ( m ) P ( i,j ,...,j m − ) . (2)Thus, the edges can be defined as unordered sets of el-ements, as each node in the edge receives input of thesame form from all the other nodes.Rewriting Eq. (1) for homogeneous node dynamics andhomogeneous edge dynamics for each hyperedge of order m leads to the dynamical evolution equation:˙ x i = F ( x i ) + d (cid:88) m =2 σ ( m ) (cid:88) e ∈E ( m ) [ I ( m ) ] i,e G ( m ) ( x i , x e \ i )= F ( x i ) + d (cid:88) m =2 σ ( m ) ( m − N (cid:88) j =1 ... N (cid:88) j m − =1 [ A ( m ) ] i,j ,...,j m − · G ( m ) ( x i , ..., x j m − ) , (3)where d is the maximum edge order present in the hyper-graph. Here, σ ( m ) denotes the strength of the m th ordercoupling.Cluster synchronization is manifested in groups ofnodes following the same trajectory over time, x i ( t ) = ... = x i L ( t ), where the groups are not fully synchro-nized with one another. We call each group of syn-chronized nodes a “node cluster” and denote them by (a) (b) (c) v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v FIG. 1. Examples of admissible cluster synchronization patterns on an example hypergraph. (a): two clusters. (b) Fourclusters, (c) Eight clusters. The nodes are colored according to their node clusters. C , C . . . , C K assuming K distinct groups exist. Ex-amples of node clusters are given in Fig. 1. Theset of dynamic trajectories followed by the nodes ineach cluster, the “node cluster trajectories”, can be ex-pressed as s ( t ) , ..., s K ( t ). The trajectories are time-dependent, but we assume the time-dependence is im-plicit and use the notation s , ..., s K for compactnessin the rest of the manuscript. Likewise, we consider“edge clusters” and “edge cluster trajectories”. A hy-peredge of order m can be characterized by the nodeclusters to which the m nodes it connects together be-long. All the hyperedges of order m with the same setof node clusters is called an edge cluster and denoted by C ( m )1 , C ( m )2 , . . . , C ( m ) K m . The edge cluster trajectories aredenoted by s C ( m )1 , s C ( m )2 , . . . , s C ( m ) Km , where s C ( m ) j is the setof dynamic trajectories followed by the nodes involved inthe j th edge cluster. A concrete example of node andedge clusters will be given in Fig. 2. The clusters withtheir corresponding trajectories will be used to facilitatestability calculations.In the following section, we show how to use the prop-erties of the incidence matrices I ( m ) (or adjacency andLaplacian tensors A ( m ) and L ( m ) ) to determine the ad-missible patterns of cluster synchronization. Then, inSec. IV, we develop the stability analysis of these pat-terns. III. BALANCED EQUIVALENCE RELATIONSAND SYNCHRONIZATION
Patterns of synchronization on hypergraphs can be de-duced from balanced equivalence relations (equitable andexternal equitable partitions), as we will demonstrate inSection III A. In Sec. III B we show the distinction be-tween admissible states on a hypergraph and on its lower-dimensional projection onto a network with multiple edgetypes. Yet, even if a pattern of cluster synchronizationis dynamically admissible, it is not guaranteed to be ob-served in biological or engineered systems for a given pa-rameter regime, or their simulations. Thus, determin- ing the stability of cluster synchronization is a crucialstep that can help relate models to real world phenom-ena. The stability calculation can become less tractableas the system size increases, thus making dimensional-ity reduction very useful. In Section IV A, we reviewstability calculations for networks of dyadic interactions.In Section IV B, we generalize these results to dynami-cal systems on hypergraphs. We also demonstrate how togeneralize the results in Ref. [29] to block diagonalize theJacobian matrix in the context of stability calculationson undirected hypergraphs.
A. Patterns of synchronization from externalequitable partition
For dynamical systems on networks with purely dyadicinteractions, equitable partitions can be used to deter-mine the synchronized clusters [31, 32] as well as otherpatterns of synchronization [33]. Equitable partitions di-vide the network into cells, where each node in a cell C i receives the same input from any cell C j including thenodes within its own cell, i = j . Each cell of the networkdefines a cluster of nodes that could be synchronized.In case of noninvasive coupling, as is the focus of thismanuscript, the conditions above only have to hold for i (cid:54) = j (in which case the partition is called external equi-table partition), since the terms representing the effect ofnodes within the same cluster upon one another becomeszero if evaluated on that state.The same idea holds for higher order interaction net-works, where the equitable partitions now need to bedefined in terms of the interactions of all orders. Theconditions can be relaxed to external equitable partitionif the coupling functions are noninvasive, which is thecase for the Laplacian-like coupling considered herein.More specifically, this can be understood in terms ofpartitioning the incidence matrix. The nodes can be sep-arated into non-overlapping cells of “node clusters”. Thisnode partition induces a partition of edges into “edgeclusters”, according to what combination of node clus- (a) (b)(c) v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v bv gy (d) (e) [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , , ] [ , , ] [ , , ] [ , , ] [ , , ] [ , , ] FIG. 2. An admissible pattern of synchronization into fourclusters of nodes based on an external equitable partition.(a) The hypergraph. (b) The 1st order quotient network. (c)The 2nd order quotient network, with nodes undergoing thateffective dynamics shown in larger size. [(d)-(e)] Incidencematrices for dyadic and triadic interactions respectively. Dotsrepresent ones. Row label colors represent the node clusters,column label colors represent the edge clusters induced bythe node clusters. As shown in (d) there are 6 types of dyadicedge clusters and in (e) there are two types of triadic edgeclusters. ters those edges contain (we only need to consider theunordered set of included node clusters in the case ofLaplacian-like coupling as the edges are undirected). Thepartition is equitable if each node in a given node clustergets the same input from each edge cluster. We demon-strate this on a concrete example by considering the hy-pergraph structure shown in Figs. 1 and 2.The structure of the hypergraph in Figs. 1 and 2 isan extension of the network in Fig.1 in Ref. [29], withextra hyperedges added to represent the higher order in-teractions, and extra edges added to highlight that strictsymmetry conditions are not necessary for our frame-work. We need to specify certain aspects of the dynamicsfor the system as well and for this we impose conditionson the coupling functions. Namely, we assume both G (2) and G (3) are both Laplacian-like (and therefore permu-tationally symmetric and noninvasive).In general, many distinct partitions, each one corre-sponding to a different pattern of cluster synchronization,can be admissible for a given hypergraph (e.g., three dif-ferent patterns of cluster synchronization are presentedin Fig. 1, with two, four and eight clusters respectively). To illustrate the partition into node and edge clus-ters, we focus on the four cluster synchronization pat-tern shown both in figure Fig. 1(b) and Fig. 2(a). Wecan divide the nodes into four non-overlapping cells (i.e.,synchronized clusters) which we label by their numberfor convenience in mathematical formulas, C , C , C , C ,or equivalently by their color for convenience when re-ferring to the figure, C g , C y , C b , C v , corresponding togreen, yellow, black, and violet. C = C g = { , , } , C = C y = { , , } , C = C b = { , , , , , } ,and C = C v = { , , } . This partition of nodesinto node clusters induces a partition of edges into edgeclusters based on the node clusters those edges span.The 6 distinct dyadic order edge clusters that exists areshown by their colors in the column labels of Fig. 2(d).Specifically, the edge clusters are C (2)1 = C (2) gg = { [1 , , [2 , , [3 , } , C (2)2 = C (2) gy = { [1 , , [2 , , [3 , } , C (2)3 = C (2) yb = { [4 , , [4 , , [5 , , [5 , , [6 , , [6 , } , C (2)4 = C (2) yv = { [4 , , [5 , , [6 , } , C (2)5 = C (2) bb = { [7 , , [8 , } , and C (2)6 = C (2) vv = { [13 , } . Thetwo distinct triadic order edge clusters are shown bytheir colors in the column labels of Fig. 2(e), where C (3)1 = C (3) gyv = { [1 , , , [2 , , , [3 , , } and C (3)2 = C (3) ybb = { [4 , , , [5 , , , [6 , , } . These cells forman external equitable partition, meaning each node in agiven node cluster gets the same external input as all theother nodes in that cluster. Therefore, this partition intonode clusters corresponds to an admissible state. Addi-tionally, we use the notation C j ∈ C ( m ) k to denote thenode clusters that are included in the k th edge cluster,for instance C (3)1 = { C g , C y , C v } = { C , C , C } , writtenrespectively in terms of colors and then index number.The condition for cluster synchronization in networkswith Laplacian-like coupling can be written as: (cid:88) j ∈ C ( m ) k I ( m ) ij = (cid:88) j ∈ C ( m ) k I ( m ) i (cid:48) j , (4)for i, i (cid:48) ∈ C l , where we are summing over the m -th orderhyperedges j that are in edge cluster C ( m ) k , and the termscoming from edge clusters C ( m ) k containing only nodes in C l can be ignored.The effective dynamics of cluster synchronized stateswith dyadic interactions can be described via the use ofquotient networks. Analogously, we can form a quotienthypergraph (which is a hypergraph that represents ef-fective interactions between nodes of different clusters)based on the partition above. The structure of that effec-tive quotient hypergraph for the state shown in Fig. 2(a)is contained in the effective incidence matrices, I ( m )eff , de-fined as: I (2)eff = [ b y ] [ y v ] [ y g ] b y v g , I (3)eff = [ bb y ] [ g y v ] b y v g , (5)and obtained from: I ( m )eff = P n I ( m ) ( P ( m ) e ) T , (6)where P n ( K × N ) and P ( m ) e ( K m × N ) are the node andedge projection matrices respectively, and I ( m ) is the m thorder incidence matrix. The projection matrix is definedby [ P n ] i,j = 1 if node i belongs to node cluster C j , and[ P ( m ) e ] i,j = 1 if the m th order edge i belongs to the m thorder edge cluster C ( m ) j . Here, we take into account thatthe dyadic coupling is noninvasive by ignoring the inter-action terms between the fully synchronized nodes (e.g.,nodes 8 and 11 in Fig. 2 (a)), and the undirectedness ofthe triadic coupling by considering, e.g., [ byb ] and [ ybb ]to be the same edge pattern.The quotient hypergraph (illustrated in Fig. 2 (b-c))can be used to read out the time evolution of each node.For instance, every node in the yellow cluster evolvesaccording to:˙ x y = F ( x y ) + G (2) ( x y − x g )+ G (2) ( x y − x v ) + 2 G (2) ( x y − x b )+ G (3) (2 x y − x b − x b ) + G (3) (2 x y − x v − x g ) , (7)with analogous equations describing the time evolution ofcompletely synchronized nodes belonging to other clus-ters. Here, G ( m ) is expressed taking into account theLaplacian-like coupling assumption.The admissibility and stability properties of variouscluster synchronization patterns can also be formulatedin terms of the adjacency tensor. To achieve that, weneed to introduce additional pieces of notation. Let { E , ..., E K } be the diagonal matrices corresponding todifferent clusters (assuming K distinct synchronized clus-ters are present). These matrices have zero off-diagonalentries and ones on the positions on the diagonal corre-sponding to [ E k ] ii = 1 if i ∈ C k . Similarly, we can define N × ... × N (cid:124) (cid:123)(cid:122) (cid:125) m times tensors corresponding to different edge clus-ters, denoted by E C ( m ) k , where [ E C ( m ) k ] ij ,...,j m − = 1 if[ i, j , .., j m − ] ∈ C ( m ) k , and elements equal to zero oth-erwise. For mathematical convenience it is necessary attimes to specify that a particular cluster be indexed firstin the tensor. Thus, we use the notation E C ( m ) k,l to indi-cate that the l th cluster is being indexed first.Using this notation, the fact that the m th order con-tributions to the dynamical evolution of nodes i and i (cid:48) in the same cluster, ( i, i (cid:48) ∈ C l ), are equal is expressed via:[ A ( m ) ◦ E C ( m ) k,l ] ij ...j m − j ...j m − = [ A ( m ) ◦ E C ( m ) k,l ] i (cid:48) j ...j m − j ...j m − , (8)where stands for a tensor of all ones, and the tensorsummation notation is used.If the conditions of Eq. (8) (or, equivalently, Eq. (4))hold for all node clusters, the synchronization patternis valid for any coupling type in hypergraphs includingdirected edges (similar to equitable partitions of networkswith dyadic interactions). The conditions can be relaxedfor different kinds of noninvasive couplings. Specifically,for Laplacian and Laplacian-like coupling, the conditionscan be relaxed by ignoring the input of edges containingonly the nodes in the same cluster.Here, we discussed the cluster synchronization condi-tions for hypergraphs from the perspective of the inci-dence matrix as well as the adjacency tensor. These con-ditions are mathematically equivalent, but have differentbenefits in terms of visualization and computation. InSection III B, we discuss how these conditions are relatedto balanced equivalence relations on a reduced-order pro-jected network with purely dyadic interactions. B. Relation to the projected network with multipleedge types
Finding the synchronized clusters of the hypergraphcan be assisted by considering its projection onto anetwork with multiple types of edges. Algorithms tofind equitable partitions of dyadic graphs are available[31, 34, 35]. The algorithm of Ref.[31], for instance, canbe applied to directed and undirected networks with oneor more edge types. Methods for finding such balancedrelations for hypergraphs are less explored in synchro-nization literature. Thus, it can be useful to consider pro-jected hypergraphs amenable to previously establishedmethods to aid in searching for potential admissible clus-ter patterns.We define a projected network, in which each hyper-edge of order m is replaced by a clique of m nodes withall-to-all coupling [25]. We additionally assign a separateedge weight α m to cliques arising from each hyperedgeorder. The elements of the projected adjacency matrix A P can be written as:[ A P ] ij = d (cid:88) m =2 (cid:88) e ∈E ( m ) i ,i,j ∈ e α m , (9)or, equivalently:[ A P ] ij = d (cid:88) m =2 α m A ijj ...j m − j ...j m − . (10) v v v v v v v v v v v v v v v v v v v v v v v v ! admissible not admissible FIG. 3. Example of a pattern not admissible for the hyper-graph (subfig. (a)), but admissible for its projected network(subfig. (b)).
Every equitable partition of a hypergraph is an eq-uitable partition of the projected network, since dyadicinteractions can be considered a subset of polyadic in-teractions that are additive. On the other hand, not ev-ery pattern of cluster synchronization admissible on theprojected graph is an admissible cluster synchronizationpattern of the original hypergraph. For instance, considerthe example in Fig. 3, where the pattern of synchroniza-tion with nodes 1 and 6 synchronized is admissible on theprojected network, Fig. 3(b), but not the hypergraph,Fig. 3(a). Therefore, after obtaining admissible synchro-nization patterns on a projected graph, it is necessary totest if conditions in Section III A are satisfied to deter-mine the admissibility on the hypergraph provided we donot take into account various types of noninvasiveness.An additional issue arises when the coupling functionsare noninvasive, i.e. when the coupling function vanishesif all (or a fraction of) nodes on a hyperedge are fullysynchronized. That can be easily taken into account inmodifying Eq. (8), but is nontrivial to take into accountin a partitioning algorithm for cluster synchronization onthe projected network.We also note that unlike the case of full synchro-nization in Ref.[25], analyzing stability of synchronizedclusters requires explicitly considering structures beyondthe projected adjacency matrix. The details of that areshown in Section IV B.
IV. STABILITY OF CLUSTERSYNCHRONIZATION PATTERNS
To illustrate the hypergraph cluster synchronizationstability formalism, we first briefly summarize existing re-sults for systems with dyadic interactions in Section IV A.
A. Stability of synchronization patterns innetworks with dyadic interactions
For networks with dyadic interactions, the dimension-ality of the stability calculation can be reduced by blockdiagonalizing the full system’s Jacobian matrix and de-termining the maximum transverse Lyapunov exponentsseparately for each lower dimensional transverse block.The simplification can be achieved both from the sym-metry perspective [12, 17, 36] and from balanced equiv-alence relations [29, 37]. These methods of reducing thedimensionality of the problem can be modified to be ap-plicable to networks with higher order interactions. Inthis section, we focus on doing so for the more generalcase of cluster synchronization arising from the equitablepartition.First, we present the variational equation that deter-mines the linear stability of the cluster synchronized statefor dyadic interactions. Recall the notation that s k isthe trajectory of the k th node cluster and s C ( m ) k is thetrajectory of the k th edge cluster for edges of order m .Likewise, we use K to denote the total number of nodeclusters and K m to denote the total number of edge clus-ters of order m . The variational equation for Laplaciancoupling (where the coupling function can be expressedas G ( x i , x j ) = H ( x j ) − H ( x i )) is: δ ˙ x = (cid:16) K (cid:88) k =1 E k ⊗ JF ( s k )+ σ K (cid:88) k =1 LE k ⊗ JH ( s k ) (cid:17) δx, (11)where k indexes over all the clusters, J is the Jacobianoperator, and ⊗ represents the tensor product. Addi-tionally, L = A − D is the Laplacian matrix, and thedegree matrix D is a diagonal matrix with elements onthe diagonal calculated as D ii = (cid:80) k A ik . The simulta-neous block diagonalization of the terms in Eq. (11) canbe performed using the method from Ref.[29] that relieson simultaneously block diagonalizing the set of matri-ces { E , ..., E K , L } (or { E , ..., E K , A } for adjacency cou-pling).The Laplacian case can be generalized to G ( x i , x j ) = G ( x i − x j ), the form that is more similar to higher orderLaplacian-like coupling. Then, the variational equationsare: δ ˙ x = (cid:16) K (cid:88) k =1 E k ⊗ JF ( s k )+ σ K (cid:88) k =1 E k K (cid:88) l =1 (cid:0) AE l − D l (cid:1) ⊗ JG ( s k , s l ) (cid:17) δx, (12)where the diagonal matrix [ D l ] ii = (cid:80) j [ AE l ] ij indicatesthe number of edges from cluster C l into each node i ,and: JG ( s k , s l ) p,q = ∂G p ( x i , x j ) ∂ [ x i ] q (cid:12)(cid:12)(cid:12)(cid:12) x i = s k ,x j = s l , (13)where p and q index over the dimensions of the statespace of x i . We denote L l = AE l − D l . Then, the set ofmatrices that needs to be simultaneously block diagonal-ized is { E , ..., E K , L , ..., L K } .In the following subsection, we extend this to systemswith higher order interactions. B. Stability of synchronization patterns for higherorder interactions
Eq. (12) acquires additional terms in presence of higherorder coupling. Here, we discuss how all the added termsin the Jacobian can be diagonalized for a given clus-ter synchronization pattern using the incidence matri-ces, or, alternatively, adjacency tensors, for the case ofLaplacian-like coupling.First, we define a Laplacian corresponding to a specificedge synchronization pattern (i.e., if there are multipleadmissible m th order edge synchronization patterns, weconsider the k th such one) as follows: A ( m ) k = I ( m ) k [ I ( m ) k ] T − D ( m ) k , L ( m ) k = A ( m ) k − m D ( m ) k , (14)where I ( m ) k is an N × | C ( m ) k | (here, | C ( m ) k | denotes thenumber of unique elements in the edge cluster C ( m ) k ) ma-trix obtained by stacking the columns of I ( m ) correspond-ing to the hyperedges in the k th cluster of order m . Forinstance, for the bby I ( m ) k is obtained by keeping thelast 3 columns. Additionally, D ( m ) k is a diagonal matrixof node degrees corresponding to the number of edgeswith that synchronization pattern.Then, the variational equations evaluated at that clus-ter synchronized state are: δ ˙ x = (cid:18) K (cid:88) k =1 E k ⊗ JF ( s k ) + d (cid:88) m =2 σ ( m ) · (15) (cid:16) K m (cid:88) k =1 (cid:88) l ∈{ C ( m ) k } E l L ( m ) k ⊗ JG ( m ) ( s l , s C ( m ) k \ l ) (cid:17)(cid:19) δx, where { C ( m ) k } is a set of unique node clusters includedin the k th edge cluster, (e.g., in Fig. 2, { C (3) ybb } = { y, b } ).Additionally, s C ( m ) k \ l is set of all the trajectories of nodesincluded in edge cluster C ( m ) k , excluding those nodes in node cluster l . The partial derivatives are computed as: JG ( m ) ( s l , s C ( m ) k \ l ) p,q = ∂G ( m ) p (cid:80) j ∈ C ( m ) k \ l x j − ( m − x l ∂ [ x j ] q (cid:12)(cid:12)(cid:12)(cid:12) x j = s j ,x l = s l = ∂G p ( z ) ( m ) ∂z q (cid:12)(cid:12)(cid:12)(cid:12) z = (cid:80) j ∈ C ( m ) k \ l s j − ( m − s l . (16)Thus, to block diagonalize the entire system of Laplacian-like coupled oscillators, it is sufficient to simultaneouslyblock diagonalize the following matrices: { E , ..., E M , L (2)1 , ..., L (2) K , ..., L ( d )1 , ..., L ( d ) K d } . (17)The form of the Laplacians describing a specific patternof cluster synchronization is similar to that of the gen-eralized Laplacian [38]. It is sufficient to block diago-nalize the generalized Laplacian if the clusters are notintertwined [28]. It is also sufficient in the case of hyper-graphs where each node is part of a single type of hyper-edge. However, in the most general case, it is necessaryto consider the entire set in Eq. (17).Alternatively, the set of matrices that needs to be si-multaneously block diagonalized can be defined in termsof the adjacency tensor: E l A ( m ) k = [ A ( m ) ◦ E C ( m ) k ] ij ...j m − j ...j m − ,E l L ( m ) k = E l A ( m ) k − ( m − D ( m ) k , D ( m ) k = n (cid:88) j =1 [ E l A ( m ) k ] ij . (18)The following set of matrices needs to be simultaneouslyblock diagonalized: { E , ..., E M } ∪ M (cid:91) m =1 K m (cid:91) k =1 (cid:91) l ∈{ C ( m ) k } { E l L ( m ) k } . (19)This block diagonalization produces the same block di-agonalization of the Jacobian as Eq. (17).An example of simultaneous block diagonalization isshown in Fig. 4 with the upper panels correspondingto the two node cluster case shown in Fig. 1(a) andthe lower panels to the four node cluster case shown inFig. 1(b). Additionally, we impose Laplacian couplingon dyadic edges, and Laplacian-like coupling on triadicedges to match the example we present in Section IV C.We observe that in the first case (two clusters, shown onthe upper panel), it is sufficient to simultaneously blockdiagonalize the cluster indicator matrices ( E y and E b )and the generalized Laplacian, since each node partici-pates in a unique triadic edge pattern. However, in thesecond case, it is necessary to simultaneously block di-agonalize two matrices arising from triadic interactions(shown in Fig. 4(e)). (a) (b) (c) E y E b A (2) A (3) ybb (d) (e) (f) E t E y E v E b A (2) A (3) ybb A (3) vyt FIG. 4. Block diagonalizing the Jacobian evaluated on thestates on Fig. 1. Subfigures (a)-(c) relate to Fig. 1(a) andsubfigures (d)-(f) relate to Fig. 1(b). [(a,b,d,e)] Matrices usedto construct the set of matrices that need to be simultaneouslyblock diagonalized. [(c,f)] Blocks of the resulting Jacobian.
C. Example: stability calculation for a specificdynamics and multiple edge types
The results presented in Section III of the manuscriptare easily generalizable to systems with different typesof nodes and edges, such as found in multilayered net-works, e.g., Ref. [12]. First, we assume that generally,only the nodes of the same type get fully synchronized.Additionally, we require that the input to each node inthe cluster from all types of edges is the same. Thus, weneed to cluster the edges not only based on the node clus-ters, but also on the edge types. Given the appropriatecluster assignment, the Jacobian block diagonalizationis guaranteed by simultaneously block diagonalizing theset of matrices in Eq. (17), where the numbers K m stillcorrespond to the number of edge clusters of order m .More specifically, checking for state admissibility canbe achieved by checking the condition in Eq. (4), wherethe incidence matrix entries are labeled according to theedge type. The condition stating that only the nodesof the same type can synchronize can be formulated interms of a trivial “first order incidence matrix” I (1) , alabeled N × F ( x i ) = β sin ( x i + π/
4) and coupling functions G (2) ( x i , x j ) = δ ij σ (2) (cid:2) sin ( x j + π/ − sin ( x i + π/ (cid:3) and G (3) = δ ijk σ (3) sin( x i + x j − x k ). Here, the values of δ ij and δ ijk are selected from values {− , } and correspond to at-tractive and repulsive hyperedges respectively. To avoidcomplications from multistability, we pick the state inFig. 1(a) for our analysis and we make edges connect-ing nodes that are in the same cluster attractive andall the other edges repulsive. Keeping the parameter β constant, we vary σ (2) and σ (3) to determine the lin-ear stability regions for different parameter regimes. Wepresent our findings for this example system in Fig. 5.Figure 5(a-e) shows the analogous plots to Fig. 2(a-e)with the state, the quotient network, and incidence ma-trices respectively. Figure 5(f) shows the set of matricesthat need to be simultaneously block diagonalized. Fig-ure 5(g) is the linear stability plot demonstrating sensi-tive dependence on both parameters σ (2) and σ (3) , withthe changes in the stability properties of the system show-ing correspondence to different regions of the bifurcationdiagram shown in Figure 5(h).Here, we demonstrated our approach in action. How-ever, further investigation is needed to determine the roleof higher order interactions in different types of dynami-cal systems on hypergraphs. V. CONCLUSION
Interdependencies in complex nonlinear dynamical sys-tems can often not be reduced to dyadic interactions.Some of these systems can be modeled as coupled oscil-lators on a hypergraph. Since there are many couplingfunction choices and oscillator dynamics types, this setupcan be appropriate in a wide variety of scenarios, andcan not be analyzed using only the methods applicableto systems with dyadic interactions.Full synchronization is one admissible state for cou-pled oscillators on hypergraphs, but so are more intricatesynchronization patterns where clusters of oscillators arefully synchronized, but distinct clusters follow differenttrajectories. We show how to use the structure of the in-cidence matrices, or, alternatively, the adjacency tensors,to determine the admissibility of cluster synchronizationpatterns and analyze their stability properties. To do so,we need to consider not only the partition into node clus-ters, but also the partition into hyperedge clusters that isinduced by the synchronization pattern of the entire setof nodes coupled on each hyperedge. This formulation interms of node and edge clusters provides a general wayto organize the analysis of dynamical processes on hyper-graphs including taking the partial derivatives requiredfor Jacobian analysis.A very important aspect of analyzing synchronizationpatterns is their stability analysis. Simultaneous blockdiagonalization can be performed to assist stability cal- (a) (b)(c) (d) (e) v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v y [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , ] [ , , ] [ , , ] [ , , ] [ , , ] [ , , ] [ , , ] [ , , ] (a) (b)(c) (d) (e)(f) (g) (h) E y E b A (2) bb A (2) yb A (3) ybb A (3) ybb FIG. 5. Example linear stability calculation for the system discussed in Section IV C. (a) A hypergraph with attractive (blue)and repulsive (violet) coupling and with nodes that obey the same dynamical equations. (b) The quotient network for thedyadic interactions. (c) The quotient network for the triadic interactions. (d) Structure of I (2) . Black and yellow represent thedistinct clusters, blue and violet represent the distinct coupling types. (e) Structure of I (3) . (f) Matrices used in simultaneousblock diagonalization to perform the stability analysis. (g) Linear stability diagram for various values of dyadic and triadiccoupling strengths, σ (2) and σ (3) . Pink areas are linearly stable, blue areas are not linearly stable. (h) Bifurcation diagramfor three distinct values of σ shown in the same color as their corresponding horizontal lines on Subfig.(g). Horizontal axisrepresents the dyadic coupling strength, vertical axis corresponds to the states of black nodes x black at the past 100 time steps.Background colors represent the calculated linear stability for each value of σ (2) . culations related to cluster synchronization on dyadicnetworks. We formulate similar conditions for systemswith higher order interactions by extending the set ofmatrices that need to be simultaneously block diagonal-ized. These conditions are more complex, especially incases when they do not depend on network symmetriesalone. Linear stability of full synchronization and clus-ter synchronization with non-intertwined clusters can beobtained using the generalized Laplacian [27, 28]. Un-like previous work [28], our analysis is not restricted tonon-intertwined clusters.The results presented in this manuscript open up anopportunity for detailed analysis of systems of theoreticaland practical significance. Our approach allows testingwhether a given cluster synchronization pattern is admis-sible based on the hypergraph structure and the type ofcoupling functions (e.g., Laplacian-like or noninvasive),and efficiently calculating the stability of specific pat-terns using simultaneous matrix block diagonalization.One of the directions for future work is developing ef-ficient algorithms for determining all the generic admis-sible patterns of cluster synchronization. While such al-gorithms exist for systems with dyadic interactions, theyhave limited applicability to higher order systems. Forinstance, some of the patterns that can be deduced from the projected dyadic network coincide with the admissi-ble patterns on the hypergraph, while others do not, asshown in Section III B. Moreover, reducing the systemto a projected graph loses information necessary to takenoninvasiveness into account. Similarly, only taking intoaccount the dyadic interactions to determine the candi-date patterns and then checking their admissibility inpresence of higher order interactions is another possibleapproach to finding a full set of cluster synchronizationpatterns. However, it can be very inefficient. For in-stance, in cases where the network of dyadic interactionshas a very large number of admissible patterns (suchas those arising from all-to-all coupling), this step failsto significantly narrow down the set of candidate syn-chronization behaviors. Additionally, determining therole of higher order interaction structure in stabilizingand destabilizing the synchronization patterns is anotherpromising direction for future work. ACKNOWLEDGMENTS
The authors would like to thank Yuanzhao Zhang andAdilson Motter for helpful discussions about theory andcode.0
Appendix A: Other types of coupling functions
The manuscript focuses on cluster synchronization onhypergraphs with Laplacian and Laplacian-like coupling.However, the methods we describe in Sections III and IVof the manuscript can be applied to a more general setof coupling types with some modifications. One of thesecoupling types is adjacency-like coupling, where the cou-pling functions G ( m ) arising in dynamical equations donot directly depend on the state x i of node i receivingthe coupling function’s input. The other type of cou-pling is noninvasive coupling, where the coupling func-tion G ( m ) on a hyperedge vanishes if all (or any fractionof) the nodes of the hyperedge are fully synchronized[27]. Finally, in all of these cases, the coupling can bedirected, which is distinct from the case we discussed inthe manuscript.The case of undirected adjacency-like coupling is verysimilar to Laplacian-like coupling. However, cluster syn-chronization patterns arise from equitable partitions inthis case. This means that to determine the admissibil-ity of a specific synchronization pattern, the in-clusteredges need to be taken into account in addition to theedges between distinct clusters in Eq. (4). For each state,the Jacobian can be block diagonalized by simultaneouslyblock diagonalizing the following matrices: { E , ..., E M , A (2)1 , ..., A (2) K , ..., A ( d )1 , ..., A ( d ) K d } , (A1)defined in Section IV B.Noninvasive coupling generalizes Laplacian-likecoupling. Cluster synchronization patterns forLaplacian-like coupling arise from partitions simi-lar to external equitable partitions. For instance, if G ( m ) ( x i , x j , ..., x j m − ) = m − (cid:81) k =1 ( x j k − x i ), the contributionof coupling functions on hyperedges where several nodesbelong to the same cluster to each of these nodesvanishes. Thus, some of the hyperedge contributions canbe ignored even if not all the nodes on these hyperedgesare fully synchronized.So far, we have considered undirected coupling. Incase of undirected coupling, the presence of the hy-peredge { i , ..., i k , ..., i m } providing input to node i via the coupling function G ( m ) , s.t. x i = ... + G ( m ) ( x i , ..., x i k , ..., x i m ), implies that hyperedge affects x i k via the same coupling function, s.t. x i k = ... + G ( m ) ( x i k , ..., x i , ..., x i m ). Additionally, the couplingfunction has to be invariant with respect to permu-tations of the elements corresponding to the input-providing node, namely, G ( m ) ( x i , ..., x i k , ..., x i k (cid:48) , ... ) = G ( m ) ( x i , ..., x i (cid:48) k , ..., x i k , ... ). If some of these conditionsare violated, the hypergraph becomes directed. The par-tition of directed hyperedges now has to take the orderingof the nodes of the hyperedge into account, as the hyper-edges containing the nodes belonging to the same clus-ters are not necessarily dynamically equivalent under thedirected coupling conditions. However, the equitable or external equitable partitions still provide cluster synchro-nization patterns under the appropriate edge clustering.Block diagonalizing the Jacobian for directed hyper-graphs is more challenging. For instance, the algorithmpresented in [28, 29] is only applicable to the undirectedcase. However, in a more restrictive case where the clus-ter synchronization state arises from symmetries (orbitalpartition of the hypergraph), the symmetry based anal-ysis presented in Appendices B and C is directly appli-cable. Appendix B: Patterns of synchronization fromhypergraph symmetries (orbital partitions)
Some of the patterns of synchronization in networks ofidentical oscillators can often be deduced from the sym-metries of the system [36, 40]. Partitions based on thesesymmetries, referred to as orbital partitions, constitutea special case of equitable partitions. More specifically,cluster synchronization in networks can occur as a resultof permutational symmetries, contained in the automor-phism group of the network adjacency matrix. Recently,symmetry has been used to obtain the patterns of clustersynchronization in complex networks, including multi-layer networks, and simplify the stability calculations forthese systems using group representation theory [12, 17].Here, we extend these admissibility and stability resultsto use symmetry considerations to study patterns of clus-ter synchronization on hypergraphs.First, we state the results for systems with dyadic in-teractions. The automorphism group of the adjacencymatrix A is formed by a set of permutation matrices P ,s.t. P A = AP . Any subgroup of that group can belinked to an admissible cluster synchronization patternvia orbital partitions. Namely, all the subsets of the net-work nodes that get mapped to themselves (and thusbelong to the same cell of the orbital partition) can becompletely synchronized [17]. The approach can be gen-eralized to systems with different types of nodes and in-teractions, e.g., multilayer networks of coupled oscillatorswhere cluster synchronization requires compatibility be-tween intra- and interlayer symmetries [12]. Similarly tomore general equitable partition methods, symmetries ofdyadic projected networks can not be immediately trans-lated to higher order interactions.Instead, to consider the hypergraph synchronizationpatterns from the symmetry perspective, we formulatethe cluster synchronization condition in terms of symme-tries of the hyperedges of each order m as: P I ( m ) = I ( m ) P ( m )edge . (B1)Here, P N × N is a permutation matrix that reorders thenodes, and [ P ( m )edge ] N m × N m correspond to the permuta-tions of the edge labels if node labels are permuted. Here, N m denotes the number of distinct hyperedges of order m in the hypergraph. These hyperedge permutation ma-1 v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v v FIG. 6. An admissible pattern of synchronization based onhypergraph symmetry. trices are defined as follows:[ P ( m )edge ] e i ,e j = [ P ( m ) ] i j ... [ P ( m ) ] i m j m , (B2)where e i = { i , ..., i m } and e j = { j , ..., j m } are the hy-peredges. Equivalently, the condition can be written interms of the adjacency tensor:[ P ( m ) ] ki [ A ( m ) ] ij ,...,j m − = [ A ( m ) ] kk ,...,k m − [ P ( m ) ] k j ... [ P ( m ) ] k m j m − . (B3)Finally, the largest common subgroup of the symmetrygroups of each edge order determines is the automor-phism group of the hypergraph. The orbits of the sub-groups of the automorphism group determine the admis-sible cluster synchronization patterns.Since the conditions in Eqs. (B1) and (B3) need to holdsimultaneously for all the higher order interactions, oneway to find the automorphism group of the higher ordersystem is to determine the symmetry group respectingthe structure of A (2) (or I (2) ) and check which of theelements of the symmetry group are also the symmetriesof higher order interactions. In this case, only standardcomputational group theory tools are necessary to findthe symmetries of the system.Fig. 6 provides an example system to illustrate howcluster synchronization arises from hypergraph symme-tries. One of the permutations that leaves the structureof the hypergraph in Fig. 6(a) invariant, in cycle nota-tion, is: P =(1 , , , , , , , , , , . (B4)The hyperedge permutation matrices derived from P are recorded in: P (2)edge = ([1 , , [2 , , [3 , , , [2 , , [3 , , , [5 , , [6 , , , [5 , , [6 , , , [5 , , [6 , , , [8 , , [9 , ,P (3)edge =([1 , , , [2 , , , [3 , , , , , [5 , , , [6 , , . (B5)Since P (2)edge and P (3)edge both leave the hypergraph struc-ture invariant, P is a symmetry of the system. The clus-ter synchronization pattern is shown on Fig. 6(a): eachnode is assigned into the same cluster as it is in Fig. 1(b),but the cluster assignment is determined by an orbitalpartition. Symmetry considerations alone would are notsufficient to analyze the state for the coupling topologyin Fig. 2.Finally, we demonstrate that the symmetries of pro-jection network do not necessarily translate to symme-tries of the original hypergraph, making considering thehyperedge permutation matrices a necessary step in find-ing the admissible cluster synchronization patterns. Forinstance, Fig. 3(b) shows the projected network of thehypergraph in Fig. 3(a). The permutation matrix per-muting nodes 1 and 6 ( P = (1 , P applied to the hyperedge[1 , ,
4] transforms it to [6 , , P . Therefore, the symme-try based conditions for cluster synchronization on thehypergraph and its projected network are distinct. Appendix C: Stability calculations for patternsarising from orbital partitions
Due to a general result from equivariant dynamical sys-tems theory, the Jacobian evaluated on the cluster syn-chronization state commutes with the elements (we de-note the actions of these elements by P ) of the symmetrygroup whose orbital partition determines the structureof that cluster synchronization state [36]. Therefore, thefollowing holds: J cs P = P J cs . (C1)where J cs is the full Jacobian of the system ( δ ˙ x = J cs δx )evaluated at the cluster synchronization state. As a re-sult, the Jacobian can be block diagonalized using thematrices that block diagonalize the symmetry group ele-ments P .The terms of Eq. (15) can be used to demonstrate whythe full Jacobian commutes with the symmetry groupaction. Using Eq. (B3) and noting that the permuta-tion matrices do not change the structure of j ,...,j m or2 E C ( m ) k,l , we show that the terms that contribute to thefull Jacobian commute with the action of the symmetrygroup: P ki [ A ( m ) ◦ E C ( m ) k,l ] ij ,...,j α ,...,j m − j ,...,j m − = [ A ( m ) ◦ E C ( m ) k,l ] kk ,...,k α ,...,k m − P k j ...P k m j m − j ,...,j m − = [ A ( m ) ◦ E C ( m ) k,l ] kj ,...,j α ,...,j m − j ,...,j m − P ij α , (C2)thus posing restrictions on the form of the Jacobian.Therefore, group representation theory can be used toblock diagonalize the Jacobian to simplify the stabilitycalculations in a similar way they are used for systems with dyadic interactions [17, 41]. Alternatively, other si-multaneous block diagonalization methods are applicableand can result in a finer block diagonal structure [29].The steps in this process may appear simpler thanthose discussed in Section III. Additionally, they apply toboth directed and undirected hypergraphs. However, thecalculation of irreducible representations that are thenused to find the transformation of J cs is computationallyexpensive. Moreover, the method is only applicable tosystems where the state arises from symmetries, and notthe larger class of systems with patterns of cluster syn-chronization arising from balanced equivalence relations. [1] Mark Newman. Networks . Oxford University Press, 2018.[2] J¨urgen Jost and Raffaella Mulas. Hypergraph Laplaceoperators for chemical reaction networks.
Advances inMathematics , 351:870–896, 2019.[3] Steffen Klamt, Utz-Uwe Haus, and Fabian Theis. Hy-pergraphs and cellular networks.
PLoS Comput Biol ,5(5):e1000385, 2009.[4] Yi Han, Bin Zhou, Jian Pei, and Yan Jia. Understandingimportance of collaborations in co-authorship networks:A supportiveness analysis approach. In
Proceedings ofthe 2009 SIAM international conference on data mining ,pages 1112–1123. SIAM, 2009.[5] Rodica Ioana Lung, No´emi Gask´o, and Mihai AlexandruSuciu. A hypergraph model for representing scientificoutput.
Scientometrics , 117(3):1361–1379, 2018.[6] Leonie Neuh¨auser, Andrew Mellor, and Renaud Lam-biotte. Multibody interactions and nonlinear consen-sus dynamics on networked systems.
Physical Review E ,101(3):032310, 2020.[7] Federico Battiston, Giulia Cencetti, Iacopo Iacopini,Vito Latora, Maxime Lucas, Alice Patania, Jean-GabrielYoung, and Giovanni Petri. Networks beyond pairwiseinteractions: structure and dynamics.
Physics Reports ,2020.[8] Leo Torres, Ann S Blevins, Danielle S Bassett, and TinaEliassi-Rad. The why, how, and when of representationsfor complex systems. arXiv preprint arXiv:2006.02870 ,2020.[9] Young Sul Cho, Takashi Nishikawa, and Adilson E Mot-ter. Stable chimeras and independently synchronizableclusters.
Physical Review Letters , 119(8):084101, 2017.[10] Lucia Valentina Gambuzza, Alessio Cardillo, AlessandroFiasconaro, Luigi Fortuna, Jesus G´omez-Gardenes, andMattia Frasca. Analysis of remote synchronization incomplex networks.
Chaos: An Interdisciplinary Journalof Nonlinear Science , 23(4):043103, 2013.[11] Xiang Wei, Jeffrey Emenheiser, Xiaoqun Wu, Jun-an Lu,and Raissa M D’Souza. Maximizing synchronizability ofduplex networks.
Chaos: An Interdisciplinary Journal ofNonlinear Science , 28(1):013110, 2018.[12] Fabio Della Rossa, Louis Pecora, Karen Blaha, AfrozaShirin, Isaac Klickstein, and Francesco Sorrentino. Sym-metries and cluster synchronization in multilayer net-works.
Nature Communications , 11(1):1–17, 2020. [13] Matteo Lodi, Fabio Della Rossa, Francesco Sorrentino,and Marco Storace. Analyzing synchronized clusters inneuron networks.
Scientific Reports , 10(1):1–14, 2020.[14] Iacopo Iacopini, Giovanni Petri, Alain Barrat, and VitoLatora. Simplicial models of social contagion.
NatureCommunications , 10(1):1–9, 2019.[15] Leonie Neuh¨auser, Michael T Schaub, Andrew Mellor,and Renaud Lambiotte. Opinion dynamics with multi-body interactions. arXiv preprint arXiv:2004.00901 ,2020.[16] Vladimir N. Belykh, Grigory V. Osipov, Valentin S.Petrov, Johan A. K. Suykens, and Joos Vande-walle. Cluster synchronization in oscillatory networks.
Chaos: An Interdisciplinary Journal of Nonlinear Sci-ence , 18(3):037106, 2008.[17] Louis M Pecora, Francesco Sorrentino, Aaron M Hager-strom, Thomas E Murphy, and Rajarshi Roy. Clustersynchronization and isolated desynchronization in com-plex networks with symmetries.
Nature Communications ,5:4079, 2014.[18] Ian Stewart, Martin Golubitsky, and Marcus Pivato.Symmetry groupoids and patterns of synchrony in cou-pled cell networks.
SIAM Journal on Applied DynamicalSystems , 2(4):609–646, 2003.[19] Peter Ashwin and Ana Rodrigues. Hopf normal formwith sn symmetry and reduction to systems of nonlin-early coupled phase oscillators.
Physica D: NonlinearPhenomena , 325:14–24, 2016.[20] Christian Bick, Peter Ashwin, and Ana Rodrigues. Chaosin generically coupled phase oscillator networks with non-pairwise interactions.
Chaos: An Interdisciplinary Jour-nal of Nonlinear Science , 26(9):094814, 2016.[21] Per Sebastian Skardal and Alex Arenas. Abrupt desyn-chronization and extensive multistability in globallycoupled oscillator simplexes.
Physical Review Letters ,122(24):248301, 2019.[22] Christian Bick, Tobias B¨ohle, and Christian Kuehn.Multi-population phase oscillator networks with higher-order interactions. arXiv preprint arXiv:2012.04943 ,2020.[23] Per Sebastian Skardal and Alex Arenas. Memory se-lection and information switching in oscillator networkswith higher-order interactions.
Journal of Physics: Com-plexity , 2020. [24] Ana P Mill´an, Joaqu´ın J Torres, and Ginestra Bianconi.Explosive higher-order kuramoto dynamics on simpli-cial complexes. Physical Review Letters , 124(21):218301,2020.[25] Guilherme Ferraz de Arruda, Michele Tizzani, andYamir Moreno. Phase transitions and stability ofdynamical processes on hypergraphs. arXiv preprintarXiv:2005.10891 , 2020.[26] Raffaella Mulas, Christian Kuehn, and J¨urgen Jost. Cou-pled dynamics on hypergraphs: Master stability of steadystates and synchronization.
Phys. Rev. E , 101:062313,Jun 2020.[27] LV Gambuzza, F Di Patti, L Gallo, S Lepri, M Romance,R Criado, M Frasca, V Latora, and S Boccaletti. Themaster stability function for synchronization in simplicialcomplexes. arXiv preprint arXiv:2004.03913 , 2020.[28] Yuanzhao Zhang, Vito Latora, and Adilson E Motter.Unified treatment of dynamical processes on generalizednetworks: Higher-order, multilayer, and temporal inter-actions. arXiv preprint arXiv:2010.00613 , 2020.[29] Yuanzhao Zhang and Adilson E. Motter. Symmetry-independent stability analysis of synchronization pat-terns.
SIAM Review , 62(4):817–836, 2020.[30] Timoteo Carletti, Duccio Fanelli, and Sara Nicoletti. Dy-namical systems on hypergraphs.
Journal of Physics:Complexity , 1(3):035006, aug 2020.[31] Hiroko Kamei and Peter JA Cock. Computation of bal-anced equivalence relations and their lattice for a coupledcell network.
SIAM Journal on Applied Dynamical Sys-tems , 12(1):352–382, 2013.[32] Michael T Schaub, Neave O’Clery, Yazan N Billeh,Jean-Charles Delvenne, Renaud Lambiotte, and Mauri-cio Barahona. Graph partitions and cluster synchroniza-tion in networks of oscillators.
Chaos: An Interdisci- plinary Journal of Nonlinear Science , 26(9):094821, 2016.[33] Anastasiya Salova and Raissa M. D’Souza. Decoupledsynchronized states in networks of linearly coupled limitcycle oscillators.
Physical Review Research , 2(4):043261,2020.[34] John W Aldis. A polynomial time algorithm to determinemaximal balanced equivalence relations.
InternationalJournal of Bifurcation and Chaos , 18(02):407–427, 2008.[35] Manuela AD Aguiar and Ana Paula S Dias. The lat-tice of synchrony subspaces of a coupled cell network:Characterization and computation algorithm.
Journal ofNonlinear Science , 24(6):949–996, 2014.[36] Martin Golubitsky and Ian Stewart.
The symmetry per-spective: from equilibrium to chaos in phase space andphysical space , volume 200. Springer Science & BusinessMedia, 2003.[37] Francesco Sorrentino, Louis M Pecora, Aaron M Hager-strom, Thomas E Murphy, and Rajarshi Roy. Completecharacterization of the stability of cluster synchroniza-tion in complex dynamical networks.
Science Advances ,2(4):e1501737, 2016.[38] Maxime Lucas, Giulia Cencetti, and Federico Battiston.Multiorder laplacian for synchronization in higher-ordernetworks.
Phys. Rev. Research , 2:033410, Sep 2020.[39] Joseph D. Hart, Don C. Schmadel, Thomas E. Murphy,and Rajarshi Roy. Experiments with arbitrary networksin time-multiplexed delay systems.
Chaos: An Interdis-ciplinary Journal of Nonlinear Science , 27(12):121103,2017.[40] Martin Golubitsky, Ian Stewart, and David G Schaeffer.
Singularities and groups in bifurcation theory , volume 2.Springer Science & Business Media, 2012.[41] Francesco Sorrentino, Louis M Pecora, and Ljiljana Tra-jkovic. Group consensus in multilayer networks.