Paths to Synchronization on Complex Networks

Paths to Synchronization on Complex Networks Jes´us G´omez-Garde˜nes,1, 2 Yamir Moreno,1 and Alex Arenas3

arXiv:cond-mat/0608314v3 [cond-mat.stat-mech] 31 Jan 2007

1

Institute for Biocomputation and Physics of Complex Systems (BIFI), University of Zaragoza, Zaragoza 50009, Spain 2 Departamento de F´ısica de la Materia Condensada, University of Zaragoza, Zaragoza E-50009, Spain 3 Departament d’Enginyeria Inform`atica i Matem`atiques, Universitat Rovira i Virgili, 43007 Tarragona, Spain (Dated: February 4, 2008) The understanding of emergent collective phenomena in natural and social systems has driven the interest of scientists from different disciplines during decades. Among these phenomena, the synchronization of a set of interacting individuals or units has been intensively studied because of its ubiquity in the natural world. In this paper, we show how for fixed coupling strengths local patterns of synchronization emerge differently in homogeneous and heterogeneous complex networks, driving the process towards a certain global synchronization degree following different paths. The dependence of the dynamics on the coupling strength and on the topology is unveiled. This study provides a new perspective and tools to understand this emerging phenomena. PACS numbers: 05.45.Xt, 89.75.Fb

In 1998 Watts and Strogatz in an effort to understand the synchronization of cricket chirps, which show a high degree of coordination over long distances as though the insects where “invisibly” connected, end up with a seminal paper about the small-world effect [1] that was the seed of the modern theory of complex networks [2, 3, 4]. Many natural and man-made networks have been, since then, successfully described within this framework. Nevertheless, the understanding of the synchronization dynamics in complex networks remains a challenge. The synchronization of non-identical interacting units occupies a privileged position among emergent collective phenomena because of its various applications in Neuroscience, Ecology, Earth Science, among others [5, 6, 7]. One of the most successful attempts to understand it is due to Kuramoto [8, 9], who analyzed a model of phase oscillators coupled through the sine of their phase differences. The Kuramoto model (KM) consists of a population of N coupled phase oscillators where the phase of the i-th unit, denoted by θi (t), evolves in time according to X dθi Λij Aij sin(θj − θi ) = ωi + dt j

i = 1, ..., N

(1)

where ωi stands for its natural frequency, Λij is the coupling strength between units and Aij is the connectivity matrix (Aij = 1 if i is linked to j and 0 otherwise). The original model studied by Kuramoto assumed mean-field interactions with Aij = 1, ∀i 6= j (all-to-all) and Λij = K/N, ∀i, j. The model can be solved in terms of an order parameter r that measures the extent of synchronization in a system of N oscillators as: reiΨ =

N 1 X iθj e N j=1

(2)

where Ψ represents an average phase of the system. The parameter 0 ≤ r ≤ 1 displays a second order phase transition in the coupling strength, being r = 0 the value of the incoherent solution, and r = 1 the value for total synchronization.

The synchronization problem has been solved in some other cases, mainly those where a mean-field approach is also valid [10]. Unfortunately, the mean-field approach requires of several constraints that are not usually fulfilled in real systems. Natural, social and technological systems show intricate patterns of connectivity between their units that are, nowadays, described as complex networks [3, 4]. The problem of synchronization in complex networks inherits the technical difficulties of the non mean-field approaches and incorporates new questions to be considered: What are the new pertinent parameters to deal with synchronization? and, What is the role of the topology in the synchronization process? Several works have partially addressed these issues by studying the stability of the synchronized state [11, 12, 13, 14, 15, 16, 17] using the Master Stability Function (MSF) formalism [18]. However, the onset of synchronization, which posses more theoretical and phenomenological challenges is much less explored, and only a few works have dealt with the study of the whole synchronization dynamics in specific scenarios [19, 20, 21, 22, 23]. The main goal of this Letter is to study the synchronizability of complex networks as a function of the coupling strengths. To do this, first, we propose and discuss a new measure of synchronization for the KM in complex networks. Second, we scrutinize and compare the synchronization patterns in Erd¨os-Renyi (ER) and scale-free (SF) networks and show that even in the incoherent solution, r = 0, the system selforganizes towards synchronization following different paths. Our study reveals that the synchronizability of these networks does depend on the coupling between units, and hence, that general statements about their synchronizability are eventually misleading. Moreover, we show that even in the incoherent solution, r = 0, the system is self-organizing towards synchronization. Let us start by mapping the KM model in finite complex networks as N X dθi Aij sin(θj − θi ) = ωi + λ dt j=1

where λ is a constant.

i = 1, ..., N

(3)

2 1000

1

SF ER

750 GC

r

0.75 0.5

250

0.25

(a)

0 1

0 150

SF ER Nc

0.75

rlink

SF ER

500

0.5

SF ER

100 50

0.25

(b)

0

0 0

0

0.05

0.1

0.15

0.2

λ FIG. 1: (color online) Evolution of a, the KM order parameter defined in Eq. (2), and b the fraction of synchronized links rlink , Eq. (4), as a function of λ. The curves separate when the incoherent solution for SF networks destabilizes. The figure clearly illustrates that the synchronizability of the networks does depend on the value of the coupling strength. Both plots are represented for Erd¨os-Renyi (ER) and scale-free (SF) networks as indicated. The size of the networks is N = 1, 000 and their average degree is hki = 6. The exponent of the SF network is γ = −3.

We study the dynamics of Eq.(3) in ER and SF networks, preserving the total number of links, Nl and nodes, N for a proper comparison [24]. We concentrate in two aspects: global and local synchronization. First, we follow the evolution of the order parameter r, as λ increases, to capture the global coherence of the synchronization in the networks. Secondly, we propose and follow the same evolution for a new parameter, rlink . This parameter measures the local construction of the synchronization patterns and allows for the exploration of how global synchronization is achieved. We define Z tr +∆t 1 X X 1 i(θi (t)−θj (t)) rlink = e dt , lim ∆t→∞ ∆t tr 2Nl i

j∈Γi

(4) that represents the fraction of all possible links that are synchronized in the network (averaged over a large enough time interval ∆t, after the system relaxes at some large time tr ), being Γi the set of neighbors of node i. We solved Eq.(3) using a 4th order Runge-Kutta method for different values of λ, with a uniform distribution of natural frequencies g(ω) in the interval [−π, π] up to achieving the stationary state. The networks are built following a model [25] that generates a one parameter family of complex networks. This parameter, α ∈ [0, 1], measures the degree of heterogeneity of the final networks. A network of size N is

0.05

0.1

0.15

0.2

λ FIG. 2: (color online) Size of largest synchronized connected component (GC) and number of synchronized connected components (N c), as a function of λ for the different topologies considered. Despite r being vanishing and hence no global synchronization is yet attained, a significant number of clusters show up. This indicates that for any λ > 0 the system self-organizes towards macroscopic synchronization. The network parameters are as in Fig. 1.

generated starting from a fully connected core of m0 nodes and a set U(0) of N − m0 unconnected nodes. At each time step, a new node (not selected before) is chosen from U(0) and linked to m other nodes. Each of the m edges is linked with probability α to a randomly chosen node (avoiding multiple and self-connections) from the whole set of N − 1 remaining nodes and with probability (1 − α) following a linear preferential attachment strategy [26]. Repeating these steps (N − m0 ) times, networks interpolating between the limiting cases of ER (α = 1) and SF (α = 0) topologies are generated [27]. In Fig. 1 we represent the evolution of both order parameters, r and rlink , as a function of the coupling strength λ. The global coherence of the synchronized state, represented by r, shows that the onset of synchronization first occurs for SF networks. A detailed finite size scaling analysis performed for both topologies shows that the critical value of the effective coupling, λc , corresponds in SF networks to λSF = 0.05(1), c and in ER networks to λER = 0.122(2), accordingly with c Fig. 1. If λ is further increased, there is a value at which r for the ER crosses over the SF curve. From this value up in λ, the ER network remains slightly more synchronized than the SF network. The behavior of rlink shows a change in synchronizability between ER and SF and provides additional information. Interestingly, the nonzero values of rlink for λ ≤ λc indicate the existence of some local synchronization patterns even in the regime of global incoherence (r ≈ 0). Right at the onset

3

FIG. 3: (color online) Synchronized clusters for several values of λ for the two different topologies studied (ER and SF). These networks are made up of 100 nodes, in order to have a sizeable picture of the system. The evolution of local synchronization patterns is always agglomerative, however, it follows two different routes. In the ER case, the growth of the GC proceeds by aggregation of small clusters of synchronized nodes, while for the SF network the central core groups the smaller clusters around it.

of synchronization for the SF network, its rlink value deviates from that of the ER. While the synchronization patterns continue to grow for the ER network at the same rate, the formation of locally synchronized structures occurs at a faster rate in the SF network. Finally, when the incoherent solution in the ER network destabilizes, the growing in its synchronization pattern increases drastically up to values of rlink comparable to those obtained in SF networks and even higher. The above results show that statements about synchronizability are dependent on the coupling strength value. Additionally, the previous discussion suggests that synchronization is attained following two different paths that depend on the underlying topology. We have studied the characteristics of the synchronization patterns along the evolution of rlink . Synchronization patterns are formed by pairs of oscillators, physically connected, whose phase difference in the stationary state tends to zero. Note that the contribution into Eq.(4) of every pair of connected be written in terms of a ma oscillators can R tr +∆t i(θ (t)−θ (t)) 1 j trix Dij = Aij lim∆t→∞ ∆t e i dt . This tr matrix is filtered using a threshold T such that the fraction of synchronized pairs equals rlink . In this way, if Dij > T oscillators i and j are considered synchronized and the synchronized patterns are extracted. In Fig. 2 we represent the number of synchronized clusters and the size of the largest one (GC) as a function of λ. The local information extracted from it is unveiling an astonishing and novel feature of the synchronization process that can not be derived from Fig. 1, and that in some sense is counterintuitive. The emergence of clusters of synchronized pairs of oscillators (links) in the networks shows that for values of λ ≤ λSF c , i.e., still in the incoherent solution r = 0, both kind of networks have developed a largest cluster of synchronized pairs of oscillators involving 50% of the nodes of the network, and an equal number of smaller synchronization clusters. From this point on, in the SF network the GC grows and the number of smaller clusters goes down, whereas for

FIG. 4: (color online) The plot shows the correlation between the likelihood that a node belongs to the GC of pairs of synchronized oscillators and its degree k as a function of the coupling strength λ. This probability ,P GC (k), is color-coded as indicated in the right panel. The figure convincingly demonstrates that highly connected nodes are those that recruit poorly connected nodes as the GC grows. Network parameters are those used in Fig. 1

the ER network the growth exploits. These results indicate that although SF networks present more coherence in terms of r and rlink , the microscopic evolution of the synchronization patterns is faster in ER networks, being these networks far more locally synchronizable than the heterogeneous ones. The observed differences in the behavior at a local scale are rooted in the growth of the GC. It turns out that for the ER networks, many different clusters of synchronized pairs of oscillators merge together to form a GC when the effective coupling is increased. The coalescence of many small clusters leads to a giant component of synchronized pairs that is almost the size of the system once the incoherent state destabilizes. This is not anymore the case for SF networks, where oscillators are incorporated to the GC practically one-by-one (forming new pairs) in terms of λ (or rlink ), but starting from a core made up of half the nodes of the network. This picture is confirmed in Fig. 3, where we have represented the evolution of local synchronization patterns in ER and SF networks for several values of λ. The ultimate reason behind these two different routes to complete synchronization is the heterogeneous character of the SF network and the role played by the hubs. In Fig. 4, we have plotted the probability that a node with degree k belongs to the GC as a function of its degree k and the coupling λ for the SF network. This probability is an increasing function of k for every λ, hence the more connected a node is, the more likely it takes part in the cluster of synchronized links. Recently [23], Zhou and Kurths have reported the study of hierarchical organization in complex networks, using the MSF and a mean-field approach in the weak coupling limit. Our results thus substantiate and generalize those about the role of hubs in the synchronization process presented in [23].

4 In summary, we have shown that synchronizability of complex networks is dependent on the effective coupling λ among oscillators. For small values of λ, SF networks outperform ER topologies, but the tendency is reverted for intermediate to large values of the coupling. On the other hand, the detailed analysis of evolution of patterns of synchronization showed that there are two radically different mechanisms to attain synchronization. In the presence of hubs, a giant component of synchronized pairs of oscillators forms and grows by recruiting nodes linked to them. On the contrary, in homogeneous structures, many small clusters first appear and then group together through a sharp merging process. These results are as far reaching as the ones obtained for percolation and epidemic spreading on top of homogeneous or heterogeneous graphs, where the radical differences of the system’s dynamics are rooted in the topology of the underlying networks, demonstrating that the same behavior may hold for nonlinear dynamical systems coupled to complex structures. However, at variance with percolation processes, here the synchronization patterns could be directly related to the growth and evolution of the network. Therefore, we have naturally incorporated a dynamics relevant for the emergence of cooperative behavior, showing that the same organizing principle may drive network evolution, i.e., if synchronization is a relevant issue, natural networks (whether they are homogeneous or heterogeneous in degree) can indeed be efficient by adaptively selecting their coupling strengths. Our study then opens new paths to clarify how synchronization is attained in complex topologies and provides new tools to analyze this ubiquitous phenomenon. We thank J.A. Acebr´on, A. D´ıaz-Guilera, C.J. P´erezVicente and V. Latora for helpful comments. J.G.G. and Y.M. are supported by MEC through a FPU grant and the Ram´on y Cajal Program, respectively. This work has been partially supported by the Spanish DGICYT Projects FIS2004-05073C04-01 and FIS2005-00337.

[1] [2] [3] [4]

Watts, D. J. & Strogatz, S. H., Nature 393, 440-442 (1998). Strogatz, S. H., Nature 410, 268-276 (2001). Newman, M. E. J., SIAM Review 45, 167-256 (2003). Boccaletti, S., Latora, V., Moreno, Y., Chavez, M. & Hwang, D.-U., Phys. Rep. 424, 175-308 (2006). [5] Winfree, A. T. The geometry of biological time (SpringerVerlag, New York, 1990). [6] Strogatz, S. H. Sync: The Emerging Science of Spontaneous Order (New York: Hyperion, 2003).

[7] Manrubia, S. C., Mikhailov, A. S. & Zanette, D. H. Emergence of Dynamical Order. Synchronization Phenomena in Complex Systems (World Scientific, Singapore, 2004). [8] Y. Kuramoto, Lect. Notes in Physics 30, 420-422 (1975). [9] Kuramoto, Y. Chemical oscillations, waves, and turbulence (Springer-Verlag, New York, 1984). [10] Acebron, J. A., Bonilla, L. L., Perez Vicente, C. J., Ritort, F. & Spigler, R., Rev. Mod. Phys. 77, 137-185 (2005). [11] Barahona, M. & Pecora, L. M., Phys. Rev. Lett. 89, 054101 (2002). [12] Nishikawa, T., Motter, A. E, Lai, Y.-C. & Hoppensteadt, F.C., Phys. Rev. Lett. 91, 014101 (2003). [13] Hong, H., Kim, B. J., Choi, M. Y. & Park, Phys. Rev. E 69, 067105 (2004). [14] Motter, A.E., Zhou, C. & Kurths, J., Phys. Rev. E 71, 016116 (2005). [15] M. Chavez, D.-U. Hwang, A. Amann, H. G. E. Hentschel, S. Boccaletti, Phys. Rev. Lett. 94, 218701 (2005) [16] Donetti, L., Hurtado, P. I. & Mu˜noz, M. A. Entangled Networks, Phys. Rev. Lett. 95, 188701 (2005). [17] Zhou, C., Motter, A. E., Kurths, J., Phys. Rev. Lett. 96, 034101 (2006). [18] Pecora, L. M. & Carroll, T. L., Phys. Rev. Lett. 80, 2109-2012 (1998). [19] Moreno,Y. & Pacheco, A. F., Europhys. Lett. 68, 603-609, (2004). [20] Oh, E., Rho, K., Hong, H. & Kahng, B., Phys. Rev. E 72, 047101 (2005). [21] Restrepo, J. G., Ott, E. & Hunt, B. R., Chaos 16, 015107 (2006). [22] Arenas, A., Diaz-Guilera, A. & Perez-Vicente, C. J., Phys. Rev. Lett. 96, 114102, (2006). [23] C. Zhou & J. Kurths, Chaos 16, 015104 (2006). [24] Note that for a proper comparison of the synchronizability of different complex networks, the global and local measures of coherence should be represented according to their respective time scales. Therefore, given two complex networks A and B with kmax = kA and kmax = kB respectively, the comparison between observables must be done for the same effective A B coupling λkA = λkB = λ. [25] G´omez-Garde˜nes, J. & Moreno, Y., Phys. Rev. E 73, 056124 (2006). [26] S. N. Dorogovtsev, J. F. F. Mendes, and A. N. Samukhin, Phys. Rev. Lett. 85, 4633 (2000). [27] The study reported henceforth has also been conducted for intermediate values of α. Starting from α = 0 and increasing its value, the networks generated are successively more homogeneous [25] with a heavy tail whose exponent is equal to (α = 0) or larger than (α > 0) γ = 3. We advance that the results for 0 < α < 1 are consistent with the picture described in the rest of the paper.