Delay Regulated Explosive Synchronization in Multiplex Networks

It is known that explosive synchronization (ES) in an isolated network of Kuramoto oscillators with inertia is significantly enhanced by the presence of time delay. Here we show that time delay in one layer of the multiplex network governs the transition to synchronization and ES in the other layers. We found that a single layer with time-delayed intra-layer coupling, depending on the values of time delay, experiences a different type of transition to synchronization, e.g., ES or continuous, and the same type of transition is incorporated simultaneously in other layer(s) as well. Hence, a suitable choice of time-delay in only one layer can lead to a desired (either ES or continuous) transition simultaneously in the delayed and other undelayed layers of multiplexed network. These results offer a platform for a better understanding of the dynamics of those complex systems which are represented by multilayered framework and contain time delays in the communication processes.


Introduction
The networks have been proven to be an important tool to investigate manifold processes such as synchronization, epidemic spreading, transport and communication in complex systems [1,2]. In recent years, a novel synchronization process exhibiting hypersensitivity or explosiveness in the natural physical and biological systems, called explosive synchronization (ES), has drawn much attention in the scientific community [3,4,5,6,7,21,9,10]. The epileptic seizures in the brain [11], chronic pain in the Fibromyalgia brain [12], the cascading failure of power grids [13], and the jamming of the Internet [14] have shown ES transitions with small initial perturbations. Further, it is known that time-delay is naturally present in communication systems and enhances For un-delayed case, all the three layers manifest a regular second-order transition, however, as a time delay is employed in a single layer (layer 1), all the three layers experience ES transition simultaneously for certain value of delay. λ depicts over-all coupling strength, and R 1 , R 2 and R 3 provide means of global synchronization as defined in section 2.
synchronization in the corresponding systems [15,16,17,18,19,20]. The effect of time delay on ES transition has been demonstrated in the Kuramoto model having the degree-frequency correlation, and in the second-order Kuramoto model [21,22].
Lately, extensive studies of ES transition on single layer networks have been carried out. However, in many complex systems several dynamics take place at the same time and to understand these systems in better manner an isolated network structure turns out to be insufficient and hence, a multilayered structural approach is required. In the multilayered approach, same set of the nodes, representing a complex system's units, are laid out in several discrete layers, each one of them representing a different type of interaction process existing among the nodes. In the multilayered networks, the synchronization of different dynamical processes evolving simultaneously in individual layer may get affected by structural and dynamical properties of other layers leading to different characteristics for the synchronization experienced by the layer. Relay, cluster and delay induced synchronization have been reported in multiplex architectures [23,24,25]. Of late, the emergence of ES has been shown in the multilayer networks comprised of Kuramoto oscillators as well as second-order Kuramoto oscillators [26,27,28].
In this paper, we discuss about the repercussions of inclusion of time delay on explosive synchronization in the multilayer networks. We represent the nodes in multiplex networks by second-order Kuramoto oscillators [29,30,31,32,33]    is introduced in one layer and its impact on synchronization in other layers are studied. It is perceived that inclusion of time delay in only one layer induces either ES or non-ES transition concurrently in all the layers depending on the value of the time-delay. It is further revealed that the occurrence of ES or non-ES transition is solely governed by the evolution of average frequency of the system. Such occurrences of ES or non-ES transition were first revealed for the single layer time delayed networks [21,22]. The important revelation coming out from our study is that by making an appropriate choice of time delay in one layer, a desired type of transition can be achieved concurrently in all the multiplexed layers. The impact of various structural and dynamical properties on the observed transition behavior has been discussed in details both numerically and analytically.

Delayed Multiplex Network Model
In the current study, we investigate how the presence of time-delay between the pair of nodes in a single layer affects transitions to synchronization in other layers of multiplex networks. Here, we consider undirected multiplex networks where nodes in each layer are represented by the second order Kuramoto oscillators with one layer subject to time delayed coupling while rest to regular coupling between each pair of nodes (figure 1). We examine the impact of time-delayed layer on transition to synchronization in other undelayed layers. The multiplex network we consider is comprised of L layers of networks having same number of the nodes N . The time-evolution dynamics of the nodes, denoted by second order Kuramoto oscillators, in the presence or absence of the time-delay in a layer are governed by respectively, where i = 1, ..., N , Ω i (0) (θ i (t)) is the natural frequency (the instantaneous phase) of ith oscillator, parameter m is mass, λ is coupling strength and τ is time-delay introduced between the phases of nodes in a layer. Superscript α denotes the single delayed layer while β = α denotes all the other L − 1 un-delayed layers. σ αβ = σ βα denotes inter-layer coupling strength between layer α and β. The multiplex network is represented by L layers encoded by a set of adjacency matrices, {A 1 , ..., A κ , ..., A L }, where element A ij κ = 1 when the nodes i and j are connected and A ij κ = 0, otherwise. Adjacency matrix (size N L × N L) of the multiplex network can be written in block matrix form as following where I and O are identity and null matrices, respectively. We quantify transition to synchronization in a layer by means of the order parameter. We define order parameters R α and R β at time t for the time-delayed and regularly coupled layer as where ψ(t) denotes average phase of the nodes. Hence, R is a measure of coherence of the collective dynamics of the nodes, i.e., the degree of synchronization of the network. R = 1 corresponds to a completely synchronized state, while R = 0 denotes a asynchronous state.

Simulation results
Here, we present and discuss in details various simulations results about the effect of time delay in a single layer on synchronous transitions in the other layers of multiplex networks. First, we compute the order parameters defined in (3)  frequencies, etc., each layer of the multiplex network is subject to the same coupling strength λ = λ α = λ β = λ γ , · · ·. We adiabatically increase the coupling strength λ from λ 0 (incoherent) to λ 0 +nδλ 0 (synchronous) state in the step of δλ 0 . We also adiabatically decrease the coupling strength from λ 0 +nδλ 0 (synchronous) to λ 0 (incoherent) in the step of δλ 0 , and compute the order parameter R for each layer at each value of λ. Before each δλ 0 step, we eliminate initial transients and integrate the system long enough (10 4 time steps) using a fourth-order Runge-Kutta method with time step dt = 0.01, to arrive at the stationary state. For all the simulations, initial phases for the nodes in the individual layer are drawn from a random uniform distribution in the range [0, 2π). We have made a specific choice for the natural frequencies, the natural frequencies of the nodes in each layer are set equal to their respective degrees, i.e., Ω i (0) = k i /k max , where k i = N j=1 A ij is the degree of a node and k max is maximum degree in the network. All the simulation results presented in this paper are for N = 100 nodes in each layer constructed with average degree k = 12 for various network topologies, otherwise mentioned elsewhere.
For numerical study, firstly, we consider a duplex network comprised of two ER random networks [34] of average degrees k 1 = k 2 = 12 and introduce a constant time delay between the phases of nodes in layer 1. We compute the order parameters for both the layers in figure 2. In the absence of time-delay, i.e., τ = 0, both the ER layers display second-order transitions and synchronize simultaneously. However, the presence of time delay significantly enhances the synchronization behavior of the multiplexed layers. For different values of the time delay, different transition behaviors, i.e., either strong or weak ES or second-order is witnessed in both the layers. For an instance in figure 2, rather a weak ES is seen for τ = 2.5 and 17, strong ES for τ = 4 and τ = 15, and very weak ES or second-order transition for τ = 0.5, τ = 7 and τ = 11.5, respectively. Moreover, both the layers synchronize simultaneously exhibiting the same  nature of transition. The onset of transition to synchronization for different values of τ is witnessed at different critical coupling strength. From these observations one can infer that introduction of a time-delay not only enhances transition in the layer that includes the time delay, but also the same nature of transition is incorporated in the other layer(s) simultaneously via inter-layer coupling. Hence in the multiplex networks, a desired (ES or continuous) transition can be incorporated simultaneously in all the layers by making a suitable choice of the time-delay.
A dependency of average frequency Ω on the time delay τ , i.e., Ω = Ω(τ ) is responsible for such interesting behavior. Therefore, to understand the underlying microscopic dynamics behind such time-delay dependent behavior observed in a multiplex network, we define effective frequency Ω i of each node and average frequency Ω of the network as where T is the total time of the averaging after eliminating the initial transients. We compute average frequencies Ω 1 and Ω 2 for the two layers as a function of time-delay τ (see figure 3). The initial (at t=0) average frequencies of the first and the second layers are Ω 0 1 = 0.575 and Ω 0 2 = 0.573, respectively. In the presence of time-delay, the average frequencies of both the layers Ω 1 and Ω 2 , simultaneously exhibit an oscillatory behavior around Ω 0 1 and Ω 0 2 ( figure 3). Starting from 0, as the value of τ is increased, Ω 1 and Ω 2 exhibit a gradual decrease. Interestingly, at a certain value of τ , Ω 1 and Ω 2 suddenly jump to higher values and again start decreasing along a new branch as τ is increased further.
From figure 3, it is apparent that the value of τ for which differences |Ω 1 − Ω 0 1 | and   are low, a second-order transition is observed in both the layers (see figure 2 for τ = 0.5 and τ = 11.5). Hence, such oscillatory dependence of Ω 1 and Ω 2 on the time delay τ causes different type of transitions in the multiplexed layers.
To further understand the microscopic dynamics of ES with time delay, we plot the order parameter R 1 along with the corresponding effective frequencies of nodes Ω i 1 and steady-state phases θ i 1 associated with the layer 1 in figure 4 as both the layers synchronize simultaneously displaying the same transitional behavior for a certain value of the time delay. It is apparent that for τ = 4, ES transition takes place as effective frequency Ω i 1 and, in turn, phases θ i 1 get locked to their respective average values together at the critical coupling strength (see figure 4). For τ = 2.5, most of the nodes are locked to their respective average frequencies and average phases at the critical coupling strength, however, rest of the nodes get locked later at a bit higher coupling strength. For τ = 0.5 (see figure 4), Ω i 1 and θ i 1 gradually start reaching to their respective stationary values near the critical coupling strength, thereby leading to the continuous transitions.

Effect of average degree of the multiplexed layers on nature of transition:
Next, we study the effect of different average degree of the network layers on synchronization transition. To investigate this, we fixed the average degree of the ER layer 1 to a value k 1 and chose an equal, a lower and a higher value for the average degree k 2 of the ER layer 2 (figure 5). For τ = 4, ES and second-order transitions corresponding to k 2 = 6 and 20, respectively, are observed at different 0.0 0.5 critical coupling strengths as compared to that of k 2 = 12. In the same fashion, weak ES transitions for τ = 11.5 corresponding to k 2 = 6 and 20 are observed at different critical coupling strengths as compared to second-order transition observed for k 2 = 12. Therefore, the variation in the average degree k 2 of the non-delayed layer shifts the critical coupling strength which, in turn, enhances or subsides the transition accordingly. Different oscillatory dependence of Ω 1 and Ω 2 on τ , corresponding to k 2 = 6 and 20 (see figure 6), accounts for this behavior.

Effect of number of multiplexed layers L on the transition:
It is also important to verify whether the same set of observations hold for the multiplex networks having more than two layers. In order to investigate this, we consider a multiplex network consisting of three ER layers and compute the order parameter for each layer for different values of the time delay. From figure 7, it is quite apparent that for each value of the time delay τ , a different path to the synchronization, i.e., ES or continuous transition takes place. It is to be noted that even though M follows the connection matrix (2), which indicates that the delayed layer 1 and the un-delayed layer 3 are not directly interconnected, yet, by means of multiplexing via layer 2, the layer 3 exhibits exactly the same transition type as that of the layer 1. Hence, this is an instance of distant coordination similar to the relay synchronization in the multiplex networks [23]. On the basis of this finding we can say that for a multiplex network of a finite number of layers, the intra-layer delay introduced in one layer is sufficient to incorporate the same transition behavior in rest of the layers. However, the oscillatory dependence of average frequency on delay for all the layers changes as the number of layers L increases, which, in turn, changes the critical coupling strength and corresponding transition type for all the layers concurrently. For an instance, transitions corresponding to τ = 2.5 and τ = 4 are observed at different coupling strengths for the 2-layered (see figure 2) and the 3-layered (see figure 7) multiplex networks.

Synchronization on a variety of multiplex networks:
Next, we investigate the transition behavior in the presence of time delay on a variety of multiplex networks constructed from selecting different network topologies. To explore this, we consider multiplex networks comprised of two small-world (SW) [35], two scalefree (SF) [36], and an ER and a SF networks having time delay introduced in one layer (in ER layer in case of ER-SF multiplexed layers). We follow the same analysis as performed for the case of two ER layers. In the absence of time delay (τ = 0) (see figure 8), the two SW layers, the two SF layers and the ER-SF layers exhibit ES, strong ES, and weak ES transitions, respectively, in contrary to the second-order transition observed in the case of two ER layers. However, depending on the value of time delay either a strong ES, weak ES or second-order transition is observed for each selected network-pair of the multiplexed layers (see figure 8). To fathom the role of average frequency of each layer in determining the transition behavior of the selected pairs of multiplexed layers, we compute the average frequencies Ω 1 and Ω 2 for the layers as a function of the time delay as shown in figure 9. For these cases too, an oscillatory dependence of average frequencies Ω 1 and Ω 2 for the two layers on time delay τ is witnessed. In the case of two SF multiplexed layers, the observed oscillatory branches of average frequencies Ω 1 and Ω 2 are stretched further along the axis τ as compared to those belonging to the other pairs of the multiplexed layers.
Therefore, it is manifested that the time delay plays an important role in determining the nature of transition in a variety of multiplexed networks. We have shown that in multiplexed networks, time delay incorporates same intra-layer transition properties in all of the layers. In the presence of time delay, difference in the average connectivity of the multiplexed layers also seems to either strengthen or subside ES transition in the layers.

Analytical treatment
In section 3, we have witnessed that the departure of average frequency Ω from its initial value Ω 0 determines the nature of transition to synchronization. Additionally, it is numerically demonstrated that the average frequency Ω of a multiplexed layer depends not only on the time delay τ but also on the average degree k of all the layers. In this section we make an attempt to elucidate analytically the dependence of average frequency Ω on dynamical as well as structural properties of the multiplex networks. Therefore, to obtain average frequency, we rewrite (1) in the terms of the order parameters R α (t) and R β (t) using (3) and by defining degree of a node as k In the stationary state, R α (t) R α and R β (t) R β . Next, we consider a separate rotating frame for each layer with average phase ψ α (t) = Ω α t. The new variables in the rotating frame are then defined as φ i α = θ i α − ψ α . Sets of equations (6) get transformed to the following in the rotating frames In synchronous state, it is safe to assume that all the nodes in each layer are locked to their respective mean field ψ α , i.e., ψ α = θ i α . Therefore, summing over the set of N equations for each layer, (7) yields the following relations for the average frequencies , where For the sake of convenience, we present our analysis to a two layered multiplex network, i.e., L = 2. For a duplex network, (8) can be expressed as Further, we define ∆Ω, i.e., the difference between Ω 1 and Ω 2 as where and Ω 1 (τ, λ 1 ) = N λ 1 R 1 k 1 sin[Ω 1 τ ] is also a constant as Ω 1 (τ ) and R 1 would be constant in the stationary state for a certain τ and λ 1 > λ c 1 . In the absence of time delay, i.e., τ = 0, (10) reduces to If Ω 0 1 and Ω 0 2 are equal, that, in turn, implies ∆Ω 0 0, consequently, (11) becomes It is quite apparent from (9) or (10) that inter-dependent Ω 1 and Ω 2 depend on parameters τ , λ, σ 12 , k 1 and k 2 , which is in further concurrence with the numerical results shown in section 3 that the average degree of the layers along with the time delay determine the characteristic of transition in the multilayered framework.

Discussions
In the present work, we have discussed the influence of time delay in a single layer on the nature of phase transition in all the layers of multiplex networks. We demonstrate that path by which multiplexed layers achieve synchronous states, can be controlled by the means of the time delays. The multiplexed layers altogether choose a common path, either a strong or weak ES or a second-order one, to reach to the respective synchronous states, which is solely determined by the value of time delay introduced in a single layer. Each value of the time delay departs the average frequencies of all the layers from their initial values leading to different but nearby values, and the extent of such departure governs the type of transition which is common to all the layers. If the extent of the departure is large, this leads to ES or first-order transition in all the layers. An intermediate departure gives rise to weak ES, whereas a small departure only results in a second-order transition. Our this finding holds true for a variety of multiplex networks constructed from different network topologies. However, in the case of SF-SF multiplex network, no second-order transition is observed, only ES or weak ES transition exists. This might be due to the strong heterogeneity of the SF layers, which favors the ES transition in presence of degree-frequency correlation. In presence of time delays, the difference in average connectivities of the layers brings upon changes in the behavior of oscillatory dependency of the average frequencies on time-delay which, in turn, may result in a different type of transition. It is also demonstrated that for the 3-layered multiplex networks with no direct inter-layer connections between the remote layers, one delayed layer is capable of inducing the same type of transition in the remotely un-delayed layers too, demonstrating that multiplex networks favor remote synchronization. We have also shown analytically that the average frequency of each layer depends on the inter-layer coupling strength and the average connectivity of all the layers along with the time-delay. Lately, the evidence of ES is reported in the brains of patients suffering from fibromyalgia (FM), a condition described as an omnipresent, chronic pain [12]. The patients with FM have brain networks predisposed towards abrupt, global responses or abnormal hypersensitivity (ES) to minor changes. ES is worth investigating as it is detrimental to certain biological and physical circumstances, e.g., epileptic seizure in which things are abruptly turned on [11] and the power grid failure in which things are abruptly turned off [13]. Investigation of ES on real world complex systems by adopting multilayer (or multiplex) architecture has more potential as the multilayer approach provides a more accurate representation of many real world complex systems.
Researchers have employed multiplex perspective to identify the most important functional region of the brain. A human functional brain network consists of peripheral and central or hub regions operating at different frequency bands, which can be modeled as multiplex network in which each layer is identified by a single frequency band carrying unique topological information [37,38]. Information processing from various functional or sensory regions of the brain is a cooperative process of neurons and time delay is naturally incorporated in the information processing to form a coherent and melded perception of the outer real world. Hence, time-delayed synchronization (ES and non-ES) of neurons, while employing multiplex network perspective in the representation of various functional sensory regions, in information processing is worth investigating. Hence, our study on ES on the multiplex framework can be prolific in not only providing better understanding but also regulating the dynamics of interwoven concurrent dynamical processes occurring in both biological and physical systems, specifically the brain network.