C OMPLEX network [1] , [2] is the unifying paradigm behind many scientific disciplines ranging from physics, mathematics to computer science and biology. Besides the study related to local and global statistical properties of complex networks [3] , researchers are equally interested in paying attention to the dynamics of their interacting units [4] . Among those collective dynamics, a widely studied phenomenon is that of synchronization [5] , [6] , [7] , [8] , [9] , [10] , [11] , [12] of coupled oscillators arranged over the constituents of complex networks. Although, in recent times, there have been attempts in examining synchronizability of time-varying networks with different underlying structural arrangements, in most of these studies of synchronization, the nodes of the complex networks are spatially static. But, in real world networks, the individuals (nodes) move around and exchange information within close proximity. Such nodes in the contextual literature are well known as mobile agents [13] , [14] , [15] , [16] , [17] , [18] , [19] , [20] , [21] referring to those nodes whose movement highly influence the systems' collective dynamics. But, in general, the states of those oscillators situated on top of those agents do not affect the agents' mobility. The studies on the process of synchronization in mobile agents [22] , [23] , [24] , [25] , [26] , [27] , [28] have, so far, primarily concentrated on those particular cases when the mobile agents interact with each other only when they are closed enough in space. But, in all of these studies, the network architecture, i.e., how the oscillators are connected with each other, depends only on the attractive interaction. So, a basic query arises that how is the flow of information gets affected when different types of interactions are present in the movement. For instance, an attempt to answer this question can be made while incorporating repulsive coupling with the attractive one in such a networked system. Such co-existing dynamical interactions can be observed in a variety of systems and have been studied previously in different network setups yielding diverse collective states [29] , [30] , [31] , [32] , [33] .
Synergy of attractive-repulsive coupling can lead to several interesting dynamical behaviors. One of such example is the emergence of solitary states in multiplex networks [34] possessing positive inter-layer and negative intra-layer interactions. Earlier, Maistrenko et al. [35] found a variety of stationary states in networks of globally coupled identical oscillators with attractive and repulsive interactions. Sathiyadevi et al. [36] , [37] found diverse chimera states due to the introduction of nonlocal repulsive coupling together with an attractive coupling in a network of coupled oscillators. Also, interplay between conformists and contrarians leads to traveling wave along with the stationary states [38] . The effect of competitive interactions on the synchronization manifold among globally coupled identical Van der Pol oscillators through their velocities is rigorously addressed in ref. [39] . Even, the presence of a tiny fraction of repulsive couplings is found to enhance the synchronization of nonidentical dynamical units that are attractively coupled in a small-world network [40] . Origination of dragon-king-like extreme events as a result of the coexistence of excitatory and inhibitory chemical synaptic couplings is reported in the ref. [41] . Besides emergent dynamical patterns, co-presence of attractive and repulsive couplings carry high significance to several biological [42] , physical [43] , ecological [44] and even social scenarios [45] .
Recently, Nag Chowdhury et al. [46] proposed an universal 0 − π rule to identify the bipartiteness of networks from the anti-phase synchronization states. This article describes such a system of attractively coupled oscillators with negative links as frustrated networks [47] . Such frustrated systems under the influence of attractive-repulsive co-existing interactions may induce unanticipated phenomenon, like extreme events [48] . Most of the existing articles [41] , [49] , [50] , [51] , [52] on extreme events in dynamical systems and networks define such an event as rare event. Mathematically, a rare event is defined on a probability space χ equipped with a σ-algebra B of events such that P (A) is small, where P is the probability measure designed to quantify the probability of the occurrence of any event A ∈ B. But, extreme events can be rare events [53] , [54] , [55] , or they can often be frequent in time and space [56] , [57] , [58] , [59] . Examples of such intermittent frequent extreme events include population size of Paris [60] , the bubbling of share market [61] before a crash and many more.
Along with its small occurrences, extreme event has another threat in the form of large impact. Financial and commodity market crashes, tsunamis, hurricanes, floods, epidemic disease spread, global warming-related changes in climate and weather, warfare and related forms of violent conflict, asteroid impacts, solar flares, acts of terrorism, industrial accidents, 8.0+ Richter magnitude earthquake are few examples of such large impact events. In dynamical systems, researchers are recently starting to measure the impact of such events in terms of the amount of variation from the central tendency of that observable [41] , [49] , [51] , [62] , [63] , [64] . So, if ψ is the observable defined on the state space U to R, then the extreme events are members of the set {u ∈ U : ψ(u) > H S }, where H S = m + dσ with d is the non-zero integer. Here, m and σ are respectively the mean value and the standard deviation of all the peak values in a time series of u. But, if H S is too large, there will be very few values to model the tail of the distribution correctly as the variance is likely to be large due to only very extreme observations remaining. On the other hand, a low threshold value of H S will include too many values giving a high bias. Recently, m+8σ is found to appropriate indicating extreme event threshold for Weibull distribution [65] .
Inspired by these observations, we consider a frustrated network of limit cycle oscillators, where each node of the network is a mobile agent. Earlier all studies [22] , [23] , [25] , [27] , [28] , [65] , [66] related to mobile agents, considered only limited local interactions. Instead, we consider here a global network of mobile agents with co-existing switching interactions. Depending solely on the relative distance between the agents, attractive or repulsive interaction is activated. Specifically, we choose attractive interaction between those agents which are staying apart from each other, while we consider repulsive coupling for the agents that are sufficiently close. We consider two types of coupling schemes, namely symmetry-breaking and symmetry-preserving couplings. Under this set up, we explore the effect of coupling strengths which reveal quite interesting phenomena like cluster synchronization [67] , complete synchronization, inhomogeneous small oscillation, extreme events etc. We also analytically derive the critical coupling strength for achieving complete synchronization using time-average Laplacian matrix in absence of repulsive interaction and numerically verify the results.
We consider N mobile agents which are moving independently on a two-dimensional (2D) physical plane S = [−g, g] × [−g, g]. Initially, the agents are randomly distributed on S. The i-th mobile agent can move in any direction with velocity v cos θ i (t), v sin θ i (t)) (any kind of collision among the agents are not allowed), where v is the uniform moving velocity and θ i is drawn randomly from [0, 2π]. Higher values of velocity v implies that the agents move the whole physical space and which increases the possibilities of interactions. Thus, if p i (t), q i (t) is the position of the i-th agent at time t, then motion updating process will be maintained by the following relation,
(1)
We define the distance between any two agents by the standard Euclidean metric as
Next, we write the the dynamical equations characterizing the time-evolution of each agent (i = 1, 2, · · · , N ) in the networked system as follows,
where x i ∈ R n is the state variable and F (x i ) is the corresponding vector field of the i-th agent. Here H(x i , x j ) is chosen as the diffusive type coupling function, k A > 0 and k R < 0 respectively correspond to the attractive and repulsive coupling strengths. Moreover, B and C are the distance dependent adjacency matrices associated to the attractive and repulsive interactions such that
where α is the parameter lying within [0, 2 √ 2g]. Here 2 √ 2g is the maximum possible distance between any two agents lying in the 2-dimensional plane [−g, g]×[−g, g]. At a particular instance, the agents whose Euclidean distance is less than α, repel to each other 
Note that the numbers of attractive and repulsive neighbors both are widely spread. The simulation is carried out for t = 3 × 10 5 . (c, d) Variation of degree of both matrices with respect to time t: Just like the upper panel, we find that the neighbors of any agent are rapidly changing depending on the agent's random movement. Here, we take g = 10 and α = 10, and the uniform moving velocity v = 2.0. and repulsive interactions arise between them. Attractive coupling occurs between the agents when the distance between them is greater than α.
To illustrate this, we draw a schematic diagram with N = 3 mobile agents in the 2-dimensional plane [−g, g] × [−g, g] in For this, we calculate the number of attractive and repulsive neighbors depending on the critical distance α. In Fig. 2 , the degree distributions of moving agents for attractive and repulsive interactions are shown. We randomly select an agent which is moving with uniform velocity v = 2.0. We then find that the number of attractive neighbors are varying within the interval [5, 92] (cf. Fig. 2(a) ) and the number of repulsive neighbors of that agent is varying within [7, 94] (cf. Fig. 2(a) ). At any particular time-step, the sum of the number of attractive neighbors and the number of repulsive neighbors is N − 1 for any i-th agent. This scenario is well expressed by the lower panel of Figs. 2(c, d).
Now, on top of each mobile agent, a two-dimensional Stuart-Landau (SL) oscillator is placed where the state dynamics of each limit cycle oscillator is represented by
where ω i = ω = 3.0 is the identical intrinsic frequency for each i = 1, 2, · · · , N . We consider two types of diffusive interactions, namely symmetry-breaking coupling H(
The interesting and important part of this distance dependent coupling functions is for α = 0, the Eq. (3) becomes globally coupled network with attractive coupling only. In this case, movement of the agents does not effect the attractive adjacency matrix B and the network becomes static. But, for α > 0, depending upon the values of attractive and repulsive coupling strengths k A and k R , two different cases can be implemented: 1) k A = 0 and k R = 0 : In this case, the mobile agents will interact attractively when the relative distance d ij is greater than α and remain disconnected with each other when they are closed to each other and relative distance d ij is less than equal to α. There is no repulsion among the agents. Here, synchronization is the most desired state for a suitable value of k A . 2) k A = 0 and k R = 0 : In this case, interplay between attractive and repulsive interactions is considered. During the spatial movement of the agents, they experience both types of interaction among them. Here, different states may arise, like inhomogeneous small oscillation, synchronization, extreme event and other intermittent states.
In the following section, we will explore the dynamics of the time-varying network (3) with two types of diffusive couplings and under above mentioned distinct coupling conditions. Our main emphasis will be to identify the parameter space by varying the two coupling strengths k A and k R with fixed values of other network parameters g = α = 10.0 and v = 2.0 related to the movement of the agents.
For numerical simulations, we integrate Eq. (3) using fifth order Runge-Kutta-Fehlberg method with a integration time step δt = 0.01. At each integrating time step, θ i (t) and thus v i (t) is updated according to rule given in (1) for each i-th mobile agent (i = 1, 2, · · · , N ). Hence, both the matrices B and C change at each integrating time step leading to fast switching approximation [68] . All simulations are done for fixed initial conditions
unless stated otherwise 1 . In this Sec. (3), the symmetry-breaking
For the chosen values of the parameters as mentioned above, we find different dynamical behaviors of the network system (3). We first look at the collective dynamics of the network whenever there is no repulsion among the agents and only attractive coupling is activated, i.e., when k A = 0 and k R = 0. In this scenario, the agents with relative distance d ij (t) ≤ α = 10 remain disconnected with each other and agents will interact with each other if their distance is greater than α. We next intend to analytically derive the critical interaction strength k critical for which synchrony appears, based on the approach of constructing the time-average Laplacian matrix G = [g ij ]. For the sake of simplicity, let us start with N = 2 mobile agents. Then, depending on their relative distance, two possible Laplacian matrices can be observed. 1) Whenever d ij > α, then the agents interact with one another so that the corresponding Laplacian matrix becomes
2) If d ij ≤ α, then the Laplacian matrix due to absence of individuals' interaction takes the form
Since G 0 is a null matrix so in this case G reduces to G = pG A , where p is the probability of interaction between the two agents. Since, we are considering a planar space S of area 4g 2 sq. units, then p = 1 − Area of the interaction Area of S becomes p = 1 − πα 2 4g 2 .
1. Inhomogeneous small oscillation depends crucially on the initial conditions, and the basin of attraction of this state is small. If the initial conditions are not chosen suitably from the basin of attraction, then the trajectories may converge to a different periodic attractor. However, other observed dynamical states can be reached from any random initial conditions starting
We next calculate the average Laplacian matrix G for N = 3 agents. For instance, depending on the spatial positions of these mobile agents for a fixed value of α, eight possible configurations have been encountered for N = 3. These possible cases are described below:
2) There are three cases in which two out of the three agents interact with each other and the third one remains isolated. The corresponding Laplacian matrices are
Laplacian matrix when i-th and j-th agents are interacting with each other, and the other third agent remains isolated. Also note that, G 12 +G 13 +G 23 = G A .
3) It is also possible that two of the agents lie within the relative distance α so that they are not interacting with each other, but the third agent stays far away and interacts with both of them. The associated Laplacians are as follows
where G k is the Laplacian matrix when the k-th agent interacts with the other two agents, but those two agents do not interact with each other. Here note that G 1 + G 2 + G 3 = 2G A . 4) Majority of the agents, i.e., all three oscillators situated within the relative distance α, hence no one can establish a connection with the others. So, Laplacian matrix can be indicated by
Then the time-average matrix for N = 3 becomes,
Here p's stand for the probabilities associated to the respective network configurations for which we have p 1 = p 2 = p 3 = p a (say), and p 12 = p 13 = p 23 = p ab (say). Then G = (2p a + p ab + p A )G A and hence G = pG A , where p = 2p a + p ab + p A is the probability of interaction between any two agents. The scenarios corresponding to N ≥ 4 can be similarly tackled, all of which yields the time-average matrix as G = pG A , where G A is the Laplacian matrix of order N × N in the form,
· · · · · · · · · · · · . · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · −1 −1 −1 · · · · · · · · · N − 1
so that G is basically a rescaled all-to-all global Laplacian matrix.
As per our described model, the links among the agents get rewired at each integrating time step, so the convergence of all the attractors to a single attractor [69] can be assured if the timeaverage matrix supports synchronization [17] . As stated in the Ref. [17] , if there exists a constant T such that G(t) satisfies 1 T t+T t G(τ )dτ = G, the time-average of G(t) and the system of coupled oscillators given bẏ
with fixed interacting Laplacian matrix G = [g ij ] possesses a stable synchronization manifold, then there exists * > 0 such that for all fixed satisfying 0 < < * , the coupled system according to the time-varying Laplacian G( t ) defined bẏ
.., N (7) also sustains a stable synchrony manifold under sufficiently fast switching sequence between the network configurations.
Then, the stability of the synchronized state can be investigated by the eigen values of G. For a time-independent coupling matrix, a necessary condition for synchronization is that the master stability function (MSF) be strictly negative in each transverse direction [70] . With our choice of parameter ω and the coupling function H(x i , x j ) = (x j − x i , 0) T , the SL-system belongs to class-I MSF [70] , i.e., the synchronization manifold is stable in the interval of type [β 1 , β 2 ]. The N eigen values of G are λ 1 = 0 and λ j = pN , j = 2, 3, ..., N . Then the critical range of k A to attain complete synchrony is obtained by the inequality β 1 ≤ k A λ j , j = 2, 3, ..., N which implies, To verify the relation (8), we plot N − k A parameter space in Fig. 3 where the red solid curve is our analytically found curve k A = β1 pN for β 1 0.0111. This curve agrees well with our numerically found synchronization region (blue shaded region).
In order to numerically study the emergence of synchrony, we define the synchronization error as
where · t represents time average obtained over a long time interval (taken as 3 × 10 5 time steps here) after initial transient of 3 × 10 5 time steps. In our work, complete synchronization corresponds to the state when E becomes less than 10 −10 accounting 50 independent network realizations. Note that, we find a small fluctuations in the relation (8) for finding k critical . These fluctuations may occur due to various reasons, including insufficient realizations, approximate choice of the value of β 1 , or for the effect of the uniform moving velocity v. Even, this small fluctuation may be due to our choice of p. Although, for static agent assumption or v = 0, our chosen p = 1 − πα 2 4g 2 works absolutely correctly. But when the agents are moving, the agents do not necessarily maintain this exact relation (8) . In fact, in deriving the relation (8), we use the negativity of the maximum Lyapunov exponent in the transverse direction of the synchronization manifold, which is a necessary condition, not a sufficient one. All these lead to small fluctuations in k critical compared to our derived critical synchronization curve k critical = β1 pN . Notably, this density dependent threshold for the emergence of synchrony is found to be relevant and is of practical interest, particularly in the studies of the bacterial infection, biofilm formation and bioluminescence where quorum-sensing transition is observed in indirectly coupled systems [71] , [72] .
This result is quite important due to the fact that if we randomly distribute N agents in the plane S and then set v = 0 (i.e., when the network becomes static), complete synchronization is not guaranteed depending on their initial placements. As the connected agents may exhibit complete synchronization depending on suitable coupling strength k A , but the disconnected agents remain isolated and as a result of that they maintain isolated trajectories. But, the mobility of each agent allows them to interact even if they are long-distant, under sufficiently long simulations. These occasional interactions induce complete synchronization among them for appropriate critical coupling strength.
Such occasional minimal interaction is found to be beneficial for the emergence of complete synchronization, instead of continuous-time coupling, from the point of view of optimal interactional cost [73] , [74] . This confined interaction may be useful in case of robotic networks and wireless communication systems, where transmission signals are turned on if the agents are lying sufficiently close to each other and where we will have to inspect if synchronization can still occur.
To analyze the effect of only attractive coupling, one requires to turn-off the repulsive coupling k R . But this is not the only possible way of setting on repulsion free interaction. Instead, one may realize the scenario by setting α = 0. This specific choice of α transforms the repulsive matrix C to null matrix, and attractive matrix in the form 1 · · · · · · · · · 1 1 0 1 · · · · · · · · · 1 1 1 0 · · · · · · · · · 1 · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · · 1 1 1 · · · · · · · · · 0
Note that the matrix B now becomes a static, i. e., a timeindependent matrix. Under this set-up, one may expect complete synchronization, but for higher values of k A , we can not obtain convergence of those attractors into one attractor. Instead of that, we find incomplete synchronization in the form of cluster synchronization for k A ≥ 0.19. We define U = {x i | i is odd} and V = {x i | i is even}. After the initial transient, the x i of trajectories of the sets U and V converge respectively to two different trajectories as shown in Figs. 4(a, b) for k A = 1.0. Both the sets U and V contain 50 oscillators. The global network evolves into two equipotent subsets of oscillators in which members of the same cluster are synchronized to the same trajectory, but members of different clusters, oscillating with different small amplitudes. Notice that the oscillators that start with negative x i as per relation (5) Such partial synchronization may show up in swarms of unmanned autonomous vehicles, power grids and swarms of animals [75] . This interesting phenomenon is also hold for random initial 1] . But, in that case, the cardinality of U and V vary significantly in each realization. We also find that, if we increase k A by keeping fixed the initial conditions as per the relation (5), still we do not observe complete synchronization for k A ≥ 0.19. They still show cluster synchronization for suitable k A . However, for smaller suitable k A , the trajectories exhibit complete synchronization converging to a single periodic attractor. For N = 100 oscillators, complete synchronization can be obtained for k A < 0.19.
Now, we consider α = 10 along with k A = 1.0 and k R = −0.10, so that both the couplings can play a decisive role. Under this choice of parameters, we observe the time-series of each oscillator exhibits chaotic behavior. Again, those x i of N = 100 oscillators distribute equally on both sides of y i = 0. Similarly, y i of those oscillators are distributed equally on both sides of x i = 0. One of such scenario is shown in Fig. 5(a) . It is observed that each of these two groups is not synchronized. They maintain their chaotic small oscillations. The simulations are done using initial conditions of relation (5) .
Recently this type of oscillations is reported by Dixit et al. [76] for two coupled SL-oscillators under the influence of dynamic interaction. But, our study is quite different from them. Our network consists of mobile agents, and their relative distance d ij (t) actually decides the mutual interaction type, which is completely different from the study [76] . In fact, we found this state for large networks, which help to deny the finite-size effect. Moreover, each split group in our case consists more than 2-oscillators and those trajectories do not lead to any cluster synchronization.
To study the effect of network size N , a detailed numerical study is perceived. We find that there exists a N critical , beyond which network shows such inhomogeneous small oscillation. Numerically, we find this N critical is 51. For N ≥ 52, such chaotic small oscillation is observed. For N < 52, the attractors are converged to a single attractor and the trajectories are exhibiting periodic behavior. In Fig. 5(b) , we plot the x i 's of all N = 51 oscillators, which are oscillating synchronously. To scrutinize the effect of k R , we fix the attractive coupling k A = 1.0, and slowly increase the absolute magnitude of k R . We observe that with increment of k R , the amplitude of the oscillators increase. After a certain value of k R , the attractors of one group, which are initially separated by the axes (See Fig. 5(a) ), are crossing the axes frequently. The trajectory of i-th oscillator crosses the axes and joins the other group for sometimes, before it comes back to its original group. To illustrate this, we choose randomly an agent from N = 100 oscillators, and plot its time series in Fig. 5(c) . The corresponding attractor is also shown in Fig. 5 (d) to demonstrate the chaotic switching behavior of the attractor. Interestingly, further slight increment of k R leads to complete synchronization. At k R = −0.28, all N = 100 oscillators oscillate synchronously and re-gain their isolated periodic dynamics. They oscillate in the limit cycle regime from [−1, 1].
Further decrement of k R (for fixed k A = 1.0) destabilizes the synchronization manifold and after a certain value of k R , the error dynamics becomes intermittent. Specifically, for most of the time the synchronization error E (cf. Eq. (9)) remains near zero (lies within the interval [10 −10 , 10 −6 ]), but occasionally this E becomes non-zero. For a narrow range of k R with fixed k A = 1.0, we find this intermittent error trajectories exhibiting extremely large values. To reveal this issue, we plot variation of synchronization error E in Fig. 6 , which is accumulated over a sufficiently long time interval. For this entire data, we calculate the sample mean m = 0.0011 and the standard deviation σ = 0.0199. We then define a extreme event threshold H S = m + 8σ 0.1606, which is plotted over the error dynamics with horizontal dashed line in Fig. 6 . In this particular accumulated data of E, we find 2570 data points having value of E, which is greater than H S . Thus, the probability of {E > H S } is 0.0008567. However, this probability varies with each realization, and for our choice of H S = m + 8σ, Cantelli's inequality [77] , [78] yields the upper bound for each realization as
In the existing literatures [79] , [80] , this type of intermittency is known as on-off intermittency. This type of intermittent error, the trajectories experiencing non-uniform, uncertain extensive expeditions from the synchronous manifold and can be considered as an appropriate candidate for extreme events. Earlier, local instability of synchronization manifold due to on-off intermittency causing extreme events in a pair of coupled chaotic electronic circuits is reported [81] . Also, Nag Chowdhury et al. [65] found a similar mechanism for the generation of extreme events in a mobile network of chaotic oscillators. To further resolve whether the spikes in Fig. 6 qualify for extreme events or not, we determine the local maximum values E n of E and plot the histogram for the event sizes E n in Fig.  7 . Note that the length of the collected time series is sufficient enough so that due to statistical regularity, inclusion of new sample does not affect the structure of the histogram. Pisarchik et al. [82] pointed out a special characteristic of extreme events, which is the L-shaped probability density function (PDF) of E n in the semilog scale. In spirit of all these facts, we define extreme events in this article as follows, 1) These short lasting events are recurrent, aperiodic and of different amplitudes higher than the significant height H S = m + 8σ. In order to further eliminate the small amplitude events, another restriction is maintained, H S > 0.01 is taken along with H S = m+8σ [65] , [82] . 2) They appear much more often than Gaussian statistics and they are unpredictable almost surely with respect to time, where an event is defined as almost surely if the set of possible exceptions may be non-empty but has zero Lebesgue measure [83] .
3) The appearance of these large events would have very small probability (though the exact measure of small needs to be quantified) [48] .
Clearly, the histogram in Fig. (7) is an approximate representation of the distribution f (x) of the numerical data. In the next, we define a return interval k. Mathematically, a return interval k occurs if X i > H S and X i+k > H S , but X j ≤ H S for i < j < i + k, where X i is the value of the observable X at the ith step. Suppose, W H S is the number of total return intervals. So, we have k i for i = 1, 2, · · · , W H S such that not necessarily all the k i are distinct. Clearly, if k i is short then a cluster of accumulated extreme events can be observed. In contrast, a large value of return interval definitely portrays few occurrences of extreme events. Also, let k 1 + k 2 + · · · + k W H S = W , then the mean return interval is
This relation (11) clearly indicates that the mean return interval R H S is relatively short if the occurrence of extreme events is frequent (i.e., W H S is high). Similarly, if W H S is small then the mean return interval is relatively high. But, this number of return intervals W H S is dependent on various factors and thus it varies on realizations. Hence, relation (11) gives different average return interval R H S based on different realizations. Using time series analogous of Kacs Lemma [84] , [85] , [86] , R H S can be described in terms of the tail of the normalized distribution density f (x) as
Thus, there exists a one-one relation between the chosen threshold H S and the mean return interval R H S , which is solely determined by the normalized distribution f (x).
For a complete understanding of interplay between attractive and repulsive interactions (i.e., for α = 0), a two dimensional k A −k R parameter space is presented in Fig. 8 . In this figure, we vary
Initially, in absence of repulsive coupling, i.e., for k R = 0.0, complete synchronization is observed for suitable choices of k A where E < 10 −10 . As per this Fig. 8 , whenever k A reaches to 0.5, then the introduction of small k R < 0.0 may destroy the limit cycle behavior of the oscillators and their synchronized behavior is also lost. The reflected region of inhomogeneous small oscillation enlarges with further increment of positive coupling k A . Interestingly, for a fixed k A ≥ 0.5, the decrement of k R helps to restore the periodic behavior of the attractors. With suitable choices of k R for fixed k A ≥ 0.5, this transition from inhomogeneous small oscillation to complete synchronized state occurs through the path of attractor switching, as already shown in Fig. 5 . The attractors are quenched and with further decrement of k R for appropriate k A , the trajectories switch between the co-existing attractors and finally regain their respective periodic synchronized attractor. When k R is further decreased, the synchronization error E becomes intermittent and reveals occasional away journey from the synchronization manifold. This regime is depicted through the red and blue regions in the Fig. (8) . In order to distinguish between these red and blue regions, we calculate the significant height H S . When H S lies within (10 −6 , 0.01], we label those states as the intermittent region (in red). For H S > 0.01, we define those states as extreme events which obey additional few clauses mentioned earlier. While in complete synchronization, all The yellow and gray regions respectively correspond to the states E < 10 −10 0 (i.e., complete synchrony) and the inhomogeneous small oscillation. The blue area stands for the extreme events whereas the red zone represents the state where the error E becomes intermittent but H S ≤ 0.01. Finally, the area in cyan corresponds to the desynchronization state. Other parameters are N = 100, g = α = 10 and v = 2.0. Here, 100 independent realizations are used to obtain this parameter space.
the attractors converge into a single attractor. But, during the appearance of extreme events, most of the time all the trajectories evolve over a converged single trajectory, but occasionally due to the local repulsion (of suitable strength k R ) few trajectories leave the converged attractor and follow their own path for a relatively short time. Hence, during those time, the error trajectory leaves the equilibrium point E = 0.0 (the synchronization manifold) and exhibit a large excursion exceeding H S = m + 8σ. This choice of H S is influenced by the Ref. [65] , where the authors gave analytical logic behind the choice of such extreme event indicating threshold. The second variables y i of all N = 100 oscillators are plotted over a short time-interval and the corresponding error dynamics are also presented through an arrow in Fig. 8 . As expected, if k R is further decreased while keeping k A fixed, desynchrony states are observed.
In this section, we deal with the scenario in which the local dynamical units, i.e., the SL-oscillators are coupled through both the variables so that the coupling becomes (rotational) symmetrypreserving. Our main emphasis will be to examine how this change affects the collective dynamics of the networked system. We depict the dynamical outcomes through Fig. 9 in which the phase diagram in the k R − k A parameter plane is portrayed. As can be seen from the figure, the presence of only attractive coupling readily induces complete synchronization in the system. Moreover, in contrast to our earlier observation for symmetry-breaking coupling (cf. Fig.   8 ), this synchrony sustains as k A increases. This is mainly because now the coupling function does not allow to break the symmetry yielding inhomogeneous small oscillation. On the other hand, switching the repulsive coupling strength k R on, one observes occasional jumps in the error dynamics (for suitable strengths of couplings) that indicates the appearance of intermittent states. More importantly, further increment of k R leads to the intermittent states that satisfy the criteria of being considered as extreme events. This result demonstrates that our observations of complete synchronization or the extreme events on the considered network set-up is not limited to a specific choice of symmetry-breaking coupling function.
In this study, we have considered Stuart-Landau oscillators on top of mobile agents moving in a finite region of two-dimensional plane. Due to the spatial movement of the agents, the relative distance between any two agents always varies with time. Using this spatial distance between the agents, two competing interactions are introduced. Whenever the agents stay inside a closed ball of radius α, the oscillators are repulsively coupled bidirectionally through diffusive linear coupling. On the contrary, whenever the mobile agents lie beyond the closed communication ball of radius α, they are coupled through attractive coupling. Under this set-up of the networked system, we have encountered diverse collective behaviors and provided rigorous investigation of each of them.
Firstly, we have addressed the effect of sole attractive interaction among the dynamical agents. In order to accomplish this, we Fig. 9 . Phase diagram in the k R −k A parameter plane: The yellow region corresponds to the states of complete synchronization. The blue area stands for the extreme events whereas the red zone represents the state where the error E becomes intermittent, but H S ≤ 0.01. Finally, the area in cyan corresponds to the desynchronization state. Other parameters are N = 100, g = α = 10 and v = 2.0. Here, symmetry-preserving vector coupling among the local dynamical units is considered in which 100 independent realizations are used to obtain this parameter space.
first set the repulsive strength to zero for which complete synchronization is found for α = 10. This result is quite interesting as for a static time-independent network, complete synchronization is not possible. Mobility of the agents play the decisive role to converge all the attractors to a single one. Analytically, the critical attractive coupling strength for synchrony is calculated, which exhibit excellent agreement with the numerical results. To investigate further the sole influence of attractive coupling, the repulsive matrix is set to be a null matrix by keeping α = 0 fixed. This different strategy leads to complete synchronization for suitable choice of k A and cluster synchronization beyond a certain value of k A .
Next we dealt with our prime issue of coexisting attraction and repulsion. The coexistence of such competing interactions can demolish the synchronized behavior of the whole system. Under suitable choices of both coupling strength, inhomogeneous small oscillation is perceived. In this state, each oscillator maintains chaotic behavior though we set initially each oscillator into limit cycle regime. Further reduction of repulsive strength helps to regain their periodic behavior through the path of attractor switching. Moreover, further lessening leaves the error dynamics intermittent. This irregular away journey of error trajectory from the synchronization manifold gives rise to infrequent large deviated events. Using few characterizations, we confirmed these states as extreme events. The mean return interval of these extreme events is approximated using time series analogs of Kac's lemma in relation (12) . This relation indicates the dependency of average return interval R H S on the threshold H S . The larger R H S clearly corresponds to larger H S and vice-versa. The route for generating such large amplitude events is also addressed in terms of on-off intermittency. An upper bound for the probability of occurrence of extreme events is calculated in relation (10) depending on the choice of H S .
Further, to elucidate the mechanism behind inhomogeneous small oscillation of those oscillators, vector coupling is introduced. This coupling helps to preserve the rotational symmetry of SLoscillators. Although, this vector coupling yields behavior like complete synchrony and extreme events, but it eliminates the appearance of such small oscillation, which is occurred during symmetry-breaking coupling. We have been able to portray the whole scenario of all these emerging collective behaviors in a two-dimensional parameter space of the competing coupling strengths. This parameter space manifests the fact that by tuning any one of the coupling strength, one can avoid such devastating extreme events. This fact can be used as one of the influential strategy for controlling extreme events. Such controlled usage of antibiotics [87] , [88] already found to be helpful in avoiding worsening medical cost and mortality especially for life-threatening bacteria infections.
In this article, we analyzed the interplay of positive and negative interactions in coupled identical limit cycle oscillators, by considering SL system on top of each mobile agent. An important direction of future generalization may include networks of chaotic oscillators, which is where our approach might make a difference. The chaotic dynamics [89] , [90] , [91] , [92] are expected to reveal an even wider spectrum of interesting dynamical states. This remains the core interesting avenue for future work.