Synchronization and resonance on networks are some of the most remarkable collective dynamical phenomena. The network topology, or the nature and distribution of the connections within an ensemble of coupled oscillators, plays a crucial role in shaping the local and global evolution of the two phenomena. This article further explores this relationship within a compact mathematical framework and provides new contributions on certain pivotal issues, including a closed bound for the average synchronization time in arbitrary topologies; new evidences of the effect of the coupling strength on this time; exact closed expressions for the resonance frequencies in terms of the eigenvalues of the Laplacian matrix; a measure of the effectiveness of an *influencer node*’s impact on the network; and, finally, a discussion on the existence of a resonant synchronized state. Some properties of the solution of the linear swing equation are also discussed within the same setting. Numerical experiments conducted on two distinct real networks—a social network and a power grid—illustrate the significance of these results and shed light on intriguing aspects of how these processes can be interpreted within networks of this kind.

Why do pink flamingos perform a pre-breeding dance in perfect unison? How does the entire Australian coral reef manage to synchronize to release millions of gametes at the same time each year? How does it happen that a number of individuals in a social network take part with increasing involvement in collective ritual events? These questions concern the problem of how multiple agents within a network can achieve an identical synchronized dynamic behavior, possibly with an amplification of the intensity of the phenomenon over time. Prompted by the ubiquity and transversality of these processes, the paper offers some novel insights into the already immense research field on the dynamics of networks of coupled identical harmonic oscillators.

## I. INTRODUCTION

Up to several thousand great flamingos, during the pre-breeding period, gather in southern France to perform in synchrony a variety of stereotyped movements. Such group performances, involving elaborate courtship rituals and synchronized dances, are believed to play a role in enabling reproductive synchronization.^{1} In the Australia’s Great Barrier Reef, many coral species tune their reproductive cycles with incredible precision, synchronizing each year all on the same day, just after the full moon in November, to collectively release millions of gametes in one of the largest scale mass reproduction events on the planet.^{2} In the epidemiological context, synchronization processes govern the dynamics of disease spread through complex networks in anthropic environments. Examples are the dynamics of childhood diseases, such as measles, mainly determined by the annual seasonal cycle, or the pulse vaccination-induced synchronization in simple epidemic models. It is a peculiar feature that ecological communities tend to adjust to periodic environmental or climatic driving variables. Although ecologists and epidemiologists have essentially opposite goals, the mathematical laws underlying the different population dynamics are very similar and frequently involve internal synchronization or resonance with external driving factors.^{3} We find countless examples in different settings. In the biomechanical context, Brumley *et al.*^{4} study a mechanism of self-propulsion in ancestral eukaryotes, unicellular organisms equipped with tiny mobile flagella that allow them to move in a liquid medium. Correlated changes in the directions of the flagellar beating can be explained by the raise of some kind of large scale coordination and synchronization, even in a weakly coupled system as is the case for the green alga *Volvox*. In the human social setting, the emergence of synchronization processes among individuals, on a physical or mental level, can be induced by collective dynamics and participation in repetitive rituals or dances, as occurs in Zikr, a Sufi mystical dance performed by Chechen Sufis and involving hundreds of men at once. However, it is well known that synchronization phenomena do not always have positive implications. Strong synchronization may indeed be related to pathological activities as is the case for the epileptic seizures and Parkinson’s disease. Hammond *et al.*^{5} have observed an abnormally synchronized oscillatory activity at multiple levels of the basal ganglia-cortical loop in humans with Parkinson’s disease. This excessive synchronization correlates with motor deficit, and its suppression by specific therapies could suggest an effective approach to improve the living conditions of patients with severe motility disorders.

These examples show that synchronization emerges, in general, from the collaboration, competition, and interaction of multiple agents when they undergo periodic and cyclic dynamics. A first interpretive model, albeit oversimplified, reduces such systems to a set of coupled identical harmonic oscillators all synchronizing on the same frequency.

The two pioneering papers on synchronization in coupled systems^{6} and synchronization in chaotic systems^{7} have stimulated a great deal of interest in the study of complete synchronization of coupled nonlinear dynamical processes. The first paradigmatic model of diffusively coupled identical oscillators on a complex network is given in Pecora *et al.*,^{8} where the master stability function is introduced for the first time. Such a function accomplishes two goals: decrease the computational load and provide a guidance in designing coupling configurations that conform to the required stability. This approach is currently the most widely applied one to determine the stability of the fully synchronized state in terms of the eigenstructure of the connectivity matrix.

The interest in the study of synchronization in coupled systems grew considerably after the discovery of the small-world and scale-free properties of many natural and artificial complex networks. Barahona and Pecora^{9} analyze the implications of the small-world phenomenon on the synchronization of oscillator networks of arbitrary topology. Subsequent analysis provided evidence that well-connected networks are easier to globally synchronize than regular networks. For instance, Chen *et al.*^{10} found that synchronization processes typically start first from a small fraction of hub nodes and then spread to nodes with smaller degrees. This is replicated at the cluster level, where a typical synchronization process starts from partial synchronization within clusters to evolve to global complete synchronization.

Ren^{11} analyzes coupled second-order linear harmonic oscillators, described as point masses, under rather mild network connectivity assumptions, showing that their positions and velocities can be synchronized in the presence or absence of a leader. A similar approach is used in Zhan *et al.*^{12} where the vibration spectra of a classical mass spring model on complex networks are investigated; and, in Shang,^{13} where synchronization of a leader–follower system of coupled harmonic oscillators connected by dampers is studied in the presence of random noises and time delays with an interaction topology modeled by a weighted directed graph.

In general, the presence of dissipative couplings has the effect to equate the state of the nodes so that the system enters a positively invariant set, the absorbing domain, in finite time. It is well-known that if $G$ is a constant matrix whose eigenvalues have only negative real parts, then $ e t G$ acts as a uniform contraction; that is, $ | | e t G | |\u2264K e \u2212 \eta t$ for some positive real constants, $K$ and $\eta $. This fact guarantees that the solution of the equation is uniformly asymptotically stable, and the synchronization manifold remains invariant under the flow of the equation that governs the system.

For numerical estimates and model prediction purposes, Pietras and Daffertshofer^{14} stress the importance of reduction of complex oscillatory systems. A number of surveys are devoted to this topic.^{15} Arola-Fernández *et al.*^{16} study the synchronized state in a population of network-coupled heterogeneous oscillators, showing that the steady-state solution of the linearized dynamics may be written as a geometric series whose terms represent the different spatial scales of the network. They prove that this expansion converges for arbitrary frequency distributions and for both undirected and directed networks, provided that the adjacency matrix is primitive. Also, Eroglu *et al.*^{17} propose a survey on the main theory behind complete, generalized, and phase synchronization phenomena in complex networks with an approach close to the one followed in the present paper.

Investigations of the transient and steady-state dynamics of synchronous machines have significantly impacted the description of power grid networks. Bhatta *et al.*^{18} study the swing equation and its linearization on transmission networks, highlighting the role of symmetries in reducing the computational complexity of modal decomposition.^{19} Bhatta *et al.*^{20} extend this approach to multilayer networks.

A crucial issue is the time and speed of synchronization. Grabow *et al.*^{21} focus on the speed of convergence as a collective time scale for synchronizing systems on different topologies, ranging from completely ordered and grid-like to completely disordered and random topologies. They find, for instance, that, at fixed in-degree, topological randomness induces shorter synchronization times, whereas intermediate randomness in the small-world regime induces longer synchronization times.

Many other contributions have been given on specific issues, such as in Skardal *et al.*^{22} and Skardal *et al.*,^{23} where authors introduce a synchrony alignment function that encodes the interplay between the network structure and oscillator frequencies as an objective measure for the synchronization properties of a network of heterogeneous oscillators with natural frequencies or in Su *et al.*,^{24} where authors revisit the synchronization problem for coupled oscillators in a dynamic proximity network. Dörfler *et al.*^{25} prove a closed-form condition for synchronization that can be stated in terms of the Moore–Penrose pseudoinverse of the Laplacian matrix. Liu and Barabási,^{26} in the wider context of control theory, address the problem of how networks organize themselves to balance control with functionality. In particular, they review the process of local pinning synchronization, when controllers, applied to a subset of nodes, help synchronize the network. Also, in Lu *et al.*,^{27} pinning synchronization in the new setting of networks of networks is investigated. Slotine *et al.*^{28} study the role, in both nature and system design, of co-existing power leaders, to which the networks synchronize, and knowledge leaders, to whose parameters the networks adapt.

One of the most important contributions to the study of synchronization in coupled systems is undoubtedly the Kuramoto model.^{29} In this model, the underlying topology is fully connected, even if it was later adapted to the scenario of nearest neighbor interactions. Despite the simplicity and effectiveness of the Kuramoto model and all its variants, several works have highlighted some of its shortcomings. For example, Shahal *et al.*^{30} address the issue of the importance of synchronization in human networks for our civilization. They study, in an accurate experimental setting, the synchronization between violin players, arranged in a network framework with full control in terms of network connectivity, coupling strength, and delay. Their results show that players can adjust their playing period and cancel connections by ignoring frustrating signals to find a stable solution. Their observations reveal that the usual models for coupled networks, such as the Kuramoto model, cannot always be applied to human networks, as additional degrees of freedom enable new strategies and produce better solutions than are possible within such a model. Corroborating these observations, it is argued that, during social interactions, participants are continuously active, each modifying their actions in response to the changing actions of their partners. This continuous mutual adaptation results in an interactional synchrony to which all members contribute. Dumas *et al.*^{31} study the brain activity of the participants in this process and found that states of interactional synchrony correlate with the emergence of an interbrain synchronization network between right centroparietal regions. These regions have been found to play a central role in social interaction within an interindividual brainweb.

Another limiting aspect of classical models is the presence of strictly local interactions. For this reason, Estrada *et al.*^{32} study how long-range interactions can affect synchronization dynamics, proposing a model of coupled oscillators based on the d-path Laplacian matrices. Their analysis reveals that long-range interactions improve network synchronizability with an impact that depends on the original structure.

The dynamics of opinion formation within a social network is also a widely studied field. In this kind of network, each node can be assigned, for instance, a discrete or continuous score according to his/her agreement or disagreement on a given topic. By introducing a coupling between members in order to describe the dragging effect in their opinion formation process, each member in the network can change his/her opinion according to a synchronization mechanism or according to some kind of resonant phenomenon with some source/leader status.

Related to this issue, Jenkins^{33} presents a very interesting analysis of the so-called self-oscillation phenomenon, which is distinct from the more familiar phenomenon of forced resonance, whereas Tönjes *et al.*^{34} focus on networks with influencers, a group of hubs that couple strongly and connect different parts of the network. They show that subjecting influencers to an optimal noise can result in enhanced network synchronization. This effect, called coherence resonance, emerges from a synergy between a network structure and stochasticity and is highly nonlinear, vanishing as the noise is too weak or too strong.

In the light of the previous literature review, this paper aims to analyze, in a sense from first principles, some aspects of the synchronization and resonance phenomena and, in general, of the transient dynamics in networks of oscillators with a heterogeneous topological structure. A unified approach is maintained in the description of all phenomena. The main new contributions and results of the present work can be summarized as follows:

a synthetic framework in which to inscribe the description of synchronization and resonance phenomena on networks of coupled identical harmonic oscillators with arbitrary topology;

a closed expression of the average synchronization time that encodes information about the topological structure of the network;

new evidences that, for a weakly coupled system, complete and dense networks achieve the state of perfect synchronization more rapidly than, for instance, the path network and that the complete network shows a synchronization time, which is shorter for weak coupling than for strong coupling;

closed expressions for the resonance frequencies of the network as a function of the eigenvalues of the Laplacian matrix;

a detailed account of how a single node in the network can serve as a leader (influencer) or a follower (receiver) depending on its topological location when a resonance phenomenon arises inside the network;

a proof that the presence of dissipative terms pulls down the possible resonance frequencies to a single one, corresponding to the synchronization frequency;

a new perspective on the linearization of the swing equation and its solution in the proposed setting; and

some new mathematical properties of the polar decomposition of the symplectic matrix governing the synchronization phenomenon.

This paper is organized as follows. Section II introduces the basic notations and the general model. Section III has a didactic purpose: it shortly illustrates the case of networked coupled identical linear oscillators without damping. In Sec. IV, existing results are framed into the broader context described, and new ones about the mathematical properties of the involved operators and the rate of synchronization are presented. Section V is devoted to resonance phenomena and, after some recalls on the case of the single driven oscillator, provides new results related to resonant frequencies on networks when an external driving force is placed at a specific node. Section VI deals with a special case of the general equation introduced, specifically the linear swing equation and its solution. Finally, Sec. VII applies the previous findings to two real world networks, the Zachary club social network and the Syrian power grid, by illustrating a number of interesting remarks. All the proofs of the propositions in this article are given in Appendixes A and B.

## II. PRELIMINARIES

Our purpose is basically to describe the behavior of a system of $n$ coupled harmonic oscillators within a network structure. To model the local interactions between the $n$ harmonic oscillators, we use an undirected graph $G=(V,E)$ without loops, where $V$ is the set of vertices or nodes and $E\u2282V\xd7V$ is the set of edges or links. Throughout the paper, we will assume that the network is connected, i.e., that a path exists between each pair of distinct vertices. The actual configuration is completely encoded in the adjacency matrix $A={ a i j},i,j=1,\u2026,n$, where $ a i j=1$ when $(i,j)\u2208E$, that is, when nodes $i$ and $j$ show an actual local interaction, and $0$ otherwise. We will assume $ a i j= a j i$, which is a perfectly symmetric interaction within a pair of nodes. The degree of a node $i$ is defined as $ k i= \u2211 j = 1 n a i j= ( A 2 ) i i$, and $K= diag{ k i}$ is the diagonal matrix of the degrees of all the nodes. We will denote by $ \lambda i$ and $ \psi i,i=1,\u2026,n$, respectively, the eigenvalues and eigenvectors of the adjacency matrix $A$, where the eigenvalues are listed in decreasing order: $ \lambda 1\u2265\cdots \u2265 \lambda n$. The Laplacian matrix is then defined as $L=K\u2212A$. We will denote by $ \mu i$ and $ \varphi i,i=1,\u2026,n$, respectively, the eigenvalues and eigenvectors of the Laplacian matrix $L$, where again, the eigenvalues are listed in decreasing order: $ \mu 1> \mu 2>\cdots > \mu n=0$. In particular, $ \varphi n= 1 n [ 1 , 1 , \u2026 , 1 ] T= 1 n u T$, being $u\u2208 R n$ the vector of all $1$’s. Similarly, we will denote by $1=u u T$ the $n\u2212$square matrix of all $1$’s and by $I$ the identical matrix.

*identical*harmonic oscillators so that $M=I$. Hence, the system (7) becomes

*leader*. The total force on node $i$ is now

Finally, in addition to the cases described above, we will study a further variant of the system (12) in order to show how this formalism can accommodate the linearization of an important nonlinear problem. In particular, by setting $ c 1=0$, $ c 2=1$, $ c 1 \u2032=1$, $ c 2 \u2032=0$, and by giving the factor $b(t)$ in Eq. (12), the form $b(t)= [ 0 , p ] T$, where $p$ is a constant vector in $ R n$, we will gain the *linear swing equation* describing the dynamic response of a network of synchronous machines to small disturbances.

Shortly, we can summarize in Table I the cases we are interested in, where $ c 1$, $ c 2$, $ c 1 \u2032$, $ c 2 \u2032$ are the coupling terms in the matrix $G$ in Eq. (9) and the last column specifies the form of the driving term $b(t)$ in Eq. (12).

Network system . | c_{1}
. | c_{2}
. | $ c 1 \u2032$ . | $ c 2 \u2032$ . | b(t)
. |
---|---|---|---|---|---|

Coupled | 1 | 1 | 0 | 0 | 0 |

Coupled and damped | 1 | 0 | 0 | 1 | 0 |

Coupled and forced | 1 | 1 | 0 | 0 | f(t)[0, e_{h}]^{T} |

Coupled, damped, and forced | 1 | 0 | 0 | 1 | f(t)[0, e_{h}]^{T} |

Linear swing equation | 0 | 1 | 1 | 0 | [0, p]^{T} |

Network system . | c_{1}
. | c_{2}
. | $ c 1 \u2032$ . | $ c 2 \u2032$ . | b(t)
. |
---|---|---|---|---|---|

Coupled | 1 | 1 | 0 | 0 | 0 |

Coupled and damped | 1 | 0 | 0 | 1 | 0 |

Coupled and forced | 1 | 1 | 0 | 0 | f(t)[0, e_{h}]^{T} |

Coupled, damped, and forced | 1 | 0 | 0 | 1 | f(t)[0, e_{h}]^{T} |

Linear swing equation | 0 | 1 | 1 | 0 | [0, p]^{T} |

The first case in Sec. III is for illustrative purposes only, and it is used to properly introduce notations. The second case in Sec. IV will lead to synchronization phenomena on the network. The third and fourth cases in Sec. V are the core of the paper and lead to the resonance phenomena on networks. Finally, the fifth case in Sec. VI will allow us to discuss some properties of the linear swing equation and its solution.

## III. NETWORK OF COUPLED IDENTICAL HARMONIC OSCILLATORS

Let us denote $H= c 1I+ c 2L$ and $B$ the square root matrix of $H$: $ B 2=H$. Since, for any value of the positive coefficients $ c 1$ and $ c 2$, the matrix $H$ is positive definite, $B$ can be easily obtained as $B=\Phi \Lambda 1 / 2 \Phi T$, where $\Phi $ is the matrix whose columns are the eigenvectors of $H$ and $L$ and $\Lambda $ is the diagonal matrix of the positive eigenvalues of $H$.

Figure 2 illustrates the harmonic motion of the four nodes of the network in Fig. 1 given by solution $ x i(t)$ in Eq. (22). The values for the two coefficients are $ c 1= c 2=1$, and the initial conditions are $ x 0= [ 1 , 0 , 0 , 0 ] T$ and $ v 0= [ 0 , 0 , 0 , 0 ] T$. On the horizontal axis, time ranges in $[0,2\pi ]$.

Let us remark that, if we set $ c 1= c 1 \u2032= c 2 \u2032=0$ and $ c 2=1$ in system (8), we have a network of free harmonic oscillators connected by identical springs.

Therefore, there are $ \u221e 1$ values for the imaginary part $b$ of the initial coefficients. Let us notice that $ B \u2212 1$ does not exist because of the null eigenvalue $ \mu n=0$ of the Laplacian $L$ and then of the matrix $B$. This eigenvalue is associated with the global translation of all the nodes with the same velocity, i.e., to the stationary state of the whole network, which, in this case, is free and not bounded to any fixed support. Let us also observe that since $ \u2211 j B i j= \u2211 i B i j=0$, i.e., the sum of the elements along a row or a column in the matrix $B$ is null, the initial velocity vector $ v 0$ has to satisfy $ \u2211 i v 0 i=2 \u2211 i \u2211 j B i j b j=2 \u2211 j \u2211 i B i j b j=0$. Then, for the system (16) to be consistent, the initial velocity vector must satisfy the condition $ \u2211 i v 0 i=0$, equivalent to the classical linear momentum conservation.

## IV. SYNCHRONIZATION OF COUPLED IDENTICAL HARMONIC OSCILLATORS

### A. Mathematical properties

As a consequence, we have that if $\lambda $ is an eigenvalue of $G$, so are its complex conjugate $ \lambda \u22c6$, $ 1 \lambda $, and, hence, $ 1 \lambda \u22c6$. Moreover, if $\lambda $ has multiplicity $m$, then $ 1 \lambda $ has multiplicity $m$. Finally, both $G$ and $ G \u2212 1=\u2212J G TJ$ have the same set of eigenvalues. In particular, denoting by $ \lambda i G \xb1,i=1,\u2026,n$, the eigenvalues of $G$ and $ \psi i G \xb1$ the corresponding eigenvectors, Ren^{11} proves the following proposition:

^{11}for the directed case. Here, we focus on the following remarks.

Proposition 1 has the following straightforward consequence: $G$ is a Hurwitz matrix; i.e., all its eigenvalues are negative or have a negative real part: $ Re( \lambda i G \xb1)<0,\u2200i=1,\u2026,n$. Indeed, we may have the following four possibilities as in Table II.

If μ_{i} > 2 | $ \lambda i G\u2208 R$, | $ \lambda i G \xb1<0$ |

If μ_{i} = 2 | $ \lambda i G\u2208 R$, | $ \lambda i G \xb1=\u22121$ |

If 0 < μ_{i} < 2 | $ \lambda i G\u2208 C$, | $\u22121< Re( \lambda i G \xb1)=\u2212 \mu i 2<0$ |

For μ_{n} = 0 | $ \lambda n G\u2208 I$, | $ \lambda n G \xb1=\xb1i$ |

If μ_{i} > 2 | $ \lambda i G\u2208 R$, | $ \lambda i G \xb1<0$ |

If μ_{i} = 2 | $ \lambda i G\u2208 R$, | $ \lambda i G \xb1=\u22121$ |

If 0 < μ_{i} < 2 | $ \lambda i G\u2208 C$, | $\u22121< Re( \lambda i G \xb1)=\u2212 \mu i 2<0$ |

For μ_{n} = 0 | $ \lambda n G\u2208 I$, | $ \lambda n G \xb1=\xb1i$ |

The next proposition states that, in terms of asymptotic behaviors, the evolution operators $ e G t$ and $ e U t$ are equivalent.

### B. Example

Let us consider again the toy model in Fig. 1. The four eigenvalues of $L$ are $ \mu 1=4, \mu 2=3, \mu 3=1, \mu 4=0$. The eight eigenvalues of operators $G$ , $P$ , $U$ and the four angles defined in Eq. (A7) in Appendix A are listed in Table III.

λ
. | $G$ . | $P$ . | $U$ . | Angles . |
---|---|---|---|---|

$ \lambda 1 +$ | −0.268 | 4.236 | −0.894 + 0.447i | θ_{1} = 153° |

$ \lambda 2 +$ | −0.382 | 3.303 | −0.832 + 0.555i | θ_{2} = 146° |

$ \lambda 3 +$ | −0.500 + 0.866i | 1.618 | −0.447 + 0.894i | θ_{3} = 117° |

$ \lambda 4 +$ | +i | 1 | +i | θ_{4} = 90° |

$ \lambda 4 \u2212$ | −i | 1 | −i | … |

$ \lambda 3 \u2212$ | −0.500 − 0.866i | 0.618 | −0.447 − 0.894i | … |

$ \lambda 2 \u2212$ | −2.618 | 0.303 | −0.832 − 0.555i | … |

$ \lambda 1 \u2212$ | −3.732 | 0.236 | −0.894 − 0.447i | … |

λ
. | $G$ . | $P$ . | $U$ . | Angles . |
---|---|---|---|---|

$ \lambda 1 +$ | −0.268 | 4.236 | −0.894 + 0.447i | θ_{1} = 153° |

$ \lambda 2 +$ | −0.382 | 3.303 | −0.832 + 0.555i | θ_{2} = 146° |

$ \lambda 3 +$ | −0.500 + 0.866i | 1.618 | −0.447 + 0.894i | θ_{3} = 117° |

$ \lambda 4 +$ | +i | 1 | +i | θ_{4} = 90° |

$ \lambda 4 \u2212$ | −i | 1 | −i | … |

$ \lambda 3 \u2212$ | −0.500 − 0.866i | 0.618 | −0.447 − 0.894i | … |

$ \lambda 2 \u2212$ | −2.618 | 0.303 | −0.832 − 0.555i | … |

$ \lambda 1 \u2212$ | −3.732 | 0.236 | −0.894 − 0.447i | … |

In Figures 3–6, we illustrate some different synchronization processes on the same network, with different initial conditions, as specified in the caption of the figures. For each initial condition, plots represent the position and velocity of the four nodes and the phase portrait with the corresponding limiting cycle.

### C. Synchronization time

The preceding analysis gives us the opportunity to add some insights into the widely discussed topic of synchronization time and speed. What follows contributes to this topic by highlighting the role of network topology in the synchronization process in order to compare the behavior of different systems and to discriminate which systems converge more or less rapidly to the synchronized state.

In the literature, a measure of synchronization is defined as $ C(t)= 1 n \u2211 i = 1 n \Theta ( \epsilon \u2212 d i ( t ) )$, where $ d i(t)= | | y i(t)\u2212 y ~ i(t) | |$, $ y ~(t)$ is the synchronization solution in Eq. (31), $\epsilon $ is a threshold, and $ \Theta $ is the step function (see, for instance, Grabow *et al.*^{21}). This index provides the percentage of nodes whose distance from the synchronization solution is less than a given threshold. As the network becomes completely synchronized on this solution, $ C(t)\u21921$.

There is actually no univocal definition of an index that measures the speed of synchronization. For example, it could also be taken as an indicator of the time at which the maximum of the absolute difference between couples of the oscillator positions, that is, $ max i j | x i ( t ) \u2212 x j ( t ) |$, or, alternatively, the mean of all the differences between couples of positions, drops below a given threshold.

In light of the closed expression of the asymptotic solution in Eq. (31), we will refer directly to the mean absolute difference between positions of nodes and such a limiting solution.

For $n$ fixed then, as $\delta $ grows, then the synchronization time decreases.

Let us test the synchronization time in Eq. (37) on the toy networks with $n=4$ and $n=5$ nodes in Fig. 7.

For the networks A, B, and C with $n=4$, the values of the parameter $ \lambda S$ are $ \lambda S A=\u22120.2928932$, $ \lambda S B=\u22120.2679492(\u22120.3819660)$, and $ \lambda S C=\u22120.2679492(\u22120.2679492)$.^{35} Then, we expect that network A would synchronize faster than B, which in turn would synchronize faster than C.

In Fig. 8, panel (a), we plot the mean of the absolute values of differences between node position at time $t$ and the synchronization solution, as in Eq. (32), as a function of time $t$: we observe that the red line (network A) decays on average more rapidly than the blue line (network B), which decays more rapidly than the green line (network C). If we compute the mean synchronization times given in formula (37) with $\epsilon = 10 \u2212 3$, we get $ t mean A=8.421$, $ t mean B=17.273$, and $ t mean C=18.264$. The empirical values, computed as the mean of the minimum times $t$ such that the absolute value of the difference between the solution for a given node and the asymptotic solution goes below $\u03f5$, are $ t emp A=8.275$, $ t emp B=17.525$, and $ t emp C=21.925$.

Similarly, for networks D–F with $n=5$, the values of the parameter $ \lambda S$ are $ \lambda S A=\u22120.1909830$, $ \lambda S B=\u22120.2087122(\u22120.2679492)$, and $ \lambda S C=\u22120.2087122(\u22120.2087122)$. Then, we expect that network D would synchronize faster than E, which in turn would synchronize faster than F. In Fig. 8, panel (b), we plot again the mean deviation of the node positions at time $t$ and the synchronization solution, as in Eq. (32), as a function of time $t$: we observe that the red line (network D) decays on average more rapidly than the blue line (network E), which decays more rapidly than the green line (network F). Again, if we compute the mean synchronization times in formula (37) with $\epsilon = 10 \u2212 3$, we get $ t mean D=11.849$, $ t mean E=19.868$, and $ t mean F=26.510$, whereas the empirical values are $ t emp D=11.140$, $ t emp E=21.760$, and $ t emp F=26.940$.

This observation could be taken as a hint that a path will synchronize faster than a more integrated or even complete network. Things are actually more elaborate as argued in the following important remark.

So far, we have assumed a constant value equal to $1$ for the coupling coefficient $ c 2 \u2032$ between nodes. In this remark, we want to show that different networks may exhibit different rankings in the speed and time of synchronization according to the values of the coupling coefficient. Let us consider again the triplet of networks D, E, and F in Fig. 7, and let us assume as synchronization time the lower bound in formula (37) with a threshold equal to $\epsilon = 10 \u2212 6$.

As $\alpha $ decreases, more eigenvalues become complex, and if $\alpha < 2 \mu 1$, they are all complex. Therefore, what matters is the real part, which is just $\u2212 \alpha \mu i 2$. The maximum non-null value is then for all $\u2212 \alpha \mu n \u2212 1 2$. For the path graph, $ \mu n \u2212 1=2 [ 1 \u2212 cos \u2061 ( \pi / n ) ]$, and for the complete graph, $ \mu n \u2212 1=n$. Then, we have to compare $\u2212\alpha ( 1 \u2212 cos \u2061 ( \pi / n ) )$ and $\u2212 \alpha n 2$. Since $n>2 ( 1 \u2212 cos \u2061 ( \pi / n ) )$ for every $n>2$, we have that $\u2212 \alpha n 2<\u2212\alpha ( 1 \u2212 cos \u2061 ( \pi / n ) )$, which shows that, for $\alpha < 2 \mu 1$, the high connectivity of the complete network enables it to synchronize faster than other topologies. Similar results for networks A–C in Fig. 7 confirm these conclusions.

## V. RESONANT STATES IN COMPLEX NETWORKS

We move now to the third case in Table I by introducing a driving periodic force acting on or even produced by a single node. We will analyze the effect of a standard periodic sinusoidal force with frequency $\omega $. Our purpose is to illustrate some new results about resonant phenomena on networks. In order to introduce method and notations, we quickly revise some ideas about the single driven resonant oscillator.

### A. The driven harmonic oscillator

Let us observe that the limit for $\omega \u2192 \omega 0$ in Eq. (52) shows different behavior than in Eq. (48): when $\omega $ approaches $ \omega 0$, only the term $sin\u2061 \omega 0 \omega \pi $ in the denominator approaches $0$, while other terms represent an oscillation of frequency $ \omega 0$. Therefore, the resonance produced by this kind of force preserves the proper frequency $ \omega 0$, while the amplitude goes to infinity uniformly for all $t$.

### B. Network case

Let us assume a driving force $f(t)= F 0sin\u2061(\omega t)$ applied on node $h$, $h=1,\u2026,n$. The general solution and the corresponding resonant states are described by the next proposition.

*proper*frequencies of the network. The diverging term in the position components in Eq. (57) is

The first conclusion is that the amplitudes of both $ x k(t)$ and $ v k(t)$, *in general*, grow linearly with $t$ as $\omega $ approaches one of the proper frequencies of the network. However, things are slightly more involved than this, as discussed in the next remark.

Let us imagine that, for $i= \u0131 ~$, the product $ \varphi \u0131 ~(h) \varphi \u0131 ~(k)$ is equal to $0$.

Let us begin by examining the case in which $ \varphi \u0131 ~(k)=0$. Then, even if $\omega \u2192 \omega \u0131 ~$, the node $k$ cannot resonate with the source node $h$ at such a frequency $ \omega \u0131 ~$. In other terms, if we make the source node $h$ oscillate with frequency $\omega = \omega \u0131 ~$ and $ \varphi \u0131 ~(k)=0$, then node $k$ is not affected by the oscillation in $h$. It is as if node $k$ is transparent to the propagation of the oscillation of frequency $ \omega \u0131 ~$ from the source node $h$. The network as a whole cannot transmit such a frequency to node $k$, and this node does not participate to the global resonant process.

On the other hand, if it is $ \varphi \u0131 ~(h)=0$, then an oscillation of frequency $\omega = \omega \u0131 ~$ in the source node $h$ cannot be spread throughout the network to any node, and no resonance phenomenon can be induced from such a node with that frequency.

The only frequency that is able to always put in resonance the entire network is $ \omega \u0131 ~= \omega n=1$. In fact, if we look at Eq. (59), being $ \varphi n(h)\u22600$ and $ \varphi n(k)\u22600$, $\u2200h,k$, when $\omega \u21921$, any receiver node $k$ resonates with any source node $h$ oscillating at such a ground frequency.

The four possible configurations of the phases of nodes, for the different resonant driving frequencies applied in node $1$, are depicted in Fig. 10.

In particular, if the driving force in node number $1$ has frequency $ \omega 2=2$, node numbers $3$ and $4$ do not resonate since $ \varphi 2(3)= \varphi 2(4)=0$; similarly, when node number $1$ oscillates with frequency $ \omega 3= 2$, only node number $3$ is prevented from entering the resonant state.

We also said that if the driving force is placed in a node $h$ such that $ \varphi i(h)=0$, then if $\omega \u2192 \omega i$, we have $ x k(t)=0,\u2200k=1,\u2026,n$. For instance, if the force is placed in node $3$, being $ \varphi 2(3)=0$, this implies that for $\omega \u2192 \omega 2=2$, all the nodes stay at rest. The network is not able to transmit such a frequency from such a node. Other frequencies may have an effect, but there will be no resonance with frequency $ \omega 2=2$ placed in node number $3$.

In Fig. 11, we illustrate the position plots of the four nodes in the toy network when an external force $f(t)=sin\u2061(\omega t)$ is applied to node number $1$ starting from initial conditions $ x 0= [ 0 , 0 , 0 , 0 ] T$ and $ v 0= [ 0 , 0 , 0 , 0 ] T$. As $\omega $ approaches $ \omega 4=1$ from below, all the nodes are activated and enter a global resonant state, as shown in Fig. 9 panel (d). In Fig. 12, we illustrate the position plots of the four nodes under the same conditions as before but as $\omega $ approaches $ \omega 3= 2$ from below; let us notice that all nodes are activated, except node number $3$ (green line in the plot), which does not respond to node $1$ stimuli at this frequency, as shown in Fig. 9 panel (c).

### C. Synchronization with resonance

where $B= L 2 \u2212 4 I$.^{36} Let us now suppose again a driving force acting on node $h$ of the form $f(t)= F 0sin\u2061\omega t$. The general solution and the corresponding resonant states are described by the next proposition.

Notice that, with respect to Eq. (26), we dropped the index G in the eigenvalues to make the expressions lighter.

The question now is what frequency or frequencies are able to induce resonance phenomena in the network as $t$ grows to infinity, bearing in mind the fact that the presence of dissipative terms still drives synchronization between nodes. This amounts to asking what are the differences between formula (59) and formula (68).

First, let us note that the ratio before the square bracket in formula (68) is always a real value. This is due to the fact that $ \lambda i \u2212$ and $ \lambda i +$ are reciprocal and to the equality $ ( \lambda i \u2212 ) 2+ ( \lambda i + ) 2= \mu i 2\u22122$. Moreover, as previously discussed, the exponential terms $ e \lambda i \xb1 t$ inside the brackets go to zero, as $t$ grows. Therefore, in general, for $t\u2192+\u221e$, we have coefficients in the sum depending on $k$, which prevents nodes from synchronizing.

However, when $\omega \u21921$, the last term in the sum becomes dominant, and this term is common to all the nodes. Indeed, the only frequency $\omega $ that can prompt resonance phenomena is the one such that $ \omega 2+ ( \lambda i \u2212 ) 2=0$ and $ \omega 2+ ( \lambda i + ) 2=0$ and the only value that can make these quantities vanishing is the pure imaginary value $ \lambda n += \lambda n \u2212=i$.

In Fig. 13, we illustrate the position plots of the four nodes in the toy network when an external force $f(t)=sin\u2061\omega t$ is applied to node number $1$ starting from initial conditions $ x 0= [ 0 , 0 , 0 , 0 ] T$ and $ v 0= [ 0 , 0 , 0 , 0 ] T$. As $\omega $ approaches $1$ from below, all the nodes are activated and enter a global resonant state, as proven above. In Fig. 14, we illustrate the position plots of the four nodes under the same conditions as before, but as $\omega $ approaches, for instance, $ 2$ from below; let us notice that, as expected, no resonant state arises within the network, in contrast with what shown in Fig. 12 under the same conditions.

## VI. THE LINEAR SWING EQUATION IN NETWORKS

^{18,20}The swing equation plays a central role in the analysis of power systems dynamics and in a variety of different contexts, such as the description of superconducting Josephson junctions or a nonlinear second-order phase-locking loop in electronics. It governs the rotational dynamics of synchronous machines in stability studies and it is nonlinear in nature, and a closed formal solution cannot be explicitly given. Equation (71), in particular, models the dynamics of power grid networks, where nodes can be either generators (generating power) or loads (consuming power), and edges are transmission lines connecting them. The requirement for a synchronous generator is that it does not lose phase synchronization when small temporary changes occur in the output of the machine. Typical examples are small perturbations in the scheduled generation of the machine, which result in a change in its rotor angle, or small loads added to the network. The response of a power system to impacts is oscillatory, and if the oscillations are damped, the system can reach in time a steady state. Transient stability studies determine whether or not synchronization is maintained after the machine has been subjected to a transient disturbance. Each variable $ x i$ in Eq. (71) represents the corresponding node displacement and $ v i= x \u02d9 i$ its nodal velocity relative to a rotational reference value; $\gamma $, which we assume common to all the nodes, is the damping coefficient and the vector $p$ has components equal (or, more, in general, proportional) to the small power deviation generated (if positive) or consumed (if negative) by nodes. The general solution of Eq. (69) can be written as

We now limit the discussion to balanced power grids; that is, we assume that $ \u2211 i = 1 n p i=0$. Furthermore, we suppose that the power grid is connected so that the multiplicity of the null eigenvalue of $L$ is equal to $1$. Under these conditions, the system $p=L x ~$ is consistent since $ rank(L |p)= rank(L)$. The family of solutions is characterized by the fact that if $ x ~ 1$ and $ x ~ 2$ are two of them, then $ x ~ 1\u2212 x ~ 2=wu$, where $w\u2208 R$ and $u\u2208 R n$ is the vector of all $1$’s. Since $G [ u , 0 ] T= [ 0 , 0 ] T$ and $ e G t [ u , 0 ] T= [ u , 0 ] T$, then the general solution $y(t)= e G t y 0+ y ~\u2212 e G t y ~$ is the same for any $w\u2208 R$. This means that we can use any solution of the linear system $p=L x ~$ to build the general solution in Eq. (72).

Let us illustrate these ideas on the toy network in Fig. 1, postponing the study of a real power grid to Sec. VII. Let us suppose that node $3$ is a generator and nodes $1$, $2$, and $4$ are loads and that the corresponding powers are $p= [ \u2212 0.50 , \u2212 0.20 , 1.05 , \u2212 0.35 ] T$. Then, a solution of $p=L x ~$ is, for instance, $ x ~= [ \u2212 0.25 , \u2212 0.15 , 0.15 , \u2212 0.20 ] T$. By substituting $ y ~= [ x ~ , 0 ] T$ and $ y 0= [ 0 , 0 ] T$ in Eq. (71), we get the results depicted in Fig. 15 for two different values of the damping parameter $\gamma $.

## VII. NUMERICAL EXPERIMENTS

### A. Application to a social network

Let us imagine that in a social network, there is a debate around a certain topic and different members have different opinions on that topic. Take, for example, the case in which there are different opinions about the qualities of a particular sportsman or of a given member within the same network.

Each member can shape a positive or negative assessment with all possible gradations from full agreement to full disagreement, and he or she can share this opinion with his or her neighbors. We can model this static situation at a fixed time, for instance, by assigning $+1$ to perfect agreement and $\u22121$ to total disagreement and providing an intermediate graduation of the level of disagreement or agreement continuously from $\u22121$ to $+1$. The value $0$ would represent the totally neutral position.

However, each member is allowed to change his or her opinion over time, adjusting the assessment made at a given instant. Our goal here is to figure out what entrainment action some nodes may play with respect to others in this opinion formation process.

Following, for instance, Cartwright and Harary^{37} and Altafini,^{38} it is reasonable to assume that if a member has a close friend, he is much more likely to be influenced by the latter’s opinion than by others’ and therefore, he will move in accordance with his close neighbors’ opinion dynamics. In other words, he will be induced to synchronize the direction and the intensity of the velocity at which he changes his opinion according to that of his friends.

The mapping between the previous theoretical model and the social problem under examination is based on the interpretation of the continuous variables $ x i$ as a measure of the level of agreement or disagreement of each member with the subject. The opinion entrainment effect induced by a node on its neighbors—and, thus, on the network—is then modeled by the couplings described above. What we intend to focus on here is not the numerical value of each individual variable at a given instant. Rather, we are interested in the collective behavior of the set of variables $ x i$ over time. Such behavior makes clear, as will be seen, the effect of the topological structure of the network in shaping processes of synchronization of members’ ideas and opinions.

A celebrated example of a social network in network theory is provided by the Zachary Club network, see Fig. 16, in which $34$ members found themselves having to decide whether to agree with the position of the instructor ( $+1$) or with the position of the president ( $\u22121$). We can start by assigning an initial “position” vector in such a way that node number $1$ in red, the instructor, has position $+1$ and node number $34$ in blue, the president, has position $\u22121$. At the beginning, at time $t=0$, all the other nodes are supposed to be totally neutral since they still have to build an opinion, and therefore, they are assigned a “position” equal to $0$.

We start by assuming that the social system evolves according to Eq. (28) with a matrix $G$ as in Eq. (24). The plot in Fig. 17, panel (a), illustrates the position vs time of two different groups of nodes with initial conditions equal to $+1$ for node number $1$ (the instructor), $\u22121$ for node number $34$ (the president), and $0$ for all the other nodes. Specifically, we choose four nodes for each one of the two communities illustrated in Fig. 16 and colored in red and blue, respectively, for the group of the instructor and the group of the president. The four nodes in red are node numbers $4,11,13$, and 18, and the four nodes in blue are node numbers $23,24,30$, and 33. As can be seen, nodes in the two different groups tend to synchronize each other from the very beginning. As time goes on, both converge to zero due to the balanced initial conditions that see neither leader prevail over the other. What is interesting and unexpected is the fact that each group seems to assume a position, which is opposite to that of their respective leader. However, it could be easily explained by saying that if node $1$ is initially in position $+1$, then it can move only in the negative verse toward $0$, dragging in the negative verse also nodes that are under its main influence and similarly for node $34$.

The plot in Fig. 17, panel (b), illustrates the position vs time of two groups, including the two respective leaders. In this case, initial conditions see all the nodes having the same initial position equal to $0$ (including nodes $1$ and $34$), but we imagine now that at time $t=0$, node $1$ starts to move with velocity $+1$ (in the positive verse) and node $34$ starts to move with velocity $\u22121$ (in the negative verse). In the plot, nodes in red are numbers $1,4,11$, and 18 and nodes in blue are numbers $34,30,32$, and 33; dashed lines represent the two leader. As can be seen, nodes in the two different groups tend to synchronize each other, to synchronize with their leader and to stay out-of-phase with members of the other group for all the time, until both converge to $0$ due to the balanced initial conditions.

Two synchronization processes with unbalanced initial conditions are illustrated in Fig. 18 for the same set of nodes, as in Fig. 17: in panel (a), node number $1$ and number $34$ have initial positions equal to $0.8$ and $\u22120.5$, respectively; in panel (b), node number $1$ and number $34$ have initial velocities equal to $0.8$ and $\u22120.5$, respectively.

We compare now the synchronization times given by Eq. (37) under different conditions. In detail, we are interested in understanding which of the two leaders is more effective in terms of speed of synchronization. Keeping in mind the discussion in the last remark in Subsec. IV C, we start with strong coupling between nodes $ c 2 \u2032=1$. If only the instructor (node number $1$) is endowed with an initial velocity (in our simulation equal to $4$) and all other nodes (including number $34$) are initially held at zero, then the time at which all nodes synchronize with the asymptotic solution in Eq. (31) below a threshold $\epsilon =0.001$ is equal to $t=91.70$ with a mean time equal to $t=5.29$. If we now assign the same initial velocity to the president (node $34$) while keeping all others (including node number $1$) at zero, then the same time becomes $t=96.20$ with a mean time equal to $t=5.47$. The two times are basically comparable with slightly higher effectiveness for the instructor than for the president, confirming the fact that the network is divided into two balanced communities. The same conclusion is validated using other initial conditions on both the velocity and the position of the two leaders.

Let us now move to the resonant phenomena. We assume now that the social system evolves according to Eq. (56) with a matrix $G$ as in Eq. (54). The resonant frequencies given in Proposition 3 range from $1$ to $4.37$, and they are listed in Table IV.

i
. | μ_{i}
. | ω_{i}
. | i
. | μ_{i}
. | ω_{i}
. | i
. | μ_{i}
. | ω_{i}
. |
---|---|---|---|---|---|---|---|---|

1 | 18.137 | $4.37455$ | 11 | 4.581 | $2.36237$ | 21 | 2.000 | $1.73205$ |

2 | 17.055 | $4.24914$ | 12 | 4.480 | $2.34094$ | 22 | 1.955 | $1.71903$ |

3 | 13.306 | $3.78234$ | 13 | 4.276 | $2.29693$ | 23 | 1.826 | $1.68109$ |

4 | 10.921 | $3.45269$ | 14 | 3.472 | $2.11476$ | 24 | 1.762 | $1.66190$ |

5 | 9.777 | $3.28287$ | 15 | 3.382 | $2.09332$ | 25 | 1.599 | $1.61223$ |

6 | 6.996 | $2.82776$ | 16 | 3.376 | $2.09193$ | 26 | 1.259 | $1.50313$ |

7 | 6.516 | $2.74145$ | 17 | 3.242 | $2.05963$ | 27 | 1.125 | $1.45774$ |

8 | 6.332 | $2.70769$ | 18 | 3.014 | $2.00349$ | 28 | 0.909 | $1.38176$ |

9 | 5.618 | $2.57255$ | 19 | 2.749 | $1.93627$ | 29 | 0.468 | $1.21183$ |

10 | 5.379 | $2.52559$ | 20 | 2.487 | $1.86738$ | 30 | 0 | $1.00000$ |

i
. | μ_{i}
. | ω_{i}
. | i
. | μ_{i}
. | ω_{i}
. | i
. | μ_{i}
. | ω_{i}
. |
---|---|---|---|---|---|---|---|---|

1 | 18.137 | $4.37455$ | 11 | 4.581 | $2.36237$ | 21 | 2.000 | $1.73205$ |

2 | 17.055 | $4.24914$ | 12 | 4.480 | $2.34094$ | 22 | 1.955 | $1.71903$ |

3 | 13.306 | $3.78234$ | 13 | 4.276 | $2.29693$ | 23 | 1.826 | $1.68109$ |

4 | 10.921 | $3.45269$ | 14 | 3.472 | $2.11476$ | 24 | 1.762 | $1.66190$ |

5 | 9.777 | $3.28287$ | 15 | 3.382 | $2.09332$ | 25 | 1.599 | $1.61223$ |

6 | 6.996 | $2.82776$ | 16 | 3.376 | $2.09193$ | 26 | 1.259 | $1.50313$ |

7 | 6.516 | $2.74145$ | 17 | 3.242 | $2.05963$ | 27 | 1.125 | $1.45774$ |

8 | 6.332 | $2.70769$ | 18 | 3.014 | $2.00349$ | 28 | 0.909 | $1.38176$ |

9 | 5.618 | $2.57255$ | 19 | 2.749 | $1.93627$ | 29 | 0.468 | $1.21183$ |

10 | 5.379 | $2.52559$ | 20 | 2.487 | $1.86738$ | 30 | 0 | $1.00000$ |

Let us imagine now that an external force with frequency $\omega $ is applied to node number $1$ or, in other words, that the first leader acts as an influencer by applying the same driving force in that point of the network. Let us consider two nodes belonging to the two different clusters previously considered: node number $11$ in the first cluster (red colored) and node number $30$ in the second cluster (blue colored). Let us consider their behavior as $\omega $ approaches two resonant frequencies, specifically $ \omega 30=1$ and $ \omega 22=1.719026$. Results are illustrated in Figs. 19 and 20.

As discussed in the remark after Proposition 3, all nodes resonate at the ground frequency $\omega =1$. Node number $11$, which belongs to the sphere of influence of the first leader, certainly resonates at the ground frequency and so does at the second frequency under consideration; on the contrary, node number $30$, which falls within the sphere of influence of the second leader, while behaving similarly to node $11$ at the ground frequency, does not respond to the second frequency from node $1$. This is easily explained if we look at the component of the corresponding eigenvector $ \varphi 22$. Indeed, for node number $11$, it is $ \varphi 22(11)=\u22120.540899038$, whereas for node number 30, it is $ \varphi 22(30)=0.002692861$. The first component is about 200 times larger than the second one, which justifies the very low influence of node 1 on node 30 at that frequency.

Let us consider now the frequency $ \omega 9=2.572554$. Since $ \varphi 9(1)=0$, we expect that the oscillation of node $1$ with that frequency does not affect the network. This is confirmed by Fig. 21, where we selected two arbitrary nodes (numbers 13 and 29) in the two different clusters, which show no amplification of their own oscillation when the frequency of node 1 approaches the indicated value.

Further confirmation comes from the frequency $ \omega 6=2.827755$. Since $ \varphi 6(1)=\u22120.0027808652$ and $ \varphi 6(4)=\u22120.8231724260$, we expect that the network globally responds much more if that oscillation starts from node $4$ than from node $1$. Figure 22 shows the different behavior of two randomly selected nodes (again numbers 13 and 29) when such a frequency is applied on the network starting from different source nodes, specifically nodes 1, 4, and 34. The amplitude of the oscillation of nodes $13$ and $29$ grows much more quickly when the source is placed in node $4$ than when it is placed in node $34$ or, worse, in node $1$. The two leaders prove to be not so effective at this specific frequency.

Finally, let us consider the case of the resonant synchronization described by proposition 4 by assuming now that the social system evolves according to Eq. (68) with a matrix $G$ as in Eq. (66). Figure 23 shows the behavior of the two leaders when the first one approaches frequency $\omega =1$ from below. All the other nodes display similar behavior and enter a synchronized resonant state with the two leaders. No other frequency is able to produce a state of synchronized resonance on all nodes in the network.

### B. Application to a power grid

^{39}A power grid consists of a network whose nodes are rotating machines and edges are transmission lines. Nodes can produce power (generating plants) or consume power (loads). The network topology of the Syrian grid is depicted in Fig. 24, panel (b). It consists of $19$ nodes with $6$ generators (nodes 1–6, in red), $13$ loads (nodes 7–19, in blue), and $24$ edges.

We aim primarily to show the qualitative behavior of the solution described in Eq. (72) for specific nodes and for different values of the involved parameters. In Fig. 25, we represent the dynamic response of the selected nodes in the Syrian power grid: nodes $3$, $4$, and $5$ correspond to the main power plants in terms of expected capacity, namely, Baniyas, Zayzoun, and Al-Thawrah plants, whereas nodes $9$, $15$, and $18$ are located in three major Syrian cities, namely, Damascus, Palmyra, and Al-Hasakah. In panel (a), we consider a distribution of nominal balanced powers such that $ p i=0.13$ for generators and $ p i=\u22120.06$ for loads with damping coefficient $\gamma =1$; in panel (b), we keep the same powers but with damping coefficient $\gamma =0.4$; and in panel (c), we assume $ p = [ 0.004 , 0.181 , 0.252 , 0.202 , 0.306 , 0.056 ]$ for the six generators and $ p i=\u22120.077$ for all the motors. Figure 25 qualitatively illustrates some important features of the transient dynamics of $ x i(t)$, such as the first peak value, the first peak time, and the value at a steady state under different conditions (see Bhatta *et al.*^{18} for a comprehensive discussion about these values). Table V collects the values of these parameters for two representative nodes in the examined simulation and under conditions corresponding to panels in Fig. 25.

For instance, as can be observed, the plant corresponding to node $4$ typically has the highest peak value, while node $3$ has the shortest peak time among the examined nodes. For the sake of brevity, we do not consider a power distribution that respects the symmetries of the network, which is such that all nodes within the same symmetry cluster have the same power. In this case, it has been shown that if the vector