Tuning synchronization of integrate-and-fire oscillators through mobility

We analyze the emergence of synchronization in a population of moving integrate-and-fire oscillators. Oscillators, while moving on a plane, interact with their nearest neighbor upon firing time. We discover a non-monotonic dependence of the synchronization time on the velocity of the agents. Moreover, we find that mechanisms that drive synchronization are different for different dynamical regimes. We report the extreme situation where an interplay between the time-scales involved in the dynamical processes completely inhibit the achievement of a coherent state. We also provide estimators for the transitions between the different regimes.

Concerning the case of time-dependent networks, up to now the interest has been focused on the need of finding new tools for characterizing topological properties of networks whose connectivity changes in time, but not as much on the interrelation between topological changes and dynamical evolution of the nodes themselves. This interrelation has been mainly considered for the asymptotic cases in which characteristic times of both dynamical processes are very different. When topological changes are very fast, a good proxy is the fast switching approximation (FSA) [8] whereas in the opposite limit, perturbation around the static case makes sense. In this paper we analyze dynamical effects in time dependent networks but understanding the network as the result of the motion of agents that interact when they are close enough [9][10][11].
When considering emergent properties of systems formed by elementary units with its own dynamics, one of the most important ones is synchronization [12]. Actually, synchronized behaviors appear in nature in groups that improve performance based on collaboration [13] and also in human actions it has shown to provide collective benefits [14,15]. The effect of changing patterns of interaction on synchronization features has been analyzed in different settings, for instance in chemotaxis [16], mobile ad hoc networks [17], wireless sensor networks [18] and the expression of segmentation clock genes [19].
Recently, a general framework of mobile oscillator networks where agents perform random walks in a twodimensional (2D) plane has been proposed [20]. This framework, that reduces to FSA when velocity is high enough, is valid for models whose evolution can be well approximated by linear dynamics. This actually holds for models such as populations of Kuramoto oscillators [21,22], whose evolution, after a short transient time, is very well described by a set of linear equations that can be solved in terms of spectral properties of the Laplacian matrix [23]. Alternative approaches based on Fokker-Planck equations have been also proposed recently [24].
Here, we focus on a dynamical system, a population of integrate and fire oscillators (IFO), where the evolution takes place in two different time scales. One for the slow evolution of the internal state variables (the phase and the orientation) and the other for the instantaneous interaction between the units (pulse coupling). During the last years it has been shown that the interaction structure plays a fundamental role in the dynamics of IFO networks [25]. Usually, IFO have been used to model neural systems but we can also find some examples of applications in other fields, as for example in economy [26].
In the present letter we consider a population of integrate-and-fire oscillators, which interact with their nearest neighbor only while freely displacing on a plane. Such a minimal interaction rule, in contrast to other approaches based on an interaction radius [27], provides a strong non-monotonic dependence of the synchronization time on the velocity of the agents. Three different regimes are identified according to the synchronization properties of the system: a slow regime, a fast switching limit and an anomalous intermediate region between them. Remarkably, in this last region a divergence of the synchronization time is observed. Depending on the studied phenomena, synchronization can be taken either as as a positive or negative feature of the collective behavior of the system. Global brain synchronization, for instance, is associated to epileptic seizures whereas local synchronization is related to some cognitive tasks [28]. Synchronization in trading activities can be harmful for the stability of financial systems [15]. Thus designing mechanisms that prevent global synchronization can become a very important tool in complex dynamical systems [29].
Our model consists on a population of N moving oscillators with identical velocity V and random orientation on a square of side length L with periodic boundary conditions. The internal phases of the agents φ ∈ (0, 1) increase uniformly with period τ , until they reach a maximum value of 1, when a firing event occurs and the phase is reset. Upon such an event at time t, the firing oscillator influences its nearest neighbor (oscillator at minimal distance, labelled nn, see Fig.  1A) producing a random reorientation of its motion and an update in its phase by a factor ǫ: The phase update is performed at frozen time until the phases of all oscillators have been updated (some agents may reach their threshold and fire upon receiving a phase update from a firing neighbor). Then the phases evolve again uniformly in time until another update is triggered.
The system is synchronized when a succession of consecutive firing events (avalanche) equal to the system size N is detected. For the sake of clarity we define the (discrete) time T , as the number of times a given oscillator (that we identify with oscillator 0 in our computer simulations) has fired. This allows us to define T sync as the number of cycles this reference oscillator takes to enter the synchronized state (i.e., the number of updates needed for an avalanche of size N to occur).
The chosen minimal interaction rule is such that the system lies far below the static percolation transition [30,31]. Therefore, global synchronization is not achievable without motion (since it is very unlikely that a giant connected component of size N exists [34]). It is the nonnull velocity of the oscillators that enables the system to reach the coherent state, therefore we could expect T sync (the average time the system needs to synchronize) to be a decreasing function of V , such that T sync → ∞ when V → 0 and T sync → T f > 0 (a constant value) when V is high enough [35]. Fig. 1C shows that, for V < V s , the synchronization time decreases as a power of V when the latter increases. Then, the decreasing slows down and T sync has a minimum at V = V m > V s . Beyond this value, the synchronization time gets larger and larger, until the system enters a region where it is unable to reach the coherent state in a finite time (gray area in Fig. 1C). For even larger values of the velocity, V > V * , it decreases abruptly, finally reaching its asymptotic value when V = V f .
Three main regions are thus identified: A 'nosynchronization zone ' in the middle that separates a 'left region ' (small V ) and a 'right region ' (high V ). In the left region, we can separate two sub-regions: on the left, at V < V s , there is what we call the slow regime, and, on the right, a transition zone. The same happens for the right region. On the left, we find a transition region, while on the right, at V > V f , the system enters the fast limit, where T sync no longer depends on V . Before entering into a more specific discussion we present a general characterization of the system behavior.
In Fig. 2, upper panels, we introduce the cumulative individual interaction network (CIN) of an oscillator (labelled 0), for two independent synchronization process at two different velocities (panel C: V = 2V f while panel D: V = V m ). A CIN represents in a visual way the role played by a given unit in the signal spreading and is constructed in the following way: whenever an oscillator i fires at oscillator 0 and oscillator 0 fires at oscillator j, a link between i and j is added. If the link already exists, its weight is increased. In the case 0 does not receive any shot, we put a link between 0 and j. A reciprocal shot is represented as a self-link. We repeat this process until the system reaches the synchronized state.
To better understand the synchronization mechanisms we also report different relevant magnitudes. In Fig. 2, lower panels (panel A: V = 2V f , panel B: V = V m ) we plot the global order parameter η(T ) = cos (2πφ i (T )) and in order to have an insight about what happens at local scale, the quantity λ(T ) = cos (2πφ 0 nn (T )) where φ 0 nn is the phase of the unit to which the reference oscillator fired at time T . These parameters measure the synchrony of our system, ranging from a uniform phase distribution of our oscillators (η(T ) ≈ λ(T ) ≈ 0) to complete synchronization (η(T ) = λ(T ) = 1). Note that the averages are calculated with respect to all oscillators' phases {φ i } in the case of the global order parameter and only with respect to the nearest neighbor's phase φ 0 nn in the case of the local one [36].
Finally, m(T ) represents the total fraction of oscillators that have been out-neighbors of the reference oscillator up to time T .
We observe strong differences between the two scenarios. When agents move fast, all nodes appear represented in the CIN since they play a global role sending and receiving signals throughout the whole system. Interactions are completely rewired at each time T , therefore m(T ) increases very rapidly and φ nn is just a random variable extracted among N − 1 possible ones. This means that λ is exactly the same as η, but with less statistics. Both quantities increase together (more or less noisily) because by means of the firing events the whole phase distribution becomes narrower.
In contrast, when the mobility of the agents is reduced, few nodes conform the CIN. Each oscillator plays a local role mediating the interactions among a small number of units that are the same all the time, no matter how long the synchronization process could be. The behaviors of η and λ appear to be uncorrelated, being λ = 1 almost all the time: each oscillator spends a long time with its neighbors, usually being able to synchronize with it before changing. At the beginning (T < 500) whenever it starts firing towards a new oscillator (black vertical lines in Fig. 2) λ experiences an abrupt decreasing while m(T ) increases as it is the first time the reference oscillator meets that neighbor. Later (T > 500), the chances to change a neighbor for another one already known, having a very similar phase, increase. The phase distribution becomes narrower and specially at local scale, among the units the oscillator of reference can meet, the dispersion is small. Consequently neighbor changes do not affect λ anymore.
The observations point out that the fast regime can be understood as an homogeneous regime while slow velocity enhance heterogeneity among units. Hence the mechanisms that allow the system to synchronize have to be different in the two cases. The system has different strategies to reach the coherent state in the left and in the right region but neither of them work in the intermediate region.
In the following paragraphs we will determine the region of the parameters space corresponding to each regime providing a quantitative estimator for both V s and V f as functions of N , ǫ, and L. Where yellow lines are so dense that they form a yellow band, it means that the oscillator has more than one in-neighbor simultaneously.
In order to understand the different regimes we need to compare the time scales of the dynamical processes that lead the system to its final synchronized state, i.e. how fast are oscillators to change neighbors effectively and how fast they are able to locally synchronize.
Let us first focus on the fast regime. As a kind of mean-field assumption we can imagine a single oscillator, i, moving on a plane where the rest of the population is fixed in the most disperse possible configuration, a square lattice such that the system is divided into (N − 1) squares of side length ℓ = L/ √ N − 1. Oscillator i changes its neighbor when it exits a square to enter one of the adjacent ones. It will do so between time T −1 and T if, at time T − 1, the component of V perpendicular to that side is larger than the distance s separating i from that same boundary. First, notice that if V > ∼ ℓ then the exit probability is equal to 1 and each unit changes its neighbor at each time step. One can conclude that ℓ is a good estimator for the value of the velocity at which the system enters the fast switching limit, i.e.
On the opposite limit, when oscillators move very slowly (V √ N /L ≪ 1), we can also compute the average time spent in changing neighbors (see the Supplemental Material [32]) This time needs to be compared with the time a characteristic cluster needs to locally synchronize [20]. Given the minimal rule of interaction we have proposed, the average size of (weakly) connected components in the resulting network of our model is S cc = (8π + 3 √ 3)/3π ≈ 3 (see [33] and [32]). For such a characteristic cluster of size 3 there is only one possible configuration [37] (see Fig. 1A) whose synchronization time is given empirically by the expression with κ = 2.0±0.1 [32]. A necessary condition for the system to reach complete local synchrony is that, on average, no topological change has to occur during a time frame T 3 . Since the average time separation between two neighbor changes of the same unit is T out , the average time separation between two topological changes in the system can be expressed as T (N ) out = T out /N ≃ πL/(4V N 3/2 ). With these two ingredients, we can easily identify the starting point of the slow regime as that set of parameters values such that T 3 ≃ T (N ) out . Hence, fixing N , ǫ and L, we obtain Both V f and V s are reported in Fig. 1 and have been checked for various values of ǫ, N , and L, showing a very good agreement. In particular, in Fig. 3 we show how our predictor mark the transition point of the slow regime. Moreover, we can identify a universal minimum for all the curves marking the value of V m .
Summarizing, we have proposed a model of moving integrate-and-fire oscillators that interact only with the nearest neighbor. This minimal interaction rule can account for some physical situations where communication between agents is minimized. We have found three very different dynamical regimes. In the region where oscillators move very slowly and oscillators tend to keep the same neighbors for a long time, local synchronization dominates and a global one proceeds through a very large number of interactions (firings and changes of neighbors). In the other limit, of very fast motion, oscillators fire at random neighbors which makes that interactions to be very effective and phases approach each other in a monotonous way and synchrony is achieved directly on the global scale. Interestingly, we find an intermediate regime where none of these mechanisms work and the system is not able to reach synchronization. The presence of a no synchronization band needs to be further studied and its robustness checked for other types of non-linear interacting oscillators in future works.
Despite synchronization is usually seen as a positive outcome of a cooperative dynamical system, under some circumstances, as it happens, for instance, in brain dynamics, global synchronization is not a healthy state. Our model represents a paradigmatic example of a system where synchronization can be prevented by tuning the rate at which the topology changes. We have determined the bounds of this region in terms of the characteristic times of the model, those related to local synchronization and to mobility.