Exploring the Interplay of Excitatory and Inhibitory Interactions in the Kuramoto Model on Circle Topologies

In the field of collective dynamics, the Kuramoto model serves as a benchmark for the investigation of synchronization phenomena. While mean-field approaches and complex networks have been widely studied, the simple topology of a circle is still relatively unexplored, especially in the context of excitatory and inhibitory interactions. In this work, we focus on the dynamics of the Kuramoto model on a circle with positive and negative connections paying attention to the existence of new attractors different from the synchronized state. Using analytical and computational methods, we find that even for identical oscillators, the introduction of inhibitory interactions significantly modifies the structure of the attractors of the system. Our results extend the current understanding of synchronization in simple topologies and open new avenues for the study of collective dynamics in physical systems.


I. INTRODUCTION
The Kuramoto model (KM), first introduced by Yoshiki Kuramoto in 1975 1,2 is a mathematical model that describes the dynamics of a system of coupled phase oscillators.It is widely used to study synchronization phenomena [3][4][5] in a variety of physical, biological and social systems [6][7][8] .The model is based on the assumption that each oscillator has a natural frequency and that the oscillators interact with each other via a diffusive nonlinear term.
Much research has been done to understand under what conditions the system will synchronize.In the simplest scenario, there is a trade-off between order and disorder.Coupling favors coherence, while the distribution of frequencies in the population disrupts it.When the intensity of the interaction reaches a certain threshold, synchronization emerges spontaneously.The model has been analyzed using many different techniques over the last 40 years, and countless variants have been proposed and discussed in the literature 6 .
Mean-field solutions of this model have been extensively studied under very different conditions.In these studies, topology did not play a relevant role.The basic assumption was to consider a fully connected system.However, the network topology has a large impact on the properties of the system 7,9,10 .They enrich the dynamics of synchronization phenomena.In contrast to simple lattice or fully connected network configurations, complex networks can exhibit features such as strong clustering, community structure and scale-free or small-world properties 11 .These topological complexities have profound effects on the system's ability to achieve synchronized states.For example, hubs can act as synchronization guides and effectively control the behavior of their connected neighbors 9,12 .The focus on KM in the context of complex networks aims to show how these complicated network architectures influence collective behavior.This makes it an indispensable tool for the study of real systems where network topology is a crucial factor 7,9 .
In situations where all oscillators are identical, one would expect the fully synchronized state to be the only attractor of the dynamics.However, in some very specific topologies, other attractors may also appear.For a sparsely connected network, it can be difficult to achieve full synchronization because isolated or poorly connected clusters act almost independently of each other.This fact usually leads to new locally stable solutions with complicated basins of attraction.In contrast, a densely connected network is more likely to converge to a single synchronized state, as each oscillator is influenced by a large fraction of the others.Some authors [13][14][15][16] have focused on the density of either random or regular networks and analyzed the necessary condition related to density that ensures the existence of a single attractor.
For this reason, the KM is often studied in its simplest form on a circle in order to understand the basic principles of synchronization.In this configuration, the oscillators are arranged along a circle, and each of them interacts with each of its nearest neighbors in a homogeneous way.This setup allows for analytical solutions and serves as an introduction to understanding how coupled phase oscillators can reach a phase locked state over time.The dynamical properties of these local solutions can exhibit a surprising degree of complexity.In this context, it is worth mentioning the works that introduced the concept of winding number to classify the phase-locked states 17,18 .
This concept is of crucial importance for the characterization of the steady-state dynamics of coupled oscillators, no matter whether they are synchronized or phase-locked states.However, determining their basins of attraction is not straightforward, as we are dealing with complex, multidimensional landscapes in which slight changes in initial conditions or network parameters can lead to completely different final states.][21] Most previous studies have focused on KMs with only positive (excitatory) connections between oscillators 22 covering different scenarios like unidirectional couplings 17,23 , longer range connections 14,18 , distribution of frequencies 24,25 , delays 26 , possibility of adding or removing nodes 27 , acyclic graphs 28 , planar graphs 29,30 , etc.However, the addition of negative (inhibitory) links also leads to a more complex and diverse range of behaviors.
In this paper, we investigate the effects of positive and negative connections on the synchronization properties of the KM.By including both excitatory and inhibitory interactions, we attempt to shed light on the interplay of these opposing forces and the resulting dynamics.We show that a circle with an arbitrary number of positive and negative connections can only be divided into two classes, despite the apparent added complexity: any circle with an even number of negative connections is equivalent to a system with all positive connections, while the circle with an odd number of negative connections is equivalent to one with a single negative link.One of the contributions of the present work is to generalize the concept of winding number.When only positive connections are present, the winding number is an integer quantity.However, we will show later that we need to introduce a new index when the number of negative connections in a circle is odd, which takes into account the symmetry of the system and ensures that the winding number can now become half-integer.KM has found application in many different areas, but one of the most successful is the dynamics of the brain as a first nonlinear approximation.In this field in particular, positive or negative links can play an important role, as they ensure that some attractors, such as complete synchronization in networks of identical oscillators, can no longer be achieved and other attractors appear that could possibly be associated with functional patterns in the system. 31et us notice that negative interactions have been also considered in the literature in a different context.For instance, in 32 the authors study a generalized KM with oscillators characterized as "conformists" and "contrarians" based on positive and negative couplings to the mean field.While our research is primarily concerned with the dynamics of positive and negative couplings, their study sheds light on related but complementary synchronization phenomena in networked systems.
The paper is organized as follows: Section II presents the fundamentals of the model and discusses the equivalence among circles varying in the number of negative links.Section III identifies the fixed points, while their stability is examined in Sect.IV.The paper concludes by thoroughly investigating the system's attractors through extensive computer simulations and summarizing the findings.

II. THE MODEL
In this paper, we consider the case of identical oscillators: where θ i represents the phase of the i-th oscillator, ω is the natural frequency, that can be taken as zero without loss of generality 6 , K is the coupling strength that fixes the timescale but for our purposes plays an irrelevant role and will be taken as 1, N is the number of oscillators, and a ij are the elements of the adjacency matrix (1 if i and j are connected, 0 otherwise).
In the particular case of a circular geometry (see Fig. 1) each oscillator is connected to its two nearest neighbors (1) where the indices of the adjacency matrix are taken mod N .The sign of the adjacency matrix elements can be either positive (excitatory) or negative (inhibitory).In this particular case and taking into account the symmetry properties of the sine function the general problem with an arbitrary number of positive or negative links can be reduced to a system with only positive links or to a system with a unique negative link as we are going to show.
Consider first that there is an even number of negative links, and let's introduce a new set of variables {φ} linearly related with the original set {θ} as follows.Choose a random node in Fig. 1, call it 0, and set φ 0 = θ 0 .Now let's move in a clockwise direction.Proceed with φ i = θ i until a negative link is found, in this case if a i−1,i < 0 make the change φ i = θ i + π, which changes the sign of the coupling.If the following links are positive keep making changes φ j = θ j + π but when getting into a new negative link a k−1,k < 0 restore the original angles φ k = θ k and proceed in this way along the ring until another negative link is found.If the number of negative links is even, this procedure can be repeated until arriving at i = N = 0 and φ 0 = θ 0 .The original equation has been converted into an all positive links KM for the new variables φ such that which means that depending on the number of negative links between oscillators 0 and i, the change of phase will be either 0 or π.However, when there is an odd number of negative links this gives rise to inconsistencies, since i = 0 cannot be changed and hence, in this case, the problem is equivalent to a unique negative link with the same change of variable provided in the previous equation, except for the last node (N ) that cannot be changed, and the last negative link is kept as being negative.
In conclusion, whichever is the number of negative links we have to consider only two possible cases (either 0 or 1 negative link), and in the next sections we will discuss how this affects the possible stationary solutions or the attractors of the dynamics.

III. FIXED POINTS
The fixed points correspond to the stationary solutions of the set of equations (1) at each node.All the equations are formally identical, and the only difference is whether a node is in between 2 positive links or between one positive and one negative link. 33

A. Local analysis
Let us define ∆ i = φ i − φ i−1 .The solution of the equation for node i depends on the sign of the links around that node.

Node with ++ links
In this case, we have Notice that solution a) involves the three phases, but solution b) involves only the extreme phases, but not φ i .

Node with +− links
In this case, we have which, again, has two solutions a) where the sign is taken to guarantee that the phase differences are in the interval (−π, π].b) which implies → φ i+1 − φ i−1 = 0. Solution a) involves again the three phases, but solution b) involves only the extreme phases, but not φ i .

B. Global analysis
The local analysis provides information about the relation between the phases of 3 consecutive nodes.Global constraints allow us to find the possible values for the phase differences.

All nodes have solution of type a)
We have a solution such that at each node, starting by i = 1, we assign a phase difference ∆ i = φ i − φ i−1 = α when the link between the nodes i−1 and i is positive and a phase difference ∆ i = φ i − φ i−1 = α − sign(α)π when this link is negative.Which values of α are possible?
The sum of phase differences along the circle has to add up to an integer multiple of 2π, this is the global constraint.If we call n the number of negative links, we have where q ∈ N is defined usually as the winding number (see 18,21 ).We can divide by π, arrange the expression, and define the lag number, ℓ, as which we will use as specification of the state of the dynamical system (attractor) from now on.In principle, q can take any integer value and the same for ℓ.Actually, ℓ accounts for the symmetry of the ring.When all links are positive ℓ is even, but when there is one negative link then it is odd.

All nodes have solution of type b)
Solutions of type b) are such that the phases of the neighboring nodes around a given one are related, and leave the phase at the reference node free.In this case, a solution of 2 alternate subnetworks of independent phases is possible.The two subnetworks are formed by alternate nodes such that • a i+1,i+2 = 1 and 0 otherwise.However, if N is an odd number this decomposition is not possible and hence a global solution of type b) is not possible.Let us assume by the moment that this solution is possible and wait for its stability in the next section.

Mixing of solutions of the two types
Let us assume first a global solution of type a) with a single solution for one node of type b) as shown in Fig. 2.

FIG. 2. Contribution of a single solution of type b) when all links are positive
In this case, we have α on the left and on the right imposed by the further nodes (labelled a)).Taking into account the type a) solution at both nodes it means that the two links have the same phase difference α, but these two phase differences have to sum π because of type b) solution, hence we conclude that α has to be π/2.We have only displayed the case of two consecutive positive links, but when we have negative links the conclusion that α = π/2 remains.A nontrivial extension to the case of two consecutive b) type solutions is shown in figure 3, also for positive links.In this case, it easily follows that between the two nodes with type b) solution the difference has to be β = π − α.On the other hand, when one of the links between the two b) type solutions is negative it implies that the phase difference is −α.In the next section we will discuss the stability of this type of fixed points.
Once the construction of the possible fixed points is finished, we have to look at the linear stability of the different solutions.

IV. LINEAR STABILITY OF THE FIXED POINTS
Let's make a small perturbation in the original equations (1) and linearize around the fixed point values of the phases: Linearizing for small {ε i } the previous equation, it allows us to write it in terms of the phase differences A. Type a) solutions Taking into account the value of the phase difference and how it depends on the sign we have Regardless of the sign of the links between node i and its neighbors, the previous equation ( 5) reads The set of equations can be rewritten in matrix form as Since the values of the elements are the same as for a totally positive cycle, we can make use of previous works that have shown that this matrix has a well known structure.Applying Gershgorin theorem (see appendix A) to our matrix, it is immediate to see that the eigenvalues are bounded in the interval [−4, 0] with a factor cos α multiplying it.It then means that when α < π/2 the fixed point is stable, and unstable otherwise.This restricts the available values of ℓ for the steady states to |ℓ| < N/2.At this point it is important to recall that in our description we have used ℓ as the lag number, which is more appropriate than the winding number widely used in the literature.For positive links both play the same role but for an odd number of negative links the stationary solutions have a slightly different interpretation.

B. 1 single b) solution
For the situation that is shown in Fig. 2 we arrived to the conclusion that the only possible value of α is π/2 which, according to the previous arguments, makes the eigenvalues to be zero.This holds for any number of isolated type b) solutions and hence these fixed points are not stable.

C. + link between two b) type solutions
On the other hand, the possibility of two consecutive nodes, as in Fig. 3, having a solution of type b) is still open.Let's look again in detail at the different contributions.
In this case, when i corresponds to the leftmost b) type node solution ∆ * i+1 = π − α and ∆ * i = α which, in the linearized equations, leads to When i is the rightmost b) type node then ∆ * i+1 = α and ∆ * i = π − α, which implies a global change in sign If the link between the two b) type solutions is negative, it is compensated by the change π − α → −α and the matrix elements are identical.
Summarizing, the matrix corresponding to the linearization around these compound solution is Applying Gershgorin theorem, in this case the eigenvalues lie in the range [−4, 2], also with a factor cos α meaning that there can be some positive eigenvalues making this configuration to be unstable.
Actually, this matrix can be considered to be a Laplacian matrix in a circle that corresponds to an adjacency matrix with all elements to be +1 except a single one that is -1.Corollary 3.7 of 34 states that the number of positive eigenvalues of a matrix like the previous one is equal to the number of negative links in the Laplacian interpretation.Hence, it shows that type b) solutions are always unstable.
If more solutions of type b) are considered then we can show that for an odd number of consecutive nodes having this type of solution the only possibility is the α = π/2 which makes again the solution to be metastable.For an even number of consecutive type b) solution what we have are boxes of the same type as discussed above and again apply the corollary to show that there are as many positive eigenvalues as boxes.

V. ATTRACTORS OF THE DYNAMICS
According to our findings in the previous section, we expect only attractors such that there is a constant phase difference between neighboring nodes, either α if the link is positive or α − π if the link is negative, being the values of α restricted to a discrete set α = πℓ/N and smaller than π/2 to guarantee the stability.Let us take into account as well that we have performed a change of angular variables ({θ i } → {φ i }) in such a way that the problem with an arbitrary number of negative links is reduced to a simpler situation with either 0 or 1 negative links.
In the first case, only positive links, the lag number ℓ = N α/π is equal to 2q, which means that it is even, and the restrictions of α limit the possible values of ℓ < N/2.On the other hand, for 1 negative link, ℓ = 2q + sign(ℓ) < N/2 and it is an odd number.In practice, we have two families of topologies and for each family we have different sets of attractors: for an even number of negative links only attractors characterized by an even lag number are possible, while for an odd number of negative links only attractors with an odd lag number are allowed.
We have performed simulations looking for the final stationary state of the evolution of populations of oscillators on a circle.Starting from random initial conditions, we look at the fraction of the initial states that end up in a particular final state (attractor), characterized by ℓ.In previous works 21 , the circle with only positive links has already been analyzed in detail.
Here, we generalize these results to a circle with an arbitrary number of negative links.We can corroborate by numerical simulations that an arbitrary even number of negative links is equivalent to the case where all the connections are positive.On the other hand, an arbitrary odd number of negative links is equivalent to the case with only one. 35In Fig. 4 the frequency of occurrence of a final state denoted by ℓ shows that, for different numbers of negative links n, the behavior of the system differs only when the parity of n changes.The simulations are performed by evolving 6 • 10 4 random initial conditions for each (N = 100, n) system.We solved the equations (1) with Runge-Kutta method implemented by the Julia library DifferentialEquations.jl 36and made use of other libraries for setting up the problem 37,38 .
Given this equivalence, we will restrict the discussion to the case n = {0, 1}.To characterize the basin of attraction of the Kuramoto ring, we now focus on the statistical properties of the final states.
To study how the distribution of attractors scales with N , we have computed, with the method described above, the evolution of 10 4 different initial conditions for each network configuration.For the symmetries of the system, we expect that E[ℓ] = 0.The interesting quantity is therefore the variance ℓ 2 (N,n) .While in Fig. 4, one cannot see a relation between P (ℓ) for n even and n odd, one can notice in Fig. 5 that for N > 8, remains the same regardless of the presence of negative links.Moreover, the variance scales linearly with the number of nodes N .We can, therefore, fit a linear model and find with κ = 0.12807±0.0001and the intercept is compatible with zero.Looking at even moments higher than two, we can see that they are compatible with a Gaussian distribution centered in zero (See fig.8 in appendix) in agreement with what was presented in 21 for q.We find, therefore, that for a fixed number of nodes, the attractor frequency is distributed along a Gaussian curve centred FIG. 4. The probability of the attractor characterized by the lag number ℓ, P (ℓ), is reported on the vertical axis.For a ring of N = 100 nodes and n negative links we computed 6 • 10 4 realizations.For the sake of clarity, we neglect bins with less than 10 data points.For n even, the only accessible attractors are the ones with ℓ even, and their frequency is independent of n (except statistical fluctuations).Conversely, n odd corresponds to ℓ odd and P (ℓ) is independent of n as well.For all configurations besides ℓ = 0 we have two almost overlapping points because, due to the symmetries of the system, P (ℓ) ∼ P (−ℓ).
in the synchronized state ℓ = 0 achievable only with even (and zero) number of negative nodes.The curve is the same for both even and odd numbers of negative links.In Fig. 6, we overlap the Gaussian function to the occurrence of the basins of attraction.κ is the coefficient defined in eq. ( 7), N the number of nodes and T is the size of the sample space, namely twice the number of realizations for each configuration to take into account the two available parities.It should be emphasized that in contrast to a Gaussian distribution, here, ℓ cannot assume any value within the set of real numbers; it must be an integer, its parity is determined by the parity of n and ℓ is strictly bounded by |ℓ| < N/2.

VI. CONCLUSIONS
In this study, we have delved into the dynamics of the KM on simple circle topologies, focusing particularly on the appearance of attractors different from the synchronized state.We have expanded the traditional KM by incorporating both excitatory and inhibitory interactions among oscillators.Our analytical findings suggest that the underlying network structure profoundly influences the system's collective behavior, generating a rich tapestry of stable solutions and their corresponding basins of attraction, which are intricately dependent on the nature of the inter-oscillator interactions.Interestingly, our work shows that even in a simple setting with identical oscillators, the introduction of inhibitory links creates a far more complex dynamical landscape than one might expect.These results challenge the conventional wisdom that assumes simplified dynamics in homogeneous systems and underscore the pivotal role of inhibitory interactions in shaping the global dynamics.Our findings open avenues for future studies aimed at understanding the consequences of these interactions in more complex networks and in systems with non-identical oscillators.
In conclusion, our study highlights the importance of considering both excitatory and inhibitory interactions when studying synchronization phenomena in networks.It invites further exploration of more complex topological structures and provides a foundation for investigating how various forms of local interactions collectively contribute to global dynamics.In the same line as the results in 21,39,40 , we have found that the basins of attraction of the different attractors form a very intertwined structure where borders are very difficult to analyze and transitions between different attractors will be very hard to predict.
In the future, our research on the dynamics of Kuramoto oscillators on a circle with positive and negative links can be extended to explore synchronization phenomena in more complex network structures, including multiplex networks 41 , higher-order interactions 40,42 , and temporal networks 43 , as suggested FIG. 6.The figure illustrates the occurrence of the basin of attraction instances characterized by the lag number ℓ using T/2 = 10 4 samples with random initial conditions.We have neglected again bins with less than 10 data points.The accessibility of the final states depends on the number of nodes N and the parity of the number of negative edges: if n even (n = 1) only even ℓ available, if n odd only, odd ℓ only are available.The solid line represents a Gaussian distribution centered at zero with variance σ2 = κN .κ is the only parameter whose numerical value is computed by the linear fit represented in Fig. 5 by previous studies.Investigating these directions holds the promise of uncovering novel synchronization patterns and behaviors, enhancing our understanding of complex systems, and enabling applications in various interdisciplinary fields, from communication networks to ecological systems, ultimately advancing our ability to model and control synchronization in real-world contexts.√ N and zero mean, we know all the higher moments.We can compute using (B1) the expected σ from higher orders: σ = p ⟨ℓ p ⟩ (p−1)!! if p is even .The plot analogous to Fig. 6, but for small values of N , illustrates that the Gaussian function inadequately approximates the distribution of the basins of attraction frequency in this case.In the left panel, the distribution is displayed on a linear scale, while in the right panel, it is represented on a logarithmic scale.In this way, we can emphasise the discrepancies in both the most and least populated states.

FIG. 1 .
FIG. 1.The model of oscillators in a ring geometry with bidirectional couplings.

FIG. 3 .
FIG. 3. Contribution of two consecutive solutions of type b) when all links are positive

FIG. 5 .
FIG.5.Second moment of the lag number ℓ distribution as a function of N .For small networks, the symmetries constraint the available values for ℓ 2 in case of odd or even number of negative edges.A linear fit for N ≥ 10 shows κ = 0.12807 ± 0.0001, with intercept compatible with zero.

FIG. 7 .
FIG.7.The occurrence of 10 4 realizations of the attractors as a function of N is depicted for all positive (left) and one negative (right) coupling edges.Labels corresponding to larger values of |ℓ| have been omitted for clarity.In this context, the area encompasses both positive and negative values of ℓ.Consequently, the probability P (ℓ = 0) may be lower than P (ℓ = |ℓ|) for ℓ > 0.

FIG. 8 .
FIG.8.If we expect the distribution to be a Gaussian, with σ proportional to √ N and zero mean, we know all the higher moments.We can compute using (B1) the expected σ from FIG. 9.The plot analogous to Fig.6, but for small values of N , illustrates that the Gaussian function inadequately approximates the distribution of the basins of attraction frequency in this case.In the left panel, the distribution is displayed on a linear scale, while in the right panel, it is represented on a logarithmic scale.In this way, we can emphasise the discrepancies in both the most and least populated states.