Dynamic Behavior of a Large Ring of Coupled Self-Excited Oscillators

The collective dynamics of a complex system consisting of a large number of self-excited oscillators interacting with each other in some way is very interesting. Self-excited oscillators occur in diverse fields such as biology [1–3], optics [4], electronics [5–7], fluid flow [8,9], and buildings [10]. The global behavior of the system is dependent on, but different from, that of the individual oscillators, and new emergent dynamic phenomena arise when they are coupled. Frequently, however, the coupling is ignored, and they are studied in isolation so that the global dynamics is not appreciated.

where $Dosc$ is a time-independent, nonlinear, differential operator that leads to self-excited oscillation $y(t)$⁠, where $t$ is time. The emphasis here is on a lack of forcing by an external, periodic input, which would provide a frequency; the self-excited oscillator has no need of that. Examples of such operators are integrate-and-fire [11], hysteretic thermostat [12,13], van der Pol [14], and many others. The periodic variable $y(t)$ is a problem-dependent quantity that may be a current, displacement, light intensity, etc., but which will simply be referred to as motion here. One of the simplest connecting configurations possible is that of a ring [15–18]. Let us look at $N$ of these oscillators in this form, as shown in Fig. 1, and let each be connected to its immediate neighbor on either side. The numbering system of the ring geometry means that oscillators $i=0$ and $i=N+1$ are actually $i=N$ and $i=1$⁠, respectively. Then

with $i=1,2,…,N$⁠, where the function $f$ takes into account the interaction of oscillator $i$ with its neighbors $i-1$ and $i+1$⁠. If all solutions of Eq. (2) are assumed to be the same, i.e., $yi=y(t)$⁠, then it can be seen by substitution that it is a solution, as long as it is also assumed that $f(y,y,y)=0$⁠. All the oscillators move in phase, there is zero interaction, and the frequency is the same as that of a single oscillator given by Eq. (1).

for $i=1,2,…,N$⁠, with $a>0$ and appropriate initial conditions.

Pioneering work on this problem was done by Lubkin and Rand [1], Linkens [15], and Endo and Mori [16]. The present authors have worked on a ring of four oscillators [14,21]. In that context, they have shown that if $b<-0.25$ the oscillators are unstable [21]. Phase shift depends on the interactive coupling; as the coupling parameter $b$ is raised, the phase shift decreases and disappears beyond a certain finite value. Others have also worked on different aspects of a ring of four coupled van der Pol oscillators [22–24].

This work will build upon previous work on interacting van der Pol rings and extend them in several directions: studying the effect of singularity of the coupling matrix on the ring dynamics, consideration of a large number of oscillators, and highly nonlinear behavior of the oscillator in the form of a large $a$ parameter. The existence and linear stability of time-independent solutions is studied in detail. Also investigated are the emergent behavior for the ring as a whole, the waves that occur, their dispersion relationship, and phase and group velocities.

These are the simplest solutions to the problem; their existence and linear stability are investigated.

where $I$ is the $N×N$ identity matrix, and $B$ is given by Eq. (7b).

and so forth. In addition, $σ2N=detA$⁠.

Polynomial stability requires that all the roots have negative real parts. The necessary and sufficient condition for this is that all coefficients of the first column of the Routh–Hurwitz (RH) array have the same sign. Table 1 shows the derived expressions for four coefficients of the RH array corresponding to the characteristic polynomial of Eq. (13). The first element of the RH array is positive whereas the second one is negative. The third coefficient is always positive for any value of $N$⁠. The fourth coefficient is negative for $N≥3$⁠. Two sign changes are detected for these four coefficients; consequently, the characteristic polynomial of $A$ has at least two right-half plane roots. This means that the RH stability condition is not satisfied and, therefore, the trivial solution is unstable.

for $j=1,2,…,N-1$⁠. Some of these critical values of $b$ are given in Table 2, where only values that are different from each other are shown. It can be seen that for even values of $N$⁠, there are $N/2$ different critical values, while for odd there are $(N-1)/2$⁠.

An interesting behavior of coupled oscillators is oscillation death, which occurs through a saddle-node bifurcation mechanism allowing the emergence of new fixed points. A new stable steady state with nonzero amplitude is created by the coupling [25,26]. When $B$ is singular, Eq. (6) can have an infinite number of steady-state solutions, and oscillation death can arise. The right column of Table 2 shows the $rank(B)$ as a function of $N$ and critical $b$⁠. For even values of $N$⁠, $rank(B)=N-1$ for $b=-0.25$ and $rank(B)=N-2$ for other values of $b$⁠. For odd values of $N$⁠, $rank(B)=N-2$⁠, irrespective of the value of $b$⁠. Then, when $B$ is singular, this matrix has one or two linearly dependent rows or columns, which can be ignored, and the dynamical system can be analyzed in a reduced dimensional space.

The individual waves travel at the phase velocity, while the information travels at the group velocity. The phase velocity decreases with the wave number $k$⁠, but the group velocity increases. Furthermore, as $N→∞$⁠, $νp→1/k$ and $νg→0$ if $b$ is held constant. The reason for this is that the angular distance between oscillators $2π/N→0$⁠, and hence the interaction term in Eq. (18)$4π2b/N2→0$ also. If, however, it is assumed that $b~N2$ such that $4π2b/N2=C$⁠, where $C$ is a constant, then both $νp=k-1(1+Ck2)1/2$ and $νg=Ck(1+Ck2)-1/2$ are independent of $N$⁠. They become independent of the wave number $k$ also for large $k$⁠, and strangely enough, become equal.

Equation (5) is solved numerically using a fourth-order Runge–Kutta method. A nonzero $a$ creates self-excitation, and no external forcing is needed, but at the same time, a small enough $a$ satisfies the linearized analysis. Numerical solutions shown are for $a=0.01$ and an arbitrarily chosen $b=1$⁠. The time step used for integration is $10-4$⁠, which is found to be small enough for numerical stability and convergence. It is observed that $N=100$ is sufficiently large and that any increase in the number of oscillators does not change the global behavior of the system.

The simplest version of Eq. (19) is with zero wave number. The motion is then independent of space and only oscillates in time. Full synchronization in which all oscillators move in phase is exhibited by the ring. Equation (20) shows that the frequency does not depend on $N$⁠.

The temporal motions of five equally-spaced oscillators around the ring are shown in Fig. 2(a) for $k=1$⁠. The phases are seen to be equally spaced also. Figure 2(b) shows three oscillators that are neighbors. The wave peaks first for $i=40$⁠, then $i=41$⁠, and finally $i=42$⁠, indicating a wave motion around the ring in a counterclockwise direction. In the biological literature, especially for cilia, this is often referred to as a metachronal wave [27–29]. Similarly, there is another solution possible, $k=-1$⁠, in which the $(i+1)th$ oscillator leads the $ith$ so that it would be an identical wave but in the clockwise direction. A standing wave is formed by identical waves running in opposite directions. Other possibilities are those of multiple integer or noninteger wavelengths in the ring depending on $k$⁠.

Van der Pol oscillators with $a≪1$ of Sec. 3 present quasi-linear behavior and, therefore, exhibit oscillations that are nearly sinusoidal. Nonlinear behavior is more evident for $a≥1$⁠. For single van der Pol oscillators, it is known that the frequency of oscillation decreases with increasing $a$⁠. The form of the oscillation also changes, the rise and fall in each period becoming steeper as $a$ increases. The numerical results shown here are obtained in the following way. Since initial conditions have a significant effect on the long-time results, the small $a$ initial conditions are introduced as before. Then, however, the value of $a$ is ramped up slowly with each time step until the desired value is reached, at which stage it is held constant for 1000 time units.

Figure 3(a) depicts a time series for $N=100$⁠, with initial conditions corresponding to Eq. (19) with $k=2$⁠. A modulated signal is clearly visible with a dominant carrier frequency and a modulation, which persists even if the integration is carried out much longer. The power spectral density, shown in Fig. 3(b), quantifies these frequencies as 0.1527 and 0.0076, respectively. The envelope is shown in further detail in Fig. 3(c) for neighboring oscillators. The envelope is itself a traveling wave that is moving at the group velocity.

Figures 4(a) and 4(b) show a time series and its corresponding power spectrum for $N=100$ and initial conditions corresponding to $k=2$⁠. Three frequencies of 0.0534, 0.1607, and 0.2680 are detected. The first is the fundamental, the others being the third and the fifth harmonics, respectively.

Self-excited oscillators connected to each other in the form of a ring with linear neighbor-to-neighbor coupling have been considered. Global behavior of the ring was investigated by using the example of a van der Pol operator. The effect of the singularity of the matrix $B$ in the ring dynamics has been explored. When this matrix becomes singular, the ring has an infinite number of steady states and oscillation death can arise. However, fixing the motion of some oscillators makes the system determinate and oscillating.

For $b>0$ and small $a$⁠, it has been shown that wavelike behavior exists for different initial conditions. These dispersive, metachronal waves traveling around the ring can be analytically determined and numerically confirmed. The dynamics for large $a$ are more complicated; nonlinear behaviors, such as modulation and harmonics, emerge.

It has been shown that interconnected oscillators can propagate signals. An oscillator has an intrinsic frequency, i.e., a stand-alone frequency if it were to be completely isolated. For a van der Pol oscillator, the intrinsic linear frequency for small $a$⁠, determined by the coefficients of the $y··$ and $y$ terms, has been taken to be unity here. However, as $a$ becomes larger the intrinsic nonlinear frequency is somewhat smaller. For communication purposes, this is a carrier frequency. Interaction between oscillators, however, is determined by the coupling constant $b$⁠. For nonzero $b$⁠, information travels at group velocity along the ring in the form of waves as an amplitude modulation of the carrier.

Nomenclature

$a$ =

damping parameter

$aij$ =

elements of $A$

$A$ =

matrix of differential system

$b$ =

coupling parameter

$B$ =

coupling matrix

$D$ =

$det(B)$

$Dosc$ =

oscillator operator

$f$ =

coupling function between oscillators

$I$ =

identity matrix

$k$ =

wavenumber

$N$ =

number of oscillators

$t$ =

time

$νg$ =

group velocity

$νp$ =

phase velocity

$y$ =

dependent variable

$y¯$ =

time-independent value of $y$

$y$ =

column vector of $y$'s and $z$'s

$z$ =

$y·$

Greek Symbols

$θ$ =

angular position of oscillator

$λ$ =

eigenvalues

$σ$ =

coefficients of characteristic polynomial

$ω$ =

frequency

Subscript

$i$ =

oscillator number