Emerging long orbits and self-similar temporal sequences in classical oscillators

2.7 MB·

by  Darka Labavić

Emerging long orbits and self-similar temporal sequences in classical oscillators
Darka Labavić and Hildegard Meyer-Ortmanns∗
Physics and Earth Sciences, Jacobs University, P.O. Box 750561, 28725 Bremen, Germany
We analyze repulsively coupled Kuramoto oscillators, which are exposed to a distribution of
naturalfrequencies. Thissourceofdisorderleadstoclosedorbitswithavarietyofdifferentperiods,
which can be orders of magnitude longer than periods of individual oscillators. By construction
theattractorspaceisquiterich. Thismaycauselongtransientsuntilthedeterministictrajectories
find their stationary orbits. The smaller the width of the distribution about the common natural
frequency is, the longer are the emerging time scales on average. Among the long-period orbits we
find self-similar sequences of temporary phase-locked motion on different time scales. The ratio of
time scales is determined by the ratio of widths of the distributions.
7
1
PACSnumbers: 05.45.-a,05.45.Xt,05.90.+m
0
2

n As model for classical oscillators we consider the Ku- size, geometry and coupling accordingly [11] and slightly
a ramoto model, a paradigm for synchronization phenom- perturb about the uniform natural frequency distribu-
J
ena with applications to natural and artificial systems tion. Later we compare the effects with cases of mod-
4
[1]. WhileinearliertimesmodificationsoftheKuramoto erate multistability. As we shall see, the time evolution
2
model were mostly considered for attractive (positive) of oscillator phases explores these attractor spaces quite

] couplingsbetweentheoscillators,nowadaysnegativecou- differently from a migration under the action of noise
O
plings are also discussed in the context of neural net- [11]. This paper will merely focus on a phenomenologi-
A works, representing inhibition [2], or social science, rep- cal description of the time evolution of trajectories with
. resenting contrarians [3]. When negative couplings come amazing features of long periods and self-similar tempo-
n
into the game, a notion of frustrated bonds can be de- ral patterns.
i
l fined, going back to [4] and later considered in [5]. As it We consider systems of N Kuramoto oscillators [12],
n
[ turnsout,frustrationinoscillatoryandexcitablesystems whose phases φi are governed by the equations:

can lead to analogous effects as in spin glasses. The first
1 onetomentionismultistablebehavior,wheretheattrac- dφi =ω + κ (cid:88)A sin(φ −φ ). (1)
v dt i N ij j i
9 tors are much more versatile in oscillatory systems than i j
8 in spin glasses. Also order-by-disorder phenomena can
Here ω denote the natural frequencies, κ parameterizes
8 be identified [6], where in replacement of thermal fluctu- i
6 the coupling strength, chosen to be κ = −2 through-
ations in spin glasses, disorder in oscillatory systems is
0 out the paper, that is deeply in the regime of multista-
realized as additive or multiplicative noise.
. bility. Further, N denotes the number of neighbors to
1 i
0 While first hints to a particularly rich attractor space which the i-th unit is connected, and Aij is the adja-
7 for certain grid sizes of repulsively coupled oscillators cency matrix with Aii = 0, Aij = 1 if i (cid:54)= j and units
1 were identified in [6], a closer exploration of the attrac- i and j are connected, otherwise Aij = 0. We choose
:
v tor spaces in [7, 8] showed what we would call a “dy- Aij to represent a hexagonal lattice of the size L × L
Xi namically generated Watanabe-Strogatz phenomenon". with periodic boundary conditions, so that Ni = 6 for
Startingfromaconsiderablefractionofinitialconditions all sites i. For the frequency distributions we choose a
r
a and uniform natural frequencies, the oscillator phases regular and a random implementation. For the random
arrange themselves into clusters of almost coalescing version we assume that the number density of natural
phases, where the clusters, considered as collective vari- frequenciesfollows,onaverage,theGaussiandistribution
ables, satisfy the criteria for dimensional reduction to g(ω)= √ 1 e−(ω−µ)2/(2σ2),describingfrequencyfluctu-
2πσ2
happenaspredictedbyWatanabeandStrogatz[9]. This ationsaboutµ=1withtypicallyσ =0.01. Fortheregu-
phenomenon goes along with N-3 conserved quantities lar implementation we generate frequencies according to
and continuous sets of solutions. This way the mecha- a deterministic procedure. The analytic expressions can
nism realizes a case of so-called extreme multistability be found in [13]. The choice is most easily illustrated
[10] with a continuum of attractors. Conditions, under in Fig. 1. The frequencies are uniformly spaced with re-
whichmoregenerallythisinterestingphenomenoncanbe specttog(ω)andgenerate“spatial"Gaussiansinnatural
observed, wereexploredbyM.Zaks[8]. Effectsbasedon frequencies with the largest frequencies in the center of
multistability should be particularly pronounced, when the grid (e.g., four sites in the middle of the 4×4 grid
such a continuum of attractors is available. Therefore of Fig. 1). This way the deterministic system retains
we shall choose our system of oscillators with respect to some symmetry in the frequency distribution over the
2
(a) g(ω) (b) (c) (d) (a)
0.005 random g(ω)
0.004
|Z|0.003
0.002
ω
0.001
ωmin μ ωmax 0.000 5000 10000 15000 20000
(b) t
0.0025
0.0020 regular g(ω)
Figure 1. Regular assignment of natural frequencies. Three Z|0.0015
|
differentsizesl ×l =4×4,5×4,5×5ofahexagonallattice. 0.0010
x y 0.0005
Colors represent a value of ω, orange being the largest, and
0.0000
0 5000 10000 15000 20000
purple the smallest. t
Figure 2. Transients of the order parameter for initial condi-
gridandcanbesmoothlytunedtowardstheuniformfre-
tions around the fixed point. Two trajectories in black and
quency. We shall analyze the dependence on the initial red are shown for a 4×4 hexagonal lattice for random ω’s,
conditions for φ , keeping the distribution of ω fixed. corresponding to two different realizations of g(ω) with the
i j
Toprojectontheinterestingfeaturesitturnsouttobe samemeanµandstandarddeviationσ,andregularω’s,cor-
responding to two different initial conditions, respectively.
quite convenient to study first the time evolution of the
All initial conditions are taken from a small radius around
Kuramoto order parameter, rather than the individual
ϕ = 0.5 ∀i. The parameters are µ = 1, σ = 0.01, and
i
trajectories of phase-correlated motion. It is defined as
κ=−2.
Z(t):=|Z|eiθ = 1 (cid:80)N eiφj, where a positive value of
N j=1
|Z|inthelimitofN →∞impliestheemergenceofphase
(a)
synchronizationandθisthephaseoftheglobalorderpa- 0.005 random g(ω)
0.004
rameter [12]. It is first |Z|, for which we observe striking |Z|0.003
0.002
longperiods,byordersofmagnitudelongerthanindivid- 0.001
0.000
ualperiodsofoscillators. Althoughidenticalvaluesof|Z| 60000 62000 64000 66000 68000
t
arecompatiblewith,butnotconclusiveforthesamepat- (b)0.0025 regular g(ω)
0.0020
treercnosrdofanpdharseeg-ilsotcekredthemseotlioonng, ipteirsioadsu.seLfuonlgfirpsetrtioodosl toof |Z|00..00001105
0.0005
the trajectories may be easily overlooked when following 0.0000
2000 22000 24000 26000 28000 30000
the individual time evolution of larger sets of oscillators. t
Asanextstepwethenfocusonselectedpointsoftheor-
der parameter and analyze the trajectories of individual Figure3. Periodictimeevolutionoftheorderparameterwith
oscillators to identify their period and synchronization initialconditionsaroundthefixedpoint,forthesametrajec-
pattern. As we shall see, apart from constant values for tories of Fig 2.
|Z|, corresponding to stable limit cycles with constant
frequencies and periods of the order of individual cycles,
we observe long transients to reach a stable orbit and monostableregime(κ>0, whichiscertainlyfarofffrom
long periods of these orbits, for both types of frequency thefinallyapproachedphasepattern), orfromauniform
distributions and both kind of initial conditions. As ini- distribution between 0 and 2π, so that the long duration
tial conditions we either start from values close to the depends on the attractor space itself, as we shall later
fixed point (in phase differences, which is a unique solu- confirm by a comparison to results for the 5×5-lattice.
tionforpositiveκandfaroffthedesiredphasepatternfor The following Fig. 3 displays the periodic time evolution
negativeκ)(asweusedforFig.2),orinstead,fromauni- of|Z|onceithasstabilizedafterthetransients. Depend-
formdistributionintheinterval(0,2π),displayedin[13]. ing on the realization of ωi’s and/or initial conditions,
All presented order parameters are periodic orbits with the system evolves to a different attractor. Attractors
differentmaximaanddifferentperiods. Transientsareof here are characterized by means of the time evolution
thesameorderofmagnitude,rangingfrom∼10000(bot- of|Z|thatindicatesdifferenttemporarysynchronization
tomgray(onlinered)trajectoryinFig.2,randomassign- patterns during the orbit, before we zoom into individ-
ment, to ∼17000 t.u. (bottom black curve, Fig. 2, regu- ual trajectories. All four order parameters in Fig. 3 have
lar assignment) (and similarly from ∼10000 to ∼16000 different periodic orbits.
for the second type of initial conditions, shown in [13]. The time evolution of the |Z| with initial conditions
ForFig.2wepursuedtheperiodicbehavioruntil500000 from ϕ ∈ (0,2π) show periods ranging from 1000 to
i
time units (t.u.). 10000 for both random and regular g(ω).
Unless the choice of initial conditions would acciden- In order to get hints on how representative the long-
tally hit the vicinity of an attractor for κ<0, transients period orbits and the long transients are for the systems
can take rather long, independently of whether starting we consider, we display in Fig. 4 histograms of the order
from the neighborhood of a distinguished value in the parameter periods for both lattice sizes and frequency
3
(a)100 4x4 random (b)400 12 5x5 random form cases considered there. From the results of [7, 8]
80 300 180 a clear difference between the attractor spaces for these
60 200 46 two lattice sizes can be traced back to a larger number
40 2
20 100 00 2000400060008000 of possible 4-cluster partitions for the 4×4 case, which
period[t.u.] can be assigned to the grid as patterns of phase-locked
0 0
0 50 100150200250300 0 2000 4000 6000 8000
period[t.u.×102] period[t.u.] motion. The lattice allows the phases to organize them-
(c)400 200 (d)400 selvesintofouractuallygloballycoupledclusters,eachof
4x4 regular 5x5 regular
300 150 300 150 four coalescing phases. It is this system of emerging four
100 100
200 50 200 50 globally coupled clusters that satisfies the conditions for
100 00.00.51.01.52.02.5 100 00 2 4 6 8 dimensional reduction à la Watanabe and Strogatz [9],
|Zconst×10-3 |Zconst×10-3 therefore leading to a continuum of attractors, with one
0 0
0 1000 2000 3000 4000 0 500 1000 1500
conserved quantity (see below), and a continuum of so-
period[t.u.] period[t.u.]
lutions differing by frequencies. Depending on the width
of our frequency distributions about µ, we violate the
Figure 4. Histograms of order parameter periods as a func-
conditions for the occurrence of the “Watanabe-Strogatz
tion of the period length for 500 realizations of g(ω) for ran-
phenomenon" more or less slightly, so we expect to see
domω’s,and500differentinitialconditionsforregularω’sfor
two different lattice sizes. Parameters are µ = 1, σ = 0.01, “remnants" of the different attractor spaces at uniform
κ = −2, initial conditions are taken from a vicinity of the frequencies for different widths, when we turn on devi-
fixed point at ϕi =0.5. ations from the uniform natural frequency. Before we
discuss the dependence of the transient times and peri-
ods on the strength of the deviations, which reveal some
implementations. An order parameter with zero period remnants, let us take a closer look at the long-period
corresponds to a constant order parameter. We see a orbits.
clear difference between the two lattice sizes for random When assigning phases to the unit cycle, each lattice
distribution of ω’s. In both cases there is a clear peak site is characterized by a gray scale (color online) and a
at zero, for the 4×4 lattice 17% of the realizations (for size ofa disk, so thatin principle we can disentangle six-
a fixed initial condition) evolves to an attractor with a teen even coalescing phases on the unit cycle as a sphere
constant order parameter, compared to 80% in the case with 16 non-overlapping internal rings, see [13].
for the 5×5 lattice. Beside this larger number of zero- In order to understand what the long-period orbits of
period order parameter orbits (that is, limit cycles), the theKuramotoorderparametermeanforthetrajectories,
periods are also longer for 4×4 by an order of magni- we show the time evolution of |Z| for one choice of ini-
tude. Similarly to the 5×5-grid, also the 4×5-grid (not tial conditions for regular ω’s. For selected instants of
displayed) has a large number of zero-period orbits, but time, indicated by the gray (online red) points in Fig. 5,
a few long-period orbits as well, whose periods are in we show a snapshot of the sixteen phases together with
general shorter than in the 4×4 case, but still orders of their representation on the unit cycle. While the time
magnitude longer than the individual oscillation period. evolution of the phases indicate a phase period of about
In the case of regular ω’s, the difference is less pro- 6 t.u. ≈ 2π/ω , the total period of the periodic orbit
i
nounced for 4×4 and 5×5 lattices: Most of the initial of ≈ 1500 t.u. is about 100 times longer in this case.
conditions lead to a constant order parameter, and we Moreover, in the concrete case of Fig. 5, only every sec-
findonlyafewdifferentperiodicorbits. Forazero-period ond absolute maximum corresponds to the same state,
order parameter, i.e. a constant, for the regular assign- because two neighboring absolute maxima correspond to
ment of ω’s, we show in the insets of (c) and (d) how differentphasepatterns, whichadduptothesamevalue
many different constant order parameters we find, corre- of |Z|. For all long-period orbits in our figures we have
spondingtodifferentlimitcyclesfordifferentinitialcon- checked when the return of the value of |Z| meant the
ditions,andindicatingmultistability. Forrandomassign- samepatternofphase-lockedmotion,sothatthelongpe-
ments, different order parameters, recorded for different riods of |Z| are actually representative for long periods
realizations, are not conclusive in view of multistability, of the system. As such an example the bottom row of
as the systems differ by their natural frequencies. So far Fig. 5 shows two (almost) the same states of the system,
we changed initial conditions only for one random real- a three cluster state (with three bubbles along the unit
ization and found bistable behavior: two different order cycle). Itisnumericallyhardtoexactlyhittheinstantof
parameters, one constant, and one periodic. time, when states like the three-cluster state are exactly
Inthiskindofsystemsacompleteoverviewofthepos- the same, as they are temporary patterns. In the movie
sible states is beyond the scope of this paper. More in- (see [13]) we show the time evolution over one whole pe-
sight into the versatility of multistability was obtained riod, that is 272 seconds, corresponding to 1355.8 t.u.
for uniform frequency distributions in [7, 8], and our or 13558 frames in the movie, so that the periodicity in
frequency distributions are perturbations about the uni- configurations of phases becomes visible if one is patient
4
(a) (b) 20000 100000
15000 81250
(c) Period 105000000 4632755000 Transient
0 25000
0.001 0.002 0.003 0.004 0.005 0.006
σ
(d) (e)
Figure6. Transienttimeandperiodasafunctionofthewidth
σ oftheGaussiandistributiong(ω)forregularassignmentof
natural frequencies, fixed µ = 1, κ = −2 and identical ini-
tialconditions. Periods(reddots)andtransienttimes(black
dots) are decreasing with an increase of the Gaussian width.
Figure 5. Time evolution of the phases and order parameter,
phase snapshots on the unit cycle for one choice of initial
conditions for regular ω’s. The parameters are µ = 1, σ =
conditions for regular frequency assignment and change
0.01, κ = −2 as before with initial conditions taken from a
only the standard deviation σ from 0.001 to 0.006, while
vicinity of the fixed point at ϕ =0.5.
i
the other parameters are kept fixed. Even starting with
the same initial condition does not guarantee that the
to wait until the end. system will evolve to the same type of attractor for dif-
ferent σ, just in a different transient time. To consider
A closer look at the evolution of phase patterns over
trajectories, which are comparable in both the transient
a long time such as a whole period reveals approximate
time and the “attractor type" itself, we make sure that
temporaryn-clusterstates,wherenvariesintimefrom2
thesystemactuallyevolvestothesametypeinthesense
to16. Fromourformeranalysis[11]foruniformfrequen-
that the periodic orbit has the same shape as illustrated
cies we know that a 3-cluster state was unstable, while
in Fig. 7, and compare time scales only for such trajec-
4-or 6-clusters states were stable as long as noise was
tories. In conclusion to Fig. 6, it confirms our conjecture
absent to make these states metastable and kick them
that the smaller the width, so the richer the space of at-
out of their basin of attraction. The dwell times close
tractors for the 4×4 grid is, the longer are the transient
to the approximate n-cluster states in the present sim-
times and periods.
ulations are not much longer than time intervals, over
which phase patterns cannot be uniquely classified in In Fig. 7 we show order parameter trajectories for the
terms of n-cluster states due to the spreading between largest (σ = 0.006) and smallest (σ = 0.001) Gaussian
the phases. So far our observations are compatible with widths, considered before in Fig. 6, to illustrate what we
the interpretation that the long-period orbits are hete- meanbythesametypeofattractor. Itshouldbenoticed
roclinic cycles, connecting unstable limit cycles. We do that the very fact is quite striking to observe this self-
not observe longer and longer dwell times when the tra- similarity of trajectories over time scales, which roughly
jectory comes close to a certain temporary pattern when differ by a factor of 6 (the same factor as given by the
repeating the cycles; longer dwell times were expected, ratio of widths). While the period is about 12500 t.u.
if these patterns would correspond to saddle points con- for σ = 0.001, it is about 2000 t.u. for σ = 0.006. The
nected by a heteroclinic cycle. As the phase space is same scaling of ratios of periods with ratios of widths
sixteen-dimensional, an analytical proof of this interpre- is observed for other widths as well, see [13]. A simi-
tation is out of the scope of our phenomenological de- lar shape of the trajectories means a similar sequence of
scription. However, in forthcoming work we will analyze synchronization patterns. Yet, the absolute values of |Z|
thestabilityofthe long-period orbitsunder variationsin differ at comparable peaks. A zoom into the actual syn-
the initial conditions and under the action of noise. chronization pattern (not displayed) at time instants of
To support our hypothesis that it is the rich attractor correspondingmaximashowedthatthetemporarystates
space,whichcauseslongtransientstowardsfindingasta- could indeed be classified by n-cluster states, for ex-
ble trajectory and a long time to close an orbit (possibly ample three clusters, but the tolerance interval, within
heteroclinic) between unstable limit cycles, we analyze which the phases coalesce within a cluster, varies be-
in Fig. 6 the dependence of transient time and period on tween both, the widths, and the individual phases. So
the width of the Gaussian distribution for a 4×4-grid: the self-similarity refers to temporal sequences of syn-
a tiny width means an almost uniform frequency distri- chronization patterns within adapted tolerance intervals
bution with a known rather rich attractor space, which of what makes up a cluster.
shouldgetreducedforlargerdeviations. Tohavecompa- A third type of long time scale is generated if we
rable starting configurations, we choose the same initial startwithinitialconditionsofafour-clusterstate, which
5
(a)0.0008 σ = 0.001 taining k and m collide, or l and n, I becomes infinite.
0.0006
Again,inagreementwithFig.7,thesmallerthewidthof
|Z|0.0004
0.0002 the distribution about the uniform case, the longer the
0.0000 transient of the system to reach the attractor.
400000 410000 420000 430000 440000 450000
t
(b) As we shall show in forthcoming work, the variety of
0.0025 σ = 0.006
0.0020 transienttimes,inparticularlongtransientsleadtophys-
|Z|00..00001105 icalagingofourdeterministicoscillatorysystems,aslong
0.0005
0.0000 as both time instants of the autocorrelation measure-
402000 404000 406000 408000
t ments fall into the transient period, so that stochastic
fluctuations are not necessary for aging.
Basic ingredients for our observations were antagonis-
Figure7. Self-similartrajectoriesoforderparametersintime
tic couplings, leading to a rich attractor space, together
for two different widths of g(ω). Notice the different time
scales. Theoveralltimeintervalissome40000t.u. shorterin with a distribution of natural frequencies and an assign-
the bottom figure. ment to a regular grid. Therefore we expect that our re-
sults are not restricted to our choice of model, but point
0.4 σ = 0.01 σ = 0.005 σ = 0.001 to a mechanism, which may be more generally at work
0.2
0.0 before after in biological systems, where these basic ingredients are
-0.2
I-0.4 frequently found.
-0.6
-0.8
-1.0 1 10 100 1000 104 105 Financial support from Deutsche Forschungsgemein-
t
schaft (DFG, contract ME-1332/25-1) is gratefully ac-
knowledged.
Figure 8. Decay of the formerly (ω = ω) stable 4-cluster
i
state, measured via the conserved quantity I as defined
klmn
inthetext. ThethickgraylineshowsthevalueofI for(σ =
0), black lines for σ = 0.001 (full), 0.005 (dashed), and 0.01
(dotted). The two insets show typical states of the system
∗ [e-mail]
before and after the collision of the clusters. The time axis
[1] J.A.Acebrón,L.L.Bonilla,C.J.Pérez-Vicente,F.Ritort,
is represented in a logarithmic scale to make the figure more
R. Spigler, Rev. Mod. Phys.77,137 (2005).
readable.
[2] E. Ledoux, N. Brunel, Front. Comput. Neurosci. 5 (25),
1(2011).
[3] H. Hong, S.H. Strogatz, Phys. Rev. E 84, 046202(2011).
satisfies the criteria for the Watanabe-Strogatz phe-
[4] H. Daido, Phys.Rev.Lett.68,1073(1992).
nomenon and would be stable as long as σ = 0. The
[5] P. Kaluza, H. Meyer-Ortmanns, Chaos 20, 043111
long decay time of the conserved quantity I =
klmn (2010).
sinϕk−2ϕl sinϕm−2ϕn as described in [9] and predicted to [6] F.Ionita,D.Labavic´,M.Zaks,H.Meyer-Ortmanns,Eur.
sinϕk−ϕm sinϕl−ϕn Phys. J. B.86(12), 511(2013).
2 2
be constant for ω = ω in the four-cluster state, slowly
i [7] P. Tomov,Interplay of Dynamics and Network Topology
decays. In particular we show in Fig. 8 how long it takes inSystemsofExcitableElements,Ph.Dthesis,Humboldt
the system to leave the 4-cluster state. The decay time Universität zu Berlin, (2016).
depends on the width of the Gaussian distribution for [8] M. Zaks, private communication.
regular ω’s. When the value of I(σ) starts to diverge [9] S. Watanabe and S. H. Strogatz, Phys. Rev. Lett. 70,
2391 (1993); Physica D 74, 197 (1994).
from a constant I(0), clusters move towards each other
[10] C.Hens,S.K.Dana,U.Feudel,Chaos25,053112(2015).
(whilefor σ =0thecluster differenceremainsconstant).
[11] F. Ionita, H. Meyer-Ortmanns, Phys. Rev. Lett.112,
Oncetheclusterscollide,thesystemevolvestooneofthe
094101 (2014).
available attractors. If the clusters containing oscillators [12] Y. Kuramoto, Chemical Oscillators, Waves, and Turbu-
k and l collide, or m and n, I becomes zero (right inset), lence, (Springer, New York, 1984).
definingthedecaytime. Alternatively,iftheclusterscon- [13] See Supplemental Material at [...] for further figures.