Graph partitions and cluster synchronization in networks of oscillators
Michael T. Schaub,1, 2, a) Neave O’Clery,3 Yazan N. Billeh,4 Jean-Charles Delvenne,1, 5 Renaud Lambiotte,2 and
Mauricio Barahona6, b)
1)

ICTEAM, Université catholique de Louvain, B-1348 Louvain-la-Neuve, Belgium
naXys and Department of Mathematics, University of Namur, B-5000 Namur,
Belgium
3)
Center for International Development, Harvard University, Cambridge, MA 02138,
United States of America
4)
Computation and Neural Systems Program, California Institute of Technology, Pasadena, CA 91115,
United States of America
5)
CORE, Université catholique de Louvain, B-1348 Louvain-la-Neuve, Belgium
6)
Department of Mathematics, Imperial College London, London SW7 2AZ, United Kingdom

arXiv:1608.04283v2 [physics.soc-ph] 16 Aug 2016

2)

Synchronization over networks depends strongly on the structure of the coupling between the oscillators.
When the coupling presents certain regularities, the dynamics can be coarse-grained into clusters by means
of External Equitable Partitions of the network graph and their associated quotient graphs. We exploit
this graph-theoretical concept to study the phenomenon of cluster synchronization, in which different groups
of nodes converge to distinct behaviors. We derive conditions and properties of networks in which such
clustered behavior emerges, and show that the ensuing dynamics is the result of the localization of the
eigenvectors of the associated graph Laplacians linked to the existence of invariant subspaces. The framework
is applied to both linear and non-linear models, first for the standard case of networks with positive edges,
before being generalized to the case of signed networks with both positive and negative interactions. We
illustrate our results with examples of both signed and unsigned graphs for consensus dynamics and for partial
synchronization of oscillator networks under the master stability function as well as Kuramoto oscillators.
Synchronization of coupled oscillators is ubiquitous in nature: from the rhythmic flashing of fireflies or the orchestrated chirping of crickets to the
entrainment of circadian rhythms or the coherent firing of neurons in epilepsy to the dynamics
of man-made networks, such as power grids and
computer networks. Synchronization is also related to consensus processes, such as the flocking
of birds or shoaling of fish, or opinion formation
in social networks.
Previous studies have typically focused on complete synchronization, where all agents on a network converge to the same dynamics. However,
many networks display patterns of synchronized
clusters, where different groups of agents converge to distinct behaviors. Here we use tools
from graph theory to study the phenomenon of
cluster synchronization. We show that cluster
synchronization can emerge in networks that can
be partitioned into groups according to an external equitable partition (EEP) of the graph. Our
graph-theoretical approach allows us to extend
the analysis to networks with positive and negative links, which are important to describe social interactions and inhibitory-excitatory interactions in biology. We showcase applications to
consensus dynamics and to generic synchronization of oscillators, including the classic Kuramoto
model, and discuss general applications to net-

a) Electronic
b) Electronic

mail: michael.schaub@uclouvain.be
mail: m.barahona@imperial.ac.uk

worked systems of interacting agents.

I.

INTRODUCTION

Synchronization phenomena are prevalent in networked systems in biology, physics, chemistry, as well as
in social and technological networks. The study of these
pervasive processes thus spans many disciplines leading
to a rich literature on this subject1–7 . The synchronization literature has traditionally focused on the problem of
total synchronization, initially under mean field or global
coupling1,8 and more recently studying how total synchronization relates to properties of the interaction topology and the dynamics of the individual agents1–4,9–14 .
Currently, there is a surge of interest in localized synchronization processes, where parts of the network become locally synchronized. This phenomenon may also
be referred to as partial synchronization, cluster synchronization, or polysynchrony 15–25 . Recent work has shown
that the predisposition of a network of coupled oscillators
to exhibit cluster synchronization is intimately linked
to symmetries present in the coupling24,25 . In particular, Pecora and collaborators showed how one can use
the inherent symmetry group of the network to blockdiagonalize the coupling, thereby assessing the stability of cluster synchronization under the master stability
function (MSF) formalism10,11 .
Here, we will also be concerned with the subject of
cluster synchronization of oscillators in networks with
general topologies. However, instead of using a grouptheoretic viewpoint, we will consider this problem from

2
an alternative graph-theoretical perspective. Specifically,
we derive results for cluster synchronization in networks
of oscillators using the notion of external equitable partition (EEP), a concept that has gained prominence in systems theory to study consensus processes26–29 . The use
of EEPs emphasizes the presence of an invariant subspace
in the coupling structure and leads to a coarse-grained
description of the network in terms of a quotient graph.
This approach complements the group-theoretical symmetry viewpoint in Refs.24,25 , while also encompassing
the analysis of networks of Kuramoto oscillators8,30 , a
prototypical model for phase synchronization which does
not lend itself to the MSF formalism.
In addition, we show how the EEP perspective of
cluster synchronization can be generalized to signed
networks, i.e., graphs with links of positive and negative
weights. To do this, we define the notion of signed
external equitable partition (sEEP) and demonstrate its
applicability on structurally balanced signed networks, a
classic model from the theory of social networks31,32 . In
structurally balanced signed networks, linear consensus
dynamics leads to a form of ‘bipolar consensus’33,34 , in
which nodes split into two factions, i.e., nodes inside
the same faction converge to a common value while the
other faction converges to the same value with opposite
sign. In the synchronization setting, we demonstrate
that the presence of sEEPs can induce a bipolar cluster
synchronization, in which each group of oscillators may
be divided into two ‘out-of-phase’ groups, with trajectories of equal magnitude but opposite sign. Below, we
show how these results appear for signed networks under
the MSF framework as well as for Kuramoto oscillators.
Notation: Our notation is standard. The number
of nodes (vertices) in the network is denoted by N ; the
number of edges (links) by E. We denote the adjacency
matrix of the graph by A = AT , where Aij corresponds to
the weight of the coupling between node (oscillator) i and
j. The graph Laplacian matrix is defined as L = D − A,
where D = diag(A1) is the matrix containing the total coupling
strength of each node on the diagonal, i.e.,
P
Dii = j Aij . From this definition, it is straightforward
to see that the vector of ones 1 is an eigenvector of L with
eigenvalue 0. It is well known that the Laplacian may be
decomposed as L = BW B T , where B is the node-to-edge
incidence matrix and W is a diagonal matrix containing
the (positive) weights of the edges. It therefore follows
that the Laplacian is a positive semidefinite matrix.
To simplify notation, but without loss of generality,
our exposition below is presented for unweighted graphs,
i.e., W = I. However, all our results apply to weighted
graphs by using edge weight matrices appropriately.

spectral properties and, consequently, for dynamical processes associated with the graph. EEPs extend the notion
of equitable partition (EP). An EP splits the graph into
non-overlapping cells {Ci } (groups of nodes), such that
the number of connections to cell Cj from any node v ∈ Ci
is only dependent on i, j. Stated differently, the nodes inside each cell of an EP have the same out-degree pattern
with respect to every cell. For EEPs, this requirement is
relaxed so that it needs to hold only for the number of
connections between different cells Ci , Cj (i 6= j).
Algebraically, these definitions can be represented as
follows26,28,35 . A partition of a graph with N nodes into c
cells is encoded by the N ×c indicator matrix H: Hij = 1
if node i is part of cell Cj and Hij = 0 otherwise. Hence
the columns of H are indicator vectors hi of the cells:
H ..= [h1 , . . . , hc ].

Given the Laplacian matrix L of a graph, we can write
the definition of an EEP as follows:
LH = HLπ .

EXTERNAL EQUITABLE PARTITIONS

External equitable partitions are of interest because
the existence of an EEP in a graph has implications for its

(2)

Here Lπ is the c × c Laplacian of the quotient graph
induced by H:
Lπ = (H T H)−1 H T LH = H + LH,

(3)

where the c × N matrix H + is the (left) Moore-Penrose
pseudoinverse of H. Observe that multiplying a vector
x ∈ RN by H T from the left sums up the components
within each cell, and that H T H is a diagonal matrix with
the number of nodes per cell on the diagonal. Hence H +
may be interpreted as a cell averaging operator 27 .
The quotient graph associated with an EEP is a coarsegrained version of the original graph, such that each cell
of the partition becomes a new node and the weights
between these new nodes are the out-degrees between
the cells in the original graph (see Fig. 1a). Although the
Laplacian of the original graph is symmetric, the quotient
Laplacian will be asymmetric in general. Note that, from
the definition of the Laplacian, there is always a trivial
EEP in which the whole graph is grouped into one cell,
i.e., H = 1 and Lπ = 0.
From (2)–(3), the definition of the EEP can be rewritten solely in terms of L:
LH = HH + LH = PH LH,

(4)

where PH ..= HH + is the projection operator onto the
cell subspace, i.e., it defines an orthogonal projection
onto the range of H.
The operator PH commutes with L:
LPH = HLπ H + = PH LPH = PH L,

II.

(1)

(5)

which follows from (2), (4) and the symmetry of L and
PH . Using the commutation (5), it is easy to show that:
H + L = Lπ H + ,

(6)

3
which summarises the relationship between the cell averaging operator H + and the Laplacians of the original
and quotient graphs.
Remark 1 [Equitable partitions and coupling via
adjacency matrices]:
It is instructive to consider
EEPs with respect to the stricter requirement of equitable partitions (EPs). Given the adjacency matrix A of
a graph, an equitable partition encoded by the indicator
matrix HEP must fulfil:
AHEP = HEP Aπ∗ .

(7)

Hence we can define the adjacency matrix of the EP quotient graph Aπ∗ induced by HEP as:
+
Aπ∗ = HEP
A HEP .

(8)

The adjacency matrix Aπ∗ has diagonal entries corresponding to self-loops in the quotient graph of the EP,
reflecting the number of edges between any two nodes
inside each cell. In contrast, the adjacency matrix of
an EEP cannot be uniquely defined; thus, such in-cell
information is not consistently specified. On the other
hand, the quotient graphs of both EPs and EEPs have
consistently defined Laplacian matrices, due to the
well-known invariance of the Laplacian to the addition
of self-loops in a graph, so that the quotient Laplacian is
unaffected by the internal connectivity inside each cell.
This algebraic argument clarifies why EEPs are defined
in terms of the Laplacian (2). It also follows directly
that every EP is necessarily an EEP, while the converse
is not true.
Remark 2 [Network symmetries and (external)
equitable partitions]:
Recently, Pecora et al.24,25
used the symmetry groups of a graph and their associated irreducible representations to identify possible synchronization clusters in networks of oscillators and to
assess their stability. Their group-theoretical analysis is
intimately related to the graph-theoretical perspective presented here. Indeed, the symmetry groups of the graph induce orbit partitions. Every orbit partition is an equitable
partition, yet the converse is not true: there exist EPs not
induced by any symmetry group35,36 . Recall that EEPs
are a relaxation of EPs in the sense that EEPs disregard
the connections inside each cell. Consequently, EEPs are
defined in terms of the Laplacian matrix (2), in contrast to EPs being defined in terms of the adjacency matrix (7). Interestingly, recent work of Sorrentino et al. 25
introduced ‘adjusted orbit partitions’ induced by symmetry groups of a “dynamically equivalent coupling matrix”
in which internal connections inside each cluster are ignored. Such ‘adjusted orbit partitions’ are in fact EEPs
but, as for EPs, there exist EEPs that cannot be generated by the symmetry groups of such dynamically equivalent coupling matrices. In this sense, EEPs provide a
generalised setting that includes the group-theoretical orbit partitions as a particular case.

As EEPs are a larger class of partitions than EPs and
Laplacians are of wide interest in applications, we concentrate here on networks with Laplacian coupling. All
our results can be applied straightforwardly to systems in
which the coupling is described by the adjacency matrix,
by considering EPs rather than EEPs.

III. CLUSTER SYNCHRONIZATION UNDER THE
EXTERNAL EQUITABLE PARTITION

We first use EEPs to study cluster synchronization on
standard networks, i.e., defined by connected undirected
graphs with positive weights. We start by considering results for linear consensus, and then apply the framework
to nonlinear cluster synchronization both under the MSF
formalism as well as Kuramoto networks.

A.

Dynamical implications of EEPs: the linear case

The definition of the EEP (2) can be understood as
a ‘quasi-commutation’ relation, which signals a certain
invariance of the partition encoded by H with respect to
the Laplacian L. Similarly, Eq. (6) shows that the cell
averaging operator H + exhibits a (distinct) invariance
with respect to L. In particular, Eq. (2) implies that
the associated cell indicator matrix H spans an invariant
subspace of L, whence it follows that there exists a set
of eigenvectors which is localised on the cells of the partition. Furthermore, the eigenvalues associated with the
eigenvectors spanning the invariant subspace are shared
with Lπ , the Laplacian of the quotient graph27 . If L has
degenerate eigenvalues, an eigenbasis can still be chosen
so that it is localised on the cells of the partition27 .
The properties of the EEP (2)–(6) have noteworthy
consequences for linear dynamics dictated by L, as illustrated by the case of linear consensus dynamics27 :
ẋ = −Lx,

(9)

where the N ×1 vector x describes the state of the system.
First, as shown in Fig. 1b, the EEP is consistent with
a form of invariance akin to ‘cluster consensus’. In particular, if the initial state vector is given by x = Hy for
some arbitrary y (i.e., all the nodes within cell Ci have
the same value yi ), the nodes inside the cells remain identical for all times and their dynamics is governed by the
quotient graph:
ẋ = H ẏ

where

ẏ = −Lπ y.

(10)

This follows directly from LHy = HLπ y.
Second, the dynamics of the cell-averaged states hxiCi
is governed by the quotient graph:
dhxiCi
= −Lπ hxiCi
dt

where

hxiCi ..= H + x,

(11)

4

(a)

1

2
1

2

1
2

(b)

nodes inside each cell remain identical for all times27 .
This simple insight, which follows from the impulseresponse of the linear system (12), will be useful when
analysing nonlinear synchronization protocols.
Finally, it is important to remark that since there is
always a trivial EEP spanning the complete graph (with
H = 1), all the results above and henceforth can be
trivially applied to the case of global consensus (global
synchronization) as a particular case.

Quotient Graph

Full Graph

4

1

10
8
6

B. EEPs and nonlinear cluster synchronization within the
MSF framework

4
2
0

(c)

0

0.5

1

time,

1.5

2

10
8
6
4

j

2
0

We now extend the notions introduced above to a more
general setting describing the dynamics of interconnected
nonlinear systems. This framework is known as the Master Stability Function (MSF) and has been pioneered by
Pecora and co-workers10–12 .
We consider networks of identical coupled oscillatory
nonlinear systems in which the dynamics of each node i
is described by:
X
ẋi = F(xi ) − γ
Lij G(xj ),
(13)

0

0.5

1

time,

1.5

2

Figure 1. External equitable partitions and invariant
consensus dynamics. (a) A graph with N = 8 nodes with
an external equitable partition into four cells (indicated with
colors) and its associated quotient graph. (b) The evolution of the consensus dynamics on the full graph (9) from an
initial condition x = Hy is shown with solid lines. The associated quotient dynamics (10) governing y is shown with
circles. Once all states within each cell are equal (i.e., they
are cluster-synchronized), the dynamics will remain clustersynchronized and its dynamics will be described by the quotient dynamics for all times. (c) For consensus dynamics, the
quotient graph dynamics also describes the cell-averaged dynamics (crosses) of the unsynchronized full graph dynamics
(solid lines), as given by (11).

which follows from H + Lx = Lπ H + x. Thus, the cellaveraged dynamics is governed by a lower dimensional
linear model, with dimensionality equal to the number of
cells in the EEP (see Fig. 1c).
Third, the results obtained for the autonomous dynamics with no inputs (9) can be equivalently rephrased for
the system with a bounded input u(t):
ẋ = −Lx + u(t).

(12)

In particular, similarly to (10), we also have cell invariance under inputs: if we apply an input consistent with
the cells of an EEP (i.e., u(t) = Hv(t), v(t) ∈ Rc ), the

where γ is a parameter that regulates the coupling
strength; xi ∈ Rd is the state vector of node i;
F : Rd → Rd is the intrinsic dynamics of each node; and
the coupling function G : Rd → Rd specifies how the
nodes in the network interact according to the interconnection topology described by the graph Laplacian L.
Although, as discussed above, we could consider a coupling mediated by the adjacency matrix (and associated
EPs), we concentrate here on the case of Laplacian coupling (and associated EEPs) as the more generic case of
interest in the literature .
To facilitate the subsequent discussion, we define
x ..= [xT1 , ..., xTN ]T ∈ RN d and use the Kronecker product
to rewrite (13) compactly as a N d-dimensional system of
ODEs:
ẋ = FN (x) − γ(L ⊗ Id )GN (x)

(14)

where FN (x) ..= [F(x1 )T , . . . , F(xN )T ]T ∈ RN d ,
GN (x) ..= [G(x1 )T , . . . , G(xN )T ]T ∈ RN d and Id is the
d-dimensional identity matrix.
A cluster-synchronized state consistent with an EEP
with indicator matrix H is then given by:
xs (t) = (H ⊗ Id ) ys (t)
T

where
T T

(15)
cd

ys (t) = [s1 (t) , . . . , sc (t) ] ∈ R .

1.

(16)

EEPs and invariance of cluster-synchronized states

Let a graph with Laplacian L exhibit a non trivial
EEP with c cells encoded by the indicator matrix H and

5

(a)

Full Graph

(b)

where Fc (y), Gc (y) ∈ Rc d are defined analogously to
FN , GN above, and we have the relations:
(H ⊗ Id )Fc (y) = FN ((H ⊗ Id )y)
(H ⊗ Id )Gc (y) = GN ((H ⊗ Id )y) .

Quotient Graph

In close parallel to the linear case (10), we can derive
the following result for cluster-synchronized dynamics.
Let us have an initial condition that is identical within
the cells of the EEP, i.e., x = (H ⊗Id )y for some arbitrary
y ∈ Rcd at t = 0. Then the nodes within cells of the EEP
remain identical for all time t ≥ 0, and their dynamics
can be described by the dynamics of the quotient graph:

7
1

(c)

(d)

20

20

10

10

0

0

-10

-10

-20
10

0

5

10

15

where

-20
10

0

0

0

-10

-10

-20

0

5

10

15

2

-20

5

10

15

20

ẋ = (H ⊗ Id )ẏ
ẏ = Fc (y) − γ(Lπ ⊗ Id )Gc (y).

(20)

25

This result follows from:

0

5

10

15

20

25

30

ẋ = (H ⊗ Id )ẏ = (H ⊗ Id ) [Fc (y) − γ(Lπ ⊗ Id )Gc (y)]
= FN ((H ⊗ Id )y) − γ(HLπ ⊗ Id )Gc (y)
= FN ((H ⊗ Id )y) − γ(LH ⊗ Id )Gc (y)
= FN ((H ⊗ Id )y) − γ(L ⊗ Id )(H ⊗ Id )Gc (y)
= FN ((H ⊗ Id )y) − γ(L ⊗ Id )GN ((H ⊗ Id )y) 

20

0

Here we have made use of the standard identity
(A ⊗ B)(C ⊗ D) = (AC) ⊗ (BD).

10
-2
-4

(18)
(19)

0
0

5

10

time, t

15

-10

0

5

10

15

20

25

time, t

Figure 2. Synchronization of identical coupled Rössler
oscillators. (a) Chaotic Rössler oscillators are coupled according to a star graph (N = 8). This graph has an
EEP with two cells (indicated with colors), shown with its
quotient graph. (b) Under certain conditions, the coupled oscillators (each with a three-variable dynamics x(t) =
(x1 (t), x2 (t), x3 (t))) can exhibit cluster synchronization commensurate with the EEP: spoke nodes (blue), centre node
(red). (c)-(d) Analogously to linear consensus (Fig. 1), given
an initial condition consistent with the EEP, the dynamics of
the nodes within each cell remain identical. The solid lines are
the full dynamics x(t) governed by (14) with initial condition
x0 = (H ⊗ Id )y0 ; the circles are the quotient graph dynamics y governed by (17). Figures (c) and (d) correspond to
two values of the coupling parameter γ, and in both cases
the dynamics remains cluster-synchronized within the EEP.
In (c), with γ = 0.3, the total synchronization of the quotient
graph is stable, and both x and y converge to the completely
synchronized solution. In (d), with γ = 0.03, total synchronization of the quotient graph dynamics is not linearly stable;
hence the system exhibits sustained cluster synchronization.

quotient Laplacian Lπ . The dynamics of the cell variables
y ..= [y1T , ..., ycT ]T ∈ Rc d associated with the quotient
graph is then given by
ẏ = Fc (y) − γ(Lπ ⊗ Id )Gc (y),

(17)

Example [Coupled Rössler oscillators] : Consider
a network of N = 8 oscillators where each node has a
three dimensional dynamics (d = 3) given by the chaotic
Rössler system37


−x2 − x3
F(x) = F([x1 , x2 , x3 ]T ) =  x1 + ax3 
(21)
b + x3 (x1 − c)
with parameters a = b = 0.2 and c = 7. The oscillators are coupled through the variable x1 according to the
linear function:


1 0 0
G(x) = 0 0 0 x.
0 0 0
The topology of interconnection is a star graph, which
has an EEP with two cells (c = 2): one cell comprises the central node, the other cell contains all other
nodes (Fig. 2).
Let the initial condition be x0 = (H ⊗ Id )y0 . Then the
variables of the nodes within each cell remain identical
at all times, i.e., the dynamics stays cluster-synchronized
(Fig. 2). For γ = 0.3, this poly-synchronous state (while
remaining cluster-synchronized at all times) evolves towards the globally synchronized state (Fig. 2c). In contrast, for γ = 0.03, the cluster synchronization of the
cells does not converge towards global synchrony, since
the completely synchronized state of the quotient graph
dynamics is no longer (linearly) stable (Fig. 2d and
Sec. III B 3).

6
2.

(b)

(a)

EEPs and cell-averaged synchronization dynamics

10
4

Although the invariance of cluster-synchronized EEP
states carries over to the nonlinear MSF setting, the second finding of the linear analysis, namely that the dynamics of cell averages is described by the quotient graph
dynamics (17), does not hold in general. Indeed, after
some algebraic manipulations it is easy to see that:

10

(H + ⊗ Id )ẋ = (H + ⊗ Id )[FN (x) − γ(L ⊗ Id )GN (x)]

5

= (H + ⊗ Id )FN (x) − γ(Lπ ⊗ Id )(H + ⊗ Id )GN (x).

0

Due to their nonlinearity, in general F and G do not
commute with the linear cell-averaging operation:
+

+

(H ⊗ Id )FN (x) 6= Fc ((H ⊗ Id )x)

5

2
0

0
-2
-5

0

5

10

15

3. Stability of EEP cluster-synchronized states through
the MSF formalism

The sections above lead naturally to consider the stability of EEP cluster synchronization. Following Pecora
et al24,25 , the linearized stability around any cluster synchronized state can be evaluated using the MSF framework via the variational expression:
E

20

25

0

5

10

15

20

25

0

5

10

15

20

25

2
0

-5
-10

-2
0

5

10

15

2

where DFN (x) and DGN (x) denote the Jacobians of FN
and GN for state x.
This result implies that if the cluster-synchronized
state is stable, the averaging operator will approximately
commute with both FN and GN when the state is close
to the cluster-synchronized state. As a consequence, an
appropriately chosen initial condition of the average cell
dynamics will remain close (or converge) to the quotient
dynamics, as seen in Fig. 3a. On the other hand, the
interplay of the Jacobians of G and F and the graph
structure encoded by L and Lπ can render the initial
perturbation unstable, and the state will exponentially
diverge. In that case, the quotient dynamics will not be
a good model for the cell-averaged dynamics, as seen in
Fig. 3b. We explore these points through the MSF formalism in Section III B 3, where we consider the stability
of the cluster-synchronized state (including the globally
synchronized state).

δ ẋ(t) =

15

4

FN (xs + ) ≈ FN ((H ⊗ Id )ys ) + DFN (xs ) ,
GN (xs + ) ≈ GN ((H ⊗ Id )ys ) + DGN (xs ) ,

(i)

10

1

Hence, unlike the linear case, the cell-averaged dynamics
is not strictly equivalent to the synchronization dynamics
governed by the Laplacian of the quotient graph.
However, an approximate equivalence is obtained if we
consider an  perturbation around a cluster-synchronized
state (15). To first order we then have,

c 
X

5

4

(H + ⊗ Id )GN (x) 6= Gc ((H + ⊗ Id )x).

"

0

#


(i)
⊗ DF(si ) − γ LE ⊗ DG(si ) δx(t),


i=1

(22)

0

0
-2

-1

0

5

10

15

time, t

time, t

Figure 3. Cell-averaged dynamics and synchronization
of identical Rössler coupled oscillators. The numerics
in this Figure follow Figure 2, but here we focus on the cellaveraged dynamics of the same system of coupled Rössler oscillators and its relationship with the dynamics of the quotient graph. (a) For γ = 0.3 and an initial condition close
to the synchronization manifold, the quotient dynamics (circles) evolves closely to the cell-averaged full system dynamics
(crosses), and both dynamics converge to the totally synchronized solution. Note the individual time courses of each of
the eight oscillators converging also to this solution. (b) For
γ = 0.03, however, the cell-averaged (crosses) and quotient
dynamics (circles) diverge, as is clearly visible for large times.

where si ∈ Rd is the (consistent) state of every node in
the ith cluster, as defined in (15), and E (i) are identity
matrices consigned to each cluster:
E (i) ..= diag(hi ),

(23)

as given by the cell indicator vectors (1). Using computational group theory, Pecora et al. block-diagonalize
the above expression to assess the stability of any clustersynchronized state.
As an alternative to symmetry-based arguments, the
MSF variational analysis may also be understood using EEPs and their associated indicator matrices. Here
we use the fact that eigenvectors and eigenvalues are
shared between the Laplacians of the original and quotient graphs27 . Let us denote the c eigenvectors of the
quotient Laplacian Lπ by V π = [v1π , . . . , vcπ ] with eigenvalues Λπ = diag(λπi ) such that Lπ V π = V π Λπ . The
properties of the EEP27 ensure that a subset of the eigenvectors of the full Laplacian L are directly related to the

7
eigenvectors of Lπ :
Vs = HV π ∈ RN ×c .

(24)

These are the eigenvectors that define the cluster synchronization manifold commensurate with the EEP.
The eigenvectors orthogonal (transversal) to the clustersynchronized manifold are denoted by V⊥ ∈ RN ×(N −c) .
These are the eigenmodes that drive the system out of a
cluster-synchronized state, and therefore we want these
modes to be damped. An orthogonal matrix of eigenvectors of L that diagonalizes the Laplacian:
V T LV = Λ

Λ ..= diag(λi ).

where

(25)

is thus given by
V = [Vs , V⊥ ] = [HV π , V⊥ ].

(26)

Hence the first c columns correspond to eigenvectors of
L (with eigenvalues λi = λπi , i = 1, . . . , c) that can be
mapped to Lπ , and the second block of (N − c) columns
corresponds to the transversal manifold.
Using V to diagonalize L via the coordinate transformation δχ = (V T ⊗ In )δx leads to:
δ χ̇(t) = (V T ⊗ In )

c
hX

E (i) ⊗ DF(si )

i=1

−γ

c
X

i
LE (i) ⊗ DG(si ) (V ⊗ In )δχ(t)

(27)

i=1

=

c
hX

V T E (i) V ⊗ DF(si )

i=1

−γ

c
X

i
V T LE (i) V ⊗ DG(si ) δχ(t)

(28)

i=1

"

#
c 



X
(i)
(i)
=
Q ⊗ DF(si ) − γ ΛQ ⊗ DG(si ) δχ(t),
i=1

(29)

where we have


V T LE (i) V = Λ V T E (i) V



=.. Λ Q(i) .

(30)

The structure of the matrices Q(i) means that the
modes in the cluster synchronization manifold are effectively decoupled from the modes transversal to it. To see
this, note that from V⊥T Vs = 0 and (24) it follows that the
transversal eigenvectors V⊥ lie in the orthogonal subspace
to H: H T V⊥ = 0. (This also means that every transversal mode is mean-free within each cell: H + V⊥ = 0.)
Therefore, we have the following effective decoupling between the cluster-synchronized and transversal modes:
V⊥T E (i) Vs = V⊥T E (i) HV π = V⊥T [0, . . . , hi , 0, . . .]V π = 0,
leading to
"
Q(i) = V T E (i) V =

(i)
Qs

0(N −c)×c

#

0c×(N −c)
.
(i)
Q⊥

By examining this matrix, we can obtain information
about the (local) stability of the cluster-synchronized
state (see Ref.25 for a related discussion). In order
to check the linear stability of the cluster-synchronized
manifold, it is enough to check that all the transversal
modes are damped. Yet such damping of the transversal modes alone does not specify the behavior within
the cluster-synchronized manifold, or indeed the convergence towards any of the different cluster-synchronized
states within it. Further damping within the clustersynchronized manifold would lead the dynamics to converge to an even lower-dimensional manifold, i.e., towards a particular subset of the cluster-synchronized
states. Stated differently, some of the cells in a clustersynchronized state could merge, leading to another state
with fewer cells. If damping within the manifold is
present, it can lead to convergence towards the completely synchronized state, akin to the numerics in Fig. 2c
(and in contrast to the numerics in Fig. 2d where such
convergence within the manifold is not observed).
C.

EEP cluster synchronization in Kuramoto networks

The MSF framework provides a powerful tool for the
analysis of nonlinear systems with diffusive couplings,
yet there are important classes of systems that do not
lend themselves naturally to this formulation. Examples
include systems with sinusoidal coupling between oscillators, as in models of power systems38 or the classic
Kuramoto model of coupled oscillators8,30 . The use of
EEPs can nevertheless afford us insight into cluster synchronization in these cases, too.
We consider the Kuramoto model with N oscillators
N
X
dθi
= ωi + λ
Aij sin(θj − θi ),
dτ
j=1

(31)

where θi and ωi describe the phase and intrinsic frequency of each oscillator, respectively, λ is the coupling
parameter, and Aij is the adjacency matrix encoding the
network connectivity.
To simplify our notation below, let us renormalize time
t = λτ . The dynamics of a network of Kuramoto oscillators coupled through a graph with Laplacian L = BB T ,
where B is the incidence matrix of the graph, can then
be rewritten in vector-matrix notation as14 :
1
ω − B sin(B T θ)
λ
1
= ω − BW(B T θ)B T θ
λ
1
= ω − LW(B T θ) θ,
λ

θ̇ =

(32)

(33)

where θ and ω are N -dimensional vectors, and we have
defined W(x) ..= diag(sinc(x)) = diag(sin(xi )/xi ).
This rewriting emphasizes the close relation of the Kuramoto model to Laplacian dynamics. Not only does the

8

(a)

(c)

(b)

Full Graph

Quotient Graph
3

0

-0.5

-0.5

-1

-1

-1.5

-1.5

5

0

-5

-2
0

1

(d)

0

5

time,

10

15

-10

-2
0

5

time,

10

15

0

5

time,

10

15

Figure 4. EEP clustered dynamics on a network of Kuramoto oscillators with identical intrinsic frequencies. (a)
Kuramoto oscillators (N = 8) coupled through the graph shown, which has an EEP indicated by the color of the nodes. This
partition is not an orbit partition, i.e., it is not directly induced by any symmetry group. The associated quotient graph is
shown below. (b) If the initial condition is constant within each cell, the dynamics of the full system (line) is exactly equivalent
to the lower dimensional Kuramoto dynamics on the quotient graph (circles). (c) Provided the initial condition for θ is close
to the cell-averaged state and within an open semi-circle, the linear cell- averaged dynamics (crosses) is closely aligned with
the quotient graph dynamics shown in (b). (d) If the initial condition is spread further on the circle, the linear cell averaging
is no longer aligned with the quotient graph.

linearization for small phase differences lead to the standard linear Laplacian dynamics, but the final equality
underscores the fact that the full Kuramoto model may
still be understood in terms of a weighted Laplacian dynamics with time-varying edge weights14 :
LW(B T θ) ..= BW(B T θ)B T .

(34)

It is therefore not surprising that EEPs give useful insights into invariant dynamics of Kuramoto networks.

1.

Case I: equal intrinsic frequencies

Let us consider first the case where all intrinsic frequencies are identical: ωi = ω, ∀i. In this case, we may
assume ω = 0 without loss of generality, as this is equivalent to grounding the system or defining the phases with
reference to a rotating frame14 . The resulting system
θ̇ = −B sin(B T θ) = −LW(B T θ) θ

(35)

invariant dynamics of cluster-synchronized states in the
full model (35). In other words, we need to show that
H ψ̇ = −B sin(B T Hψ),

i.e., Hψ is invariant under the full dynamics.
To establish this, we use the following fact:
A given EEP for a network remains an EEP if all
edge weights between two distinct cells are multiplied by
a factor that depends only on the two cells.
This fact is a direct consequence of the definition of an
EEP, since such scaling changes the out-degree patterns
of all nodes within a cell consistently.
Remark 3 [EEPs and structured weights]: A particular case of such a rescaling that will be useful below
can be represented algebraically as follows. Consider a
graph with Laplacian L = BB T and an EEP with indicator matrix H. Let the edge weights be scaled consistently
across cells (in the above sense) leading to the modified
Laplacian:
Lw(B T Hξ) = B diag(w(B T Hξ)) B T ,

4,14

is well known to converge
to the totally synchronized
state with identical phases.
Let the graph with Laplacian L = BB T be endowed
with an EEP with partition matrix H. We can then
define the following Kuramoto dynamics taking place on
the quotient graph of the EEP:
ψ̇ = −H + B sin(B T Hψ) = −H + LW(B T Hψ) Hψ (36)

(38)

(39)

c

where ξ ∈ R is a cell vector, and w(x) = w(−x) is a
symmetric function applied element-wise. Since the EEP
remains unchanged under this rescaling, it follows from
(5) that the projection operator associated with H also
commutes with the modified Laplacian:
PH L = LPH =⇒ PH Lw(B T Hξ) = Lw(B T Hξ) PH . (40)

where ψ is the c-dimensional vector containing the phases
associated with the quotient graph, and we use the definition (3) to factorize the quotient Laplacian appropriately:

We now use (40) to show that EEP clustersynchronized states are invariant under Kuramoto dynamics. To see this, left multiply (36) with H:

Lπ = H + LH = (H + B)(B T H).

H ψ̇ = −HH + B sin(B T Hψ)
= −PH LW(B T Hψ) Hψ = −LW(B T Hψ) PH Hψ

(37)

As for the linear case (10), we wish to show that the
cell dynamics on the quotient graph (36) describes an

= −LW(B T Hψ) Hψ = −B sin(B T Hψ),

9
where PH H = H follows from the definition of the projection operator. Note that W(x) = sinc(x) is symmetric. 
The proof shows that the full Kuramoto model follows
the quotient dynamics (36) for all times, if it ever synchronizes to a particular EEP. We illustrate this behavior
in Fig. 4b, where we use a network topology (Fig. 4a) inspired by a construction outlined by Chan and Godsil35,36
highlighting the difference between orbit partitions (generated from symmetry groups) and equitable partitions.
As shown in Fig. 4c, the cell averages are also well described by the quotient dynamics provided the initial condition is not too far away from the EEP-averaged state.
If the phases of the initial condition are outside the open
semicircle (as in Fig. 4d), a naive linear averaging does
not fully capture the convergence on the torus.

2. Case II: non-equal intrinsic frequencies commensurate
with an EEP

Full Graph

(a)

Quotient Graph

3
1

(b) 1.5

1
0.5

1

0
-0.5

0.5

-1
-1.5
0

1

2

3

4

5

6

0
-0.5
-1
0

5

0

5

(c) 1.5

time,

10

15

10

15

1

The analysis for the Kuramoto model with equal frequencies does not apply in general to a network of oscillators with non-equal intrinsic frequencies. However, similar results hold when the oscillators within each cell have
the same frequency. In particular, consider the Kuramoto
system (32) with EEP-commensurate frequencies:
θ̇ =

1
H$ − B sin(B T θ),
λ

(41)

where $ is a c-dimensional vector containing the frequencies of the cells.
In the case of heterogeneous frequencies, the model can
not reach globally identical synchronization, so the ‘most
synchronous’ behavior is the cluster-synchronized state
with identical phases within each cell. By the arguments
in Section III C 1, mutatis mutandis, it is easy to see that
the cluster-synchronized state Hψ is invariant under (41)
and governed by the quotient graph:
ψ̇ =

1
$ − H + B sin(B T Hψ).
λ

(42)

A numerical illustration of this invariance is given in Figure 5. We note that there is a close analogy here to the
scenario of the linear consensus system with an input
commensurate with the EEP (12). Indeed, the intrinsic
frequencies of the cells $ can be interpreted as constant
inputs to each of the cells.
We remark that our results for Kuramoto systems
here are concerned with the invariance of solutions
and not their stability.
As studied previously4,14 ,
the stability of the synchronous state depends on the
magnitude of the spread of the frequencies ωi along the
edges of the graph relative to the coupling parameter
λ. Hence as the coupling λ becomes smaller, and the
norm of $/λ becomes larger, the synchronized (and

0.5
0
-0.5
-1

time,

Figure 5. Cluster synchronization in Kuramoto networks with EEP-commensurate intrinsic frequencies.
(a) Kuramoto dynamics (41) over the same network as in
Figure 4a, but this time with non-identical intrinsic frequencies, yet aligned with the EEP. (b) If the initial condition is
constant within cells, the dynamics of the full system (lines)
is identical to the dynamics of the quotient graph (circles),
and the system eventually settles to a cluster-synchronized
state. Inset: the same dynamics without subtracting the
(time-dependent) mean phase. (c) If the initial phases are
within an open semicircle and close to the cell-averages, then
the quotient dynamics is a good descriptor for the dynamics
for all times.

cluster-synchronized) solutions become unstable.
Remark 4 [Kuramoto model with a phase offset]:
To gain insight into the effect of a phase offset, let us consider the Kuramoto model with equal intrinsic frequencies
and a constant phase offset discussed in Ref.39 , which can
be rewritten as:
θ̇ = −B sin(B T θ + α1)
T

(43)
T

= sin(α)B cos(B θ) − cos(α)B sin(B θ).

(44)

For α → 0, the first term vanishes, and we recover

10
the standard Kuramoto model (35) for which the EEP
analysis holds. Therefore, if θ is within (close proximity to) the cluster synchronization manifold Hψ, by our
arguments above, the system will remain in a polysynchronous, clustered state. As α increases, the magnitude of the Kuramoto coupling parameter (cos α) decreases, whereas at the same time the magnitude of the
spread of the input intrinsic frequencies (sin α) increases.
Hence, as α is increased above a threshold, we expect the
cluster synchronization manifold to lose stability, as for
the case of non-equal frequencies above. This is in line
with Nicosia et al. 39 , who observed numerically that cluster synchrony is lost above a critical value of α.

matrix:
(
[Bσ ]ie =

Many mathematical models for real-world networks
need to incorporate positive and negative interactions.
For instance, in social networks, relationships can be
friendly or hostile, or they reflect trust or distrust between individuals. Therefore, the sign of a link is a
central concept in social psychology, associated with the
emergence of conflict and tension in social systems31,32 ,
and it has gained popularity recently in the study of online social networks40 and online cooperation41 . In biological systems, the sign of an edge is also a key element,
in particular when modelling dynamical processes. For
instance, genes can either promote or repress the expression of other genes in genetic regulatory networks42 , and
neurons can excite or inhibit the firing of other neurons
in neuronal networks and thereby shape the global dynamics of the system43,44 .
A. Signed networks and structural balance: the signed
external equitable partition
1.

The signed Laplacian matrix

For a network with positive and negative interactions,
we can define the signed Laplacian matrix of the network
as follows33,45 :
Lσ = Dabs − A,

(45)

where Dabs = diag(|A|1) is the diagonal absolute degree
matrix, and A is again the adjacency matrix (which may
hereafter contain both positive and negative weights). As
for the standard Laplacian, it can be shown that the
signed Laplacian is positive semidefinite and its spectrum
contains one zero eigenvalue when the graph is connected
and structurally balanced. To see this, note that the
signed Laplacian can be expressed as:
Lσ = Bσ Wabs BσT ,

(46)

where Wabs = diag(|we |) is the absolute edge weight matrix and Bσ ∈ RN ×E is the signed node-to-edge incidence

if i is the tail of edge e,
if i is the head of edge e,

Henceforth, we assume Wabs = I without loss of generality. By using the signed Laplacian, the construction of an
EEP can be extended to signed graphs. To do so, however, we must first introduce the notion of structurally
balanced graph, which will enable us to define the notion
of a signed external equitable partition (sEEP).

2.
IV. CLUSTER SYNCHRONIZATION IN NETWORKS
WITH POSITIVE AND NEGATIVE WEIGHTS

1
−sign(e)

Structurally balanced graphs

Following Cartwright and Harary32 , a signed graph is
defined to be structurally balanced if the product of the
signs along any closed path in the network is positive.
This definition implies that only ‘consistent’ social relationships are allowed in triangles of three nodes: either all
interactions are positive, or there are exactly 2 negative
links, which may be interpreted in the sense that “the enemy of my enemy is my friend”31 . Equivalently, a signed
network is structurally balanced if it can be split into
two factions, where each faction contains only positive
interactions internally, while the connections between the
two factions are purely antagonistic (see Fig. 6). It has
been shown that many social networks are close to being
structurally balanced46 , suggesting that there might be a
dynamical process acting on such systems driving them
towards structural balance47,48 .
The following characterization of a structurally balanced graph based on the signed Laplacian was highlighted by Altafini33,34 . A network is structurally balanced if there exists a diagonal matrix Σ = diag(σ), with
±1 on the diagonal, such that the matrix:
L0 = ΣLσ Σ

(47)

contains only negative elements on the off-diagonal. In
other words, the signed Laplacian can be transformed
into the standard Laplacian of an associated graph with
only positive weights through the similarity transformation defined by Σ. The matrix Σ is called switching equivalence, signature similarity, or gauge transformation in
the literature34 . Using this characterization, one can efficiently determine whether a network is structurally balanced46 and obtain the corresponding switching equivalence matrix Σ. Note that it follows trivially that a
standard network with only positive weights is always
structurally balanced with Σ = IN . Hence Lσ is a generalization of the standard Laplacian.

3.

Signed external equitable partitions

Using the signed Laplacian, we extend the concept of
EEP to structurally balanced signed networks. Consider

11

(a)

B.

(b)

The definition of a sEEP provides us with an appropriate tool for the analysis of cluster synchronization in
structurally balanced signed networks, as we now show.
The results in this section parallel those obtained for unsigned graphs, hence we concentrate on the distinctive
features of clustered dynamics in signed networks.

Dynamics and signed external equitable partitions

2
1. The linear case: bipolar cluster synchronization in
signed consensus dynamics

0
A remarkable feature of structurally balanced networks
is that the linear signed consensus dynamics33 :

-2

ẋ = −Lσ x,

0

5

time,

10

15

Figure 6. Structurally balanced graphs and bipolar
consensus. (a) Example of a structurally balanced signed
graph (red links positive, blue links negative). Every cycle has
an even number of negative links or, equivalently, the graph
can divided into two factions given by their polarization σi
(corresponding to the green and magenta groups). Note that
each of the factions has only positive weights inside and only
negative weights between them. (b) The signed consensus
dynamics (51) on a structurally balanced graph always leads
to a bipolar consensus, in which each node agrees with the
nodes within its own faction, but has exactly the opposite
sign to any node in the other faction.

a structurally balanced signed graph with signed Laplacian Lσ and denote the Laplacian of the positive switching equivalent graph by L0 = ΣLσ Σ. Let H denote the
indicator matrix of an EEP of L0 :
L0 H = HLπσ .

(48)

Then there exists a signed indicator matrix Hσ = ΣH
that defines an invariant subspace for Lσ :
Lσ Hσ = Hσ Lπσ ,

(49)

which follows from the definition (47) and Σ2 = IN .
We define the partition πσ with indicator matrix Hσ =
ΣH as the signed external equitable partition (sEEP), and
its associated quotient graph is given by:
Lπσ = Hσ+ Lσ Hσ = H + L0 H.

(50)

Therefore cells in a sEEP contain nodes with the same
out-degree pattern in absolute value. An illustration of a
sEEP and associated quotient graph is shown in Fig. 7.
Note that the quotient graph only has positive weights.

(51)

converges to a polarized state, in which the nodes are
divided into two sets with final values that are equal
in magnitude but opposite in sign (Fig. 6). Stated
differently, the eigenvector of Lσ associated with the
zero eigenvalue has the form σ = [σ1 , . . . , σN ]T , where
σi ∈ {−1, +1}, ∀i. As shown in Ref.33 , this implies that
the system dynamics (51) converges to the final state
lim x(t) =

t→∞

σ T x0
σ,
N

(52)

and the sign pattern of the eigenvector σ corresponds precisely to the switching equivalence transformation, i.e.,
Σ = diag(σ). In the following, we will refer to Σii = σi
as the polarization of node i. Note that the the vector
σ is only defined up to an arbitrary sign, so only the
relative polarization of the nodes is relevant.
We now extend the analysis to networks endowed with
a sEEP. The presence of a sEEP (49) has similar dynamical implications to the presence of an EEP in the case of
a positive graph. The following statements can be proved
analogously to the standard (unsigned) consensus case.
First, sEEP cluster-synchronized states Hσ y(t) are invariant under the full linear dynamics (51). Hence, given
an initial condition x = Hσ y consistent with a sEEP, x(t)
remains in the sEEP state x(t) = Hσ y(t) for all times,
and the dynamics is governed by the quotient graph:
ẏ = −Lπσ y (see Figure 7b).
In contrast to standard unsigned graphs, the variable
of every node within a cell of the cluster-synchronized
state will have the same magnitude, but its sign may be
inverted depending on its polarization σi . Therefore, in
signed networks each cell maybe itself divided into two
factions whose values are of equal magnitude (as given
by the quotient dynamics), but of opposite sign, as illustrated in Figure 7 (see how node 8 has the opposite sign
to nodes 6, 7 all in the cyan cell). We use the term bipolar
cluster synchronization to account for this phenomenon
in signed networks.
Second, the signed cell-averaged dynamics hxiσCi =
+
Hσ x is determined by the dynamics of the quotient graph

12

(a)

3

(c)

(b)

Full Graph

0

2

8

5

1

-5
5

7

0

1

2

3

4

0

1

2

3

4

0

1

2
time,

3

4

1

1

2

0
6

4

0

-1
2

-1

Quotient Graph
2
2

2

0

1

1

2
time,

4

0
-2

3

1

3

Figure 7. Signed external equitable partitions and bipolar clustered consensus dynamics. (a) A signed graph (red
links are positive, blue links are negative) with a sEEP with four cells (indicated by colors). Note how one of the cells (cyan)
contains nodes with different polarizations and another cell (blue) is of negative polarization. The associated quotient graph
is also shown (bottom). (b) Similar to the standard graphs with positive weights, if the full dynamics is given by Hσ y at any
time, then the full dynamics will be exactly determined by the quotient dynamics, but potentially having the opposite sign, like
the node in the blue cell with negative polarization (whose negative trajectory is shown as a dashed line to make this apparent).
(c) The sign-adjusted cell averages hxiσCi = Hσ+ x (crosses) are also determined by the quotient dynamics. Trajectories from
a random initial condition are shown as solid lines. Note how the trajectory of the blue node has the opposite sign to its
sign-adjusted cell average, and in the cyan cell two of the three nodes (with negative polarization) converge to the sign-flipped
value.

(Fig. 7c). The signed consensus dynamics approaches the
bipolar consensus (52), and the final sign of each node
variable is determined by σ T x0 , as seen in Fig. 7c for
nodes 1 − 4 and 8 (positively polarized) and nodes 5, 6, 7
(negatively polarized).
Third, a system with inputs aligned with the sEEP
ẋ = −Lσ x + Hσ u,

(53)

exhibits a bipolar cluster-synchronized state. To see this,
consider the component of the state orthogonal to the
bipolar cluster synchronization manifold:
δ(t)

..=

(I −

Hσ Hσ+ ) x(t),

(54)

which has a dynamics δ̇ = −Lσ δ, and the orthogonal
component decays asymptotically to
σT δ0
σ = 0,
(55)
t→∞
N
since σ T δ 0 = σ T (I − Hσ Hσ+ ) x0 = (σ T − σ T ) x0 = 0.
Hence the system converges to the sEEP manifold.
lim δ(t) =

but now with the additional feature that the dynamics
can support a bipolar cluster synchronization based on
a sEEP. We do not discuss this case in detail again, instead we illustrate these findings for signed Kuramoto
networks.

3. Bipolar cluster synchronization in signed Kuramoto
networks

While standard Kuramoto networks with positive couplings have been studied extensively7 , the literature on
Kuramoto networks with both attractive (positive) and
repulsive (negative) couplings is comparatively sparse,
with only a handful of mean-field results49 .
Using the definition of the signed Laplacian, we write
the Kuramoto model on a signed graph as:
θ̇ = −Bσ sin(BσT θ) = −Bσ W(BσT θ)BσT θ.

(57)

Likewise, the Kuramoto dynamics on the quotient graph
becomes:
2. Bipolar cluster synchronization for nonlinear dynamics
with Laplacian couplings

All the results obtained in Section III B apply to nonlinear dynamics on signed networks of the form:
X
ẋi = F(xi ) − γ
[Lσ ]ij G(xj ),
(56)
j

ψ̇ = −Hσ+ B sin(B T Hσ ψ).

(58)

For structurally balanced signed networks, making the
necessary adjustments for the switching equivalence Σ,
we then reach the same conclusions as in Section III C.
In Figure 8, we provide numerical examples that replicate our findings for signed graphs. As expected, in Fig-

13

(a)

Quotient Graph

Full Graph

model converges to a final state consistent with the cells
of the sEEP, yet exhibiting out-of-phase behavior within
each cell due to the polarization of the nodes.

2
V.

2

(b)

1

0.5

0

-0.5

(c)

0

1

2

3

4

5

6

7

8

0

1

2

3

4

5

6

7

8

time,

1

0.5

0

-0.5

time,

Figure 8. Bipolar cluster-synchronization on a signed
graph of coupled Kuramoto oscillators with identical frequencies. (a) Signed graph (red links positive, blue
links negative) with sEEP indicated by colors. Also shown
is the associated quotient graph. (b) For an initial condition
aligned with the sEEP, the dynamics of the full system (line)
is exactly equivalent to the lower dimensional Kuramoto dynamics on the quotient graph (circles) up to the sign, given
by the polarization of each node. (c) For an initial condition
θ 0 not too spread out on the circle, the sign-adjusted cellaveraged dynamics (crosses) governed by the quotient graph
is closely aligned with the full dynamics. The system converges to a state where nodes within a cell have phases with
opposite signs.

ure 8b) we see that a bipolar cluster-synchronized solution remains invariant for all times. Figure 8c) shows
that for an initial condition θ 0 that is not too spread
out on the unit circle, we observe numerically that the
(sign adjusted) cell-averaged hθiσCi = Hσ+ θ is well aligned
with the quotient dynamics. For non-identical intrinsic
frequencies commensurate with the sEEP, the Kuramoto

DISCUSSION

In this paper, we have shown how to coarse-grain the
dynamics of generic synchronization processes by using
the graph-theoretical framework of external equitable
partitions. Exploiting regularities present in the underlying coupling network, EEPs give clusters of nodes
that play an equivalent dynamical role27 . The resulting
coarse-grained dynamics corresponds to cluster synchronization, in which all nodes within a cell follow the same
trajectory. Importantly, one can extend the notion of
EEPs to other types of coupling schemes, as shown by
our analysis of signed networks. In structurally balanced
signed networks, we showed that each of the cells splits
into two ‘out-of-phase’ dynamical factions, with the same
magnitude but opposite sign.
Connections with symmetry groups. We have
shown how our graph-theoretical approach complements
the use of symmetry groups for the analysis of synchronization dynamics over networks24,25,50 . As discussed
in Pecora et al. 24 , there exist efficient software tools to
compute symmetry groups in networks51,52 and hence orbit partitions24,25 . The more general problem of obtaining all EEPs for a graph appears to be computationally
more challenging25,53 . However, there exist efficient algorithms to compute EEPs centered around a node27,54 ,
which can be used to characterize the dynamical influence of particular nodes on the global dynamics of the
network26,27,29 .
Other signed coupling schemes. We have chosen
to consider couplings given by the signed Laplacian Lσ ,
as it provides a direct generalization of the standard
Laplacian and has direct connections with dynamical
properties. In particular, the positive semi-definiteness
of Lσ allows us to express the Kuramoto model in terms
of signed incidence matrices (57), thus facilitating our
proof and interpretation. However, the ideas developed
here may be applied to other signed coupling schemes,
provided the equivalent invariance condition to (2) can
be found. For instance, another interesting coupling is
given by the Laplacian L± = D − A, where A is a signed
adjacency matrix and D = diag(A1). For the network in
Fig. 7a, the EEP with respect to Laplacian L± is almost
identical to the one with respect to Lσ , except that node
8 forms its own cell. It is worth remarking, however, that
while the algebraic characterization of such invariant partitions can still be exploited, the graph-theoretical notion
of ’equitability’, related to the combinatorial count of
inter-cell degrees, can be lost. As many algorithms leverage such combinatorial properties to search for EEPs,
this may make the presence of such invariant partitions
harder to detect. Furthermore, their dynamical interpretation might be problematic in generic systems since L±

14
(and other coupling matrices) will in general be indefinite, hence impacting the dynamical stability.

ACKNOWLEDGMENTS

Relation to other synchronization notions. Because EEPs exploit the graph structure in connection
with dynamics, our approach complements other methods for the analysis of synchronization. Here, we have
concentrated on the existence and invariance of clustersynchronized states, with some discussion of their stability in the context of the MSF. Further insights could
be gained by combining our EEP-based analysis with
other methodologies, such as the analysis of potential
energy coupling landscapes or various mean-field analyses (see Refs.3,4,7 for overviews). In particular, as
EEPs are linked to the existence of invariant subspaces,
contraction-based arguments may be fruitfully applied to
global system dynamics in such networks55–58 .

We thank Karol Bacik for comments and carefully
reading the manuscript. MTS, JCD, RL acknowledge
support from: FRS-FNRS; the Belgian Network DYSCO
(Dynamical Systems, Control and Optimisation) funded
by the Interuniversity Attraction Poles Programme initiated by the Belgian State Science Policy Office; and the
ARC (Action de Recherche Concerte) on Mining and Optimization of Big Data Models funded by the WalloniaBrussels Federation. NO’C was funded by a Wellcome
Trust Doctoral Studentship at Imperial College London
during this work. YNB thanks the G. Harold & Leila Y.
Mathers Foundation. MB acknowledges support through
EPSRC grants EP/I017267/1 and EP/N014529/1.
No new data was collected in the course of this research.

A particular notion of synchronization worth mentioning is that of chimera states6,59 , in which parts of the
network act in unison while another parts appear unsynchronized. One may conjecture that a possible mechanism to reach such a state is to endow the underlying
network with an EEP comprising one large cell and a
multitude of single node cells. If, by carefully configuring the dynamics, the large cell could be made to remain
stable while the single node cells follow independent trajectories, a chimera state might be obtained.
Future work. Several other avenues of future work
appear to be worth pursuing. As the idea of signed networks and social balance is at the core of social network
theory31,32,60 , it would be important to investigate if the
bipolar cluster synchronization described here can be related to models evolving towards structural balance47,48 .
Following the insight by Hendrickx61 that signed opinion
dynamics can be understood as a 2N dimensional dynamics with positive interactions, it would also be interesting
to understand the symmetry requirements that an EEP
implies on the lifted 2N -dimensional graph, and whether
the EEP could be used to elucidate further properties.
Another extension would be to relax the strict requirements of EEPs (e.g., by allowing minor perturbations on
a graph with an EEP) in order to study how the dynamics
of the system is affected. Generalizations of EEPs that
allow different kinds of couplings (e.g., directed, timevarying62,63 , delays64 ) would also be of interest.
Finally, it is worth remarking that while we focussed
here on the dynamics of synchronization in linear (consensus) and nonlinear processes (coupled oscillators, Kuramoto), the concept of external equitable partitions is
applicable to more general scenarios where agents interact over a graph structure. The key ingredient is the
presence of a low-dimensional invariant subspace in the
coupling (spanned by the partition matrix) which can be
exploited to obtain a dynamical dimensionality reduction
leading to a coarse-grained system description. While
EEPs have been used in consensus and control, other application areas such as ecological networks or chemical
reaction networks would be worth investigating.

REFERENCES
1 S.

H. Strogatz, “From Kuramoto to Crawford: exploring the
onset of synchronization in populations of coupled oscillators,”
Physica D: Nonlinear Phenomena 143, 1–20 (2000).
2 S. Boccaletti, J. Kurths, G. Osipov, D. Valladares, and C. Zhou,
“The synchronization of chaotic systems,” Physics Reports 366,
1–101 (2002).
3 A. Arenas, A. Dı́az-Guilera, J. Kurths, Y. Moreno, and C. Zhou,
“Synchronization in complex networks,” Physics Reports 469,
93–153 (2008).
4 F. Dörfler and F. Bullo, “Synchronization in complex networks of
phase oscillators: A survey,” Automatica 50, 1539–1564 (2014).
5 T. Nishikawa and A. E. Motter, “Comparative analysis of existing
models for power-grid synchronization,” New Journal of Physics
17, 015012– (2015).
6 M. J. Panaggio and D. M. Abrams, “Chimera states: coexistence
of coherence and incoherence in networks of coupled oscillators,”
Nonlinearity 28, R67– (2015).
7 F. A. Rodrigues, T. K. D. Peron, P. Ji, and J. Kurths, “The
Kuramoto model in complex networks,” Physics Reports 610,
1–98 (2016).
8 Y. Kuramoto, Chemical oscillations, waves, and turbulence,
Vol. 19 (Springer Science & Business Media, 2012).
9 A. Jadbabaie, J. Lin, and A. S. Morse, “Coordination of groups
of mobile autonomous agents using nearest neighbor rules,” IEEE
Transactions on Automatic Control 48, 988–1001 (2003).
10 J. Heagy, T. Carroll, and L. Pecora, “Synchronous chaos in
coupled oscillator systems,” Physical Review E 50, 1874 (1994).
11 L. M. Pecora and T. L. Carroll, “Master stability functions for
synchronized coupled systems,” Physical Review Letters 80, 2109
(1998).
12 M. Barahona and L. M. Pecora, “Synchronization in Small-World
Systems,” Phys. Rev. Lett. 89, 054101 (2002).
13 L. M. Pecora and M. Barahona, “Synchronization of oscillators
in complex networks,” Chaos and Complexity Letters 1, 61–91
(2005).
14 A. Jadbabaie, N. Motee, and M. Barahona, “On the stability of the Kuramoto model of coupled nonlinear oscillators,” in
Proceedings of the American Control Conference, 2004., Vol. 5
(2004) pp. 4296 –4301.
15 K. Judd, “Networked dynamical systems with linear coupling:
synchronisation patterns, coherence and other behaviours,”
Chaos: An Interdisciplinary Journal of Nonlinear Science 23,
043112 (2013).
16 O. D’Huys, R. Vicente, T. Erneux, J. Danckaert, and I. Fis-

15
cher, “Synchronization properties of network motifs: Influence of
coupling delay and symmetry,” Chaos 18, 037116 (2008).
17 W. Lu, B. Liu, and T. Chen, “Cluster synchronization in networks of coupled nonidentical dynamical systems,” Chaos: An
Interdisciplinary Journal of Nonlinear Science 20, 013120 (2010).
18 T. Dahms, J. Lehnert, and E. Schöll, “Cluster and group synchronization in delay-coupled networks,” Physical Review E 86,
016202 (2012).
19 A.-L. Do, J. Höfener, and T. Gross, “Engineering mesoscale
structures with distinct dynamical implications,” New Journal
of Physics 14, 115022– (2012).
20 C. Fu, Z. Deng, L. Huang, and X. Wang, “Topological control
of synchronous patterns in systems of networked chaotic oscillators,” Phys. Rev. E 87, 032909– (2013).
21 D. P. Rosin, D. Rontani, D. J. Gauthier, and E. Schöll, “Control
of Synchronization Patterns in Neural-like Boolean Networks,”
Phys. Rev. Lett. 110, 104102– (2013).
22 F. Sorrentino and E. Ott, “Network synchronization of groups,”
Phys. Rev. E 76, 056114– (2007).
23 C. R. S. Williams, T. E. Murphy, R. Roy, F. Sorrentino,
T. Dahms, and E. Schöll, “Experimental Observations of Group
Synchrony in a System of Chaotic Optoelectronic Oscillators,”
Phys. Rev. Lett. 110, 064104– (2013).
24 L. M. Pecora, F. Sorrentino, A. M. Hagerstrom, T. E. Murphy,
and R. Roy, “Cluster synchronization and isolated desynchronization in complex networks with symmetries,” Nature communications 5 (2014).
25 F. Sorrentino, L. M. Pecora, A. M. Hagerstrom, T. E. Murphy,
and R. Roy, “Complete characterization of the stability of cluster
synchronization in complex dynamical networks,” Science Advances 2 (2016), 10.1126/sciadv.1501737.
26 M. Egerstedt, S. Martini, M. Cao, K. Camlibel, and A. Bicchi, “Interacting with Networks: How Does Structure Relate to
Controllability in Single-Leader, Consensus Networks?” Control
Systems, IEEE 32, 66–73 (2012).
27 N. O’Clery, Y. Yuan, G.-B. Stan, and M. Barahona, “Observability and coarse graining of consensus dynamics through the external equitable partition,” Physical Review E 88, 042805 (2013).
28 D. M. Cardoso, C. Delorme, and P. Rama, “Laplacian eigenvectors and eigenvalues and almost equitable partitions,” European
Journal of Combinatorics 28, 665–673 (2007).
29 S. Martini, M. Egerstedt, and A. Bicchi, “Controllability analysis of multi-agent systems using relaxed equitable partitions,”
International Journal of Systems, Control and Communications
2, 100–121 (2010).
30 Y. Kuramoto, “Self-entrainment of a population of coupled nonlinear oscillators,” in International Symposium on Mathematical
Problems in Theoretical Physics (Springer, 1975) pp. 420–422.
31 F. Heider, “Attitudes and cognitive organization,” The Journal
of Psychology 21, 107–112 (1946).
32 D. Cartwright and F. Harary, “Structural balance: a generalization of Heider’s theory.” Psychological Review 63, 277 (1956).
33 C. Altafini, “Consensus Problems on Networks With Antagonistic Interactions,” Automatic Control, IEEE Transactions on 58,
935–946 (2013).
34 C. Altafini, “Stability analysis of diagonally equipotent matrices,” Automatica 49, 2780 – 2785 (2013).
35 A. Chan and C. D. Godsil, “Symmetry and eigenvectors,” in
Graph Symmetry (Springer, 1997) pp. 75–106.
36 S. Kudose, “Equitable Partitions and Orbit Partitions,”
http://www.math.uchicago.edu/~may/VIGRE/VIGRE2009/
REUPapers/Kudose.pdf (2009).
37 O. E. Rössler, “An equation for continuous chaos,” Physics Letters A 57, 397–398 (1976).
38 F. Dörfler and F. Bullo, “Synchronization and transient Stability
in Power Networks and nonuniform Kuramoto oscillators,” SIAM
Journal of Control and Optimization 50, 1616–1642 (2012).
39 V. Nicosia, M. Valencia, M. Chavez, A. Dı́az-Guilera,
and
V. Latora, “Remote synchronization reveals network symmetries
and functional modules,” Physical Review Letters 110, 174102

(2013).
Leskovec, D. Huttenlocher, and J. Kleinberg, “Signed networks in social media,” in Proceedings of the SIGCHI conference
on human factors in computing systems (ACM, 2010) pp. 1361–
1370.
41 M. Szell, R. Lambiotte, and S. Thurner, “Multirelational organization of large-scale social networks in an online world,” Proceedings of the National Academy of Sciences 107, 13636–13641
(2010).
42 E. Davidson and M. Levin, “Gene regulatory networks,” Proceedings of the National Academy of Sciences 102, 4935–4935
(2005).
43 W. Gerstner, W. M. Kistler, R. Naud, and L. Paninski, Neuronal dynamics: From single neurons to networks and models of
cognition (Cambridge University Press, 2014).
44 M. T. Schaub, Y. N. Billeh, C. A. Anastassiou, C. Koch, and
M. Barahona, “Emergence of slow-switching assemblies in structured neuronal networks,” PLoS Comput Biol 11, e1004196
(2015).
45 J. Kunegis, S. Schmidt, A. Lommatzsch, J. Lerner, E. W.
De Luca, and S. Albayrak, “Spectral analysis of signed graphs
for clustering, prediction and visualization.” in SDM, Vol. 10
(SIAM, 2010) pp. 559–559.
46 G. Facchetti, G. Iacono, and C. Altafini, “Computing global
structural balance in large-scale signed social networks,” Proceedings of the National Academy of Sciences 108, 20953–20958
(2011).
47 S. A. Marvel, J. Kleinberg, R. D. Kleinberg, and S. H. Strogatz,
“Continuous-time model of structural balance,” Proceedings of
the National Academy of Sciences 108, 1771–1776 (2011).
48 V. A. Traag, P. Van Dooren, and P. De Leenheer, “Dynamical
models explaining social balance and evolution of cooperation,”
PloS one 8, e60063 (2013).
49 A. El-Ati and E. Panteley, “Phase locked synchronization for Kuramoto model, with attractive and repulsive interconnections,”
in Systems, Man, and Cybernetics (SMC), 2013 IEEE International Conference on (IEEE, 2013) pp. 1253–1258.
50 M. Barahona, E. Trias, T. P. Orlando, A. E. Duwel, H. S.
van der Zant, S. Watanabe, and S. H. Strogatz, “Resonances
of dynamical checkerboard states in Josephson arrays with selfinductance,” Physical Review B 55, 11989–11992 (1997).
51 P. T. Darga, K. A. Sakallah, and I. L. Markov, “Faster symmetry
discovery using sparsity of symmetries,” in Proceedings of the
45th annual Design Automation Conference (ACM, 2008) pp.
149–154.
52 T. G. Group, “GAP: Groups, Algorithms, and Programming,”
http://www.gap-system.org.
53 H. Kamei and P. J. Cock, “Computation of balanced equivalence
relations and their lattice for a coupled Cell network,” SIAM
Journal on Applied Dynamical Systems 12, 352–382 (2013).
54 M. Cao, S. Zhang, and M. K. Camlibel, “A Class of Uncontrollable Diffusively Coupled Multiagent Systems with Multichain
Topologies,” IEEE Transactions on Automatic Control 58, 465–
469 (2013).
55 L. Moreau, “Stability of continuous-time distributed consensus
algorithms,” in 43rd IEEE Conference on Decision and Control
(CDC) 2004, Vol. 4 (IEEE, 2004) pp. 3998–4003.
56 R. Sepulchre, “Consensus on nonlinear spaces,” Annual Reviews
in Control 35, 56–64 (2011).
57 G. Russo and J.-J. E. Slotine, “Symmetries, stability, and control in nonlinear systems and networks,” Physical Review E 84,
041929 (2011).
58 A. Mauroy and R. Sepulchre, “Contraction of monotone phasecoupled oscillators,” Systems & Control Letters 61, 1097–1102
(2012).
59 D. M. Abrams and S. H. Strogatz, “Chimera States for Coupled
Oscillators,” Phys. Rev. Lett. 93, 174102 (2004).
60 P. Doreian and A. Mrvar, “A partitioning approach to structural
balance,” Social networks 18, 149–168 (1996).
61 J. Hendrickx, “A lifting approach to models of opinion dynamics
40 J.

16
with antagonisms,” in Decision and Control (CDC), 2014 IEEE
53rd Annual Conference on (2014) pp. 2118–2123.
62 L. Scardovi and R. Sepulchre, “Synchronization in networks of
identical linear systems,” Automatica 45, 2557–2562 (2009).
63 D. A. Paley, N. E. Leonard, R. Sepulchre, D. Grünbaum, and

J. K. Parrish, “Oscillator models and collective motion,” Control
Systems, IEEE 27, 89–105 (2007).
64 A. Papachristodoulou, A. Jadbabaie, and U. Munz, “Effects of
delay in multi-agent consensus and oscillator synchronization,”
Automatic Control, IEEE Transactions on 55, 1471–1477 (2010).

