Chimera states in ensembles of excitable FitzHugh-Nagumo systems
EEPJ manuscript No. (will be inserted by the editor)
Chimera states in ensembles of excitableFitzHugh–Nagumo systems
Nadezhda Semenova , , a D´epartement d’Optique P. M. Duffieux, Institut FEMTO-ST, Universit´e Bourgogne-Franche-Comt´e CNRS UMR 6174, Besan¸con, France. Department of Physics, Saratov State University, Astrakhanskaya str. 83, 410012 Saratov,Russia.
Abstract.
An ensemble of nonlocally coupled excitable FitzHugh–Nagumosystems is studied. In the presence of noise the explored system canexhibit a special kind of chimera states called coherence-resonancechimera. As previously thought, noise plays principal role in formingthese structures. It is shown in the present paper that these regimes ap-pear because of the specific coupling between the elements. The actionof coupling involve a spatial wave regime, which occurs in ensembleof excitable nodes even if the noise is switched off. In addition, a newchimera state is obtained in an excitable regime. It is shown that thenoise makes this chimera more stable near an Andronov–Hopf bifurca-tion.
All real systems are inevitably affected by noise. The impact of noise plays a principalrole and strongly dictates the properties of the oscillatory dynamics. Noise influenceslead to deterioration of the observed effects or to their complete destruction. However,this does not always lead only to destructive effects. It is known that noise canqualitatively change the oscillatory behaviour, and induce bifurcation phenomenawhich give rise to the appearance of new regimes of functioning [1–7]. Moreover, suchtransitions can be accompanied by increasing of regularity with growth of the noiseintensity. Such cases may include, for example, coherence or stochastic resonance [2,8]. However, the affected system must have some special characteristics and propertiesto demonstrate new effect in the presence of noise.Noise induced effects in single oscillators are already quite well understood. Cur-rently, most of open questions address the noise impact in networks and ensembles[9–12]. Ensembles of identical and non-identical oscillators initially have a wide rangeof possible regimes. Thus, apart from typical synchronization, desynchronization, spa-tial incoherence or coherence, such systems can simultaneously demonstrate the co-existence of several modes. An example of such coexistence are chimera states. Itis a spatio-temporal pattern consisting of the neighboring clusters of elements withcoherent and incoherent dynamics in one network. Initially, this effect was discov-ered by Kuramoto and Battogtoh [13], and then revealed in more details by Strogatz a e-mail: [email protected] a r X i v : . [ n li n . AO ] N ov Will be inserted by the editor and Abrams [14, 15]. Since then, quite a long time has passed, and now this effecthas already been found in the ensembles of many systems of various nature: phaseoscillators, chaotic mappings, bistable elements, neuronal models, and many others[12]. Usually one of the main conditions imposed on the ensemble is non-local connec-tion, when each partial element has a finite radius of influence. Nowadays, chimerastates have already been found in networks of various topologies [16–18], as well asin multilayer networks [19, 20]The ensemble of nonlocally coupled neuron models is considered in the present pa-per. The FitzHugh–Nagumo system in excitable mode is chosen as a partial element.Initially, chimera states have already been found for an ensemble of similar systemsbut in oscillatory regime [21]. Later, it has been shown that chimeras can be detectedin the excitable mode in the presence of noise [22]. This effect was called coherence-resonance chimeras (CR chimera) because it combined the properties of chimeras andcoherence resonance and was characterized by periodic switching of the coherent andincoherent parts position. Another feature was that CR chimera appeared only whenthe noise intensity belonged to a certain interval. If it was smaller, there were nooscillations in the system, but too strong noise led to complete spatial incoherenceand desynchronization.This work is devoted to the study of coupling features being principal for realiza-tion of chimera states in ensembles of excitable oscillators. The comparative analysisof noise role and coupling peculiarities is carried out for the example of an ensembleof nonlocally coupled FitzHugh–Nagumo systems. In order to find out which cou-pling features lead to such a vulnerability of the system to noise, and whether otherspatio-temporal regimes can be obtained taking these features into account.
In the present paper the dynamics of a one-dimensional ensemble of nonlocally coupledFitzHugh–Nagumo systems is studied numerically. The ensemble is described by thefollowing system of equations: ε du i dt = u i − u i − v i + σ R i + P (cid:80) j = i − P [ b uu ( u j − u i ) + b uv ( v j − v i )] , dv i dt = u i + a + σ R i + R (cid:80) j = i − R [ b vu ( u j − u i ) + b vv ( v j − v i )] + √ Dξ i ( t ) , (1)where u i and v i are activator and inhibitor variables of i th oscillator. The numberof all nodes is N . The dynamics of partial elements depends on two parameters: ε and a . The first one controls the time scale in the system and a is the thresholdparameter. Depending on its values, the individual FHN element can demonstratean oscillatory ( | a | <
1) or excitable ( a ¿1) regime. In oscillatory regime the systemexhibits the periodic behaviour associated with the limit cycle in the phase plane.In the excitable regime all oscillations decay and there is only one stable fixed pointin the phase plane. But this regime is interesting because of a coherence resonance.This effect consists in stochastically excited spiking, which is the most regular at acertain noise intensity range [2]. The system (1) contains white Gaussian noise source ξ i ( t ) ∈ R of the intensity D . Index i indicates that the noise sources are not correlatedinside the network, and each oscillator has only its own noise source. All of them havethe same statistical characteristics and intensities.The components of Eq. 1 with sums describe the connectivity. Here a parameter σ is the coupling strength. P is the number of elements connected with each i th oscillatoron the either side. Normalizing this value by the total number of elements, we get the ill be inserted by the editor 3 coupling radius r = P/N . Under the sign of the sum there are two summands in eachof the equations. In the first equation the parameter b uu defines the contribution of u -variables of the connected neighbours. A coefficient b uv is the same for the variables v in the first equation. Parameters b vu and b vv in the second equation have the samemeaning. Thus, the parameters b uv and b vu are responsible for the cross-linking. Thetype of connection in Eq. 1 came from neuroscience. It is described in more detailin [23–25]. Taking into account the large number of parameters, it is reasonable toenter one main operator controlling all four b -parameters. To do this, the rotationalcoupling matrix is introduced: B = (cid:18) b uu b uv b vu b vv (cid:19) = (cid:18) cos φ sin φ − sin φ cos φ (cid:19) . (2)Now there is only one parameter φ , which controls the impact of cross-linking( b uv , b vu ) and self-linking ( b uu , b vv ) in the equation. Chimera states in the ensembles with nonlocal coupling have been found for thesame values of φ = π/ − . σ = 0 . r = 0 . r = 0 .
12 there are two incoherent domainsin the network, and at r = 0 .
08 there are three of them. The chimera state exists in aquite wide range of σ parameter values, but at σ < . Z k = (cid:12)(cid:12)(cid:12) δ Z (cid:88) | j − k |≤ δ Z e iΘ j (cid:12)(cid:12)(cid:12) , k = 1 , . . . N (3)where the geometric phase of the j th element is defined by Θ j = arctan ( v j /u j ) [21].The values Z k = 1 and Z k < δ Z is fixed δ Z = 20.Figure 1,c shows the instantaneous spatial profiles (snapshots) of wave profiles indifferent half-periods. This panel clearly shows that the position of the incoherentpart is first located on the edges of the ensemble, and then switches to the middle.Despite the seeming instability of the state, it is saved even at t = 10 .It was obtained that for the value of φ = π/ − . . ≤ a ≤ . D ∈ [0 . . φ = π/ − . Will be inserted by the editor
Fig. 1.
Coherence resonance chimera with corresponding space-time plot for the variable u ti (a), local order parameter Z ti (b) and snapshots (c) at time t = 996 . t = 991 . φ = π/ − . ε = 0 . a = 1 . r = 0 . σ = 0 . D = 0 . nodes, on the other hand, are spiking only because of the coupling [27]. The natureof the switching effect has not been revealed.This work is dedicated to finding the reasons of CR chimeras. It seems that theproperty of these switches should also be caused by the specifics of coupling. However,in order to ensure that the effect of the coupling is not confused with the effect ofnoise, the latter must be excluded. To do this, set the noise intensity to D = 0. In thiscase the same parameters as was for CR chimeras, lead to the complete disappearanceof any oscillations. The noise influence can often lead to a shift of bifurcation values[5, 28]. Therefore, a small variation of the parameter φ may lead to the CR chimeraspredecessor. And this regime is found for φ = π/ . u = 0, ˙ v = 0. It leadsto the solution: u = const = − a , v ( u ) = u − u /
3. So, if a = 1 .
001 there is a nullclinecrossing at the point u = − .
001 and v = u − u /
3. These values correspond tothe coordinates of the equilibrium state in one isolated FitzHugh–Nagumo system.However, the presence of a stable wave regime (Fig. 2) in the ensemble at the samevalues of parameters suggests that the equation should include the influence of thecoupling. Taking into account the components of the coupling, these equations are: v ( u ) = u − u / C u u = − a − C v , (4)where C u , C v are additional values produced by the nonlocal coupling. By convertingthe equation (4) onto the ensemble, we can consider what happens to the nullclines ill be inserted by the editor 5 Fig. 2.
Switching spatial wave with corresponding space-time plot for the variables u ti (a)and v ti (b); and snapshots (c) at time moments t = 997 . t = 992 . φ = π/ . ε = 0 . a = 1 . r = 0 . σ = 0 . D = 0. of the nodes belonging to two different stairs of observed regime (Fig. 2,c). Thecorresponding nullclines and projections to the phase plane are shown in Fig. 3.The gray lines in Fig. 3 represent nullclines for one isolated FHN system. Thecolors show the nullclines taking into account the connections (4) for the oscillators i = 60 (red) and i = 315 (blue). Panel (a) shows that changes of nullclines for the firstspiking oscillator are not essential. Initially, when the chimera states were considered,it was predicted that the incoherent part of the chimera started to spike only becauseof the noise, and the oscillators from the coherent domain make the path along thecycle only because of the coupling. That is why last oscillators make the spikingbehaviour more coherent. Here the red oscillator i = 60 should start to move alongthe cycle first, but there is almost no change in its nullclines. Moreover, now thereis no noise which could excite the oscillations. Nevertheless, it starts to move alongthe limit cycle (b). The blue oscillator belonging to the opposite stair is still nearthe equilibrium state. Significant changes in nullclines occur only when most of theoscillators are on the opposite side of the limit cycle (c). When all the oscillators havealready begun their journey through the cycle, the blue oscillator joins them. Thepanel (d) clearly shows a change in its v -nullcline. This leads to destruction of thestate of equilibrium. In the end, all oscillators come to the vicinity of the equilibriumstate. The red oscillator made the spike first, and it turns out to be on the right sideof the v -nullcline (e), and the blue one, because it was the last one, turns out to beon the left side. Next half-period, they will switch places, and the blue oscillator willbe the first to start its spiking behaviour.Despite the fact that all nullclines almost coincide in projections on a phase por-trait, a small deviation is enough to destroy the state of equilibrium. The summand C u affects the nullcline ˙ u = 0. It shifts it up or down depending on the sign. Thisaffects the coordinates of the equilibrium state on the phase plane. The summand C v shifts the isocline ˙ v = 0 left or right. It affects the existence of equilibrium point ingeneral. If v -nullcline shifts to the right, then there will be a higher chance of oscillatorgoes out and starts the spiking behaviour [29]. Will be inserted by the editor
Fig. 3.
Nullclines (4) for spatial wave shown in Fig. 2. The red points indicate the node i = 60 and its nullclines. Blue color corresponds to the oscillator i = 315 and its nullclines.The other are shown by dark-gray color. All points shows the projections of oscillators on thephase plane ˙ u , ˙ v . Gray lines represent nullclines for one isolated FitzHugh–Nagumo systemunder the same values of parameters ε = 0 . a = 1 . φ = π/ . r = 0 . σ = 0 . D = 0.ill be inserted by the editor 7 Fig. 4.
Coupling terms in first C u and second C v equations of the system (1) prepared forthe oscillator i = 60. Blue vertical lines shows the presence of oscillator near the equilibriumpoint. Parameters: φ = π/ . ε = 0 . a = 1 . r = 0 . σ = 0 . D = 0. Figure 4 shows temporal implementations of C u and C v coupling values. Findingall oscillators near the equilibrium state is accompanied by closeness of both values of C u and C v to zero. This happens when the oscillator i = 60, for which this graph isprepared, starts to spike first and last. Thus, we can assume that the oscillator thatis the first one is influenced by the sign of these values. From Fig. 4 we can see thatbefore the phase “last” both values of C u and C v of the oscillator i = 60 have veryclose to zero values. If we consider what happens before the phase “first”, we canclearly see that the summand C v has a small negative value. When C v is less thanzero, v -nullcline leads to the value u = − a + | C v | , and the nullcline shifts to the right.Let us consider what shifts the nullcline. At φ = π/ . b vv = cos ( φ ) ≈ − . b vu ≈ − . v values and a small disorder of u -variables. Therefore, the main contribution is madeby the following summands with cross-link: (cid:80) j b vu ( u j − u i ). This summand has anegative strength b vu ≈ − . C v is negative, all u -variables ofneighbouring oscillators u j must be greater than u i . This occurs for the oscillator,which comes to the state of equilibrium last. Therefore, the same oscillator becomesthe first one in the next half-period.Thus, the existence of CR chimeras is mainly caused by a special type of spatialprofile. It appears in the system without noise due to the presence of cross-link in thecoupling. In the case of a noisy system, there is a slight shift in the parameters at whichthis profile appears. The presence of an incoherent part is caused by noise exposure.The coupling leads to a special movement of nullclines, and due to the noise severaloscillators have the opportunity to start moving first. At φ = π/ . b uu = b vv = cos ( φ ) ≈ − . b uv = sin ( φ ) ≈ . b vu = − sin ( φ ) ≈ − . φ = π/ − . b uu = b vv ≈ . φ is anargument of trigonometric functions, such situation should take place when φ is closeto π/ π/
2. However, the chimera near 3 π/ Will be inserted by the editor
Fig. 5.
New chimera regime with corresponding space-time plot for the variable u ti (a),snapshot (b) at time t = 995 .
8, mean phase velocities ω i (c) and cross-correlation function(d). Parameters: φ = 5 . ε = 0 . a = 1 . r = 0 . σ = 0 . D = 0 . the positive sign of cross-linking parameters is important, which is negative when φ is about 3 π/ This section shows that chimera state in the system under study occurs not onlywhen φ is close to π/
2, but also for other values. The parameter φ was changed inthe interval φ ∈ [0; 2 π ]. It is found that at value φ = 5 . t = 995 .
8. Two breaksof spatial profile are clearly seen between top and bottom parts. The areas of spa-tial incoherence born near the breaks. Such an instantaneous images are typical forchimera states appearing in rings of nonlocally connected chaotic systems. For exam-ple, logistic maps, Henon maps, R¨ossler systems and many others [21, 30–32]. Such atype of a spatial profile in some literature is called “phase chimera”. However, in theensembles of excitable systems (for example, FitzHugh–Nagumo model) this has notbeen encountered before.Position of incoherent and coherent domains does not change in time, and alloscillators continue to demonstrate spiking behaviour (see the space-time plot inFig. 5,a). Due to the stationary position of the incoherent and coherent domains, itis possible to calculate mean phase velocity ω i and cross-correlation function Ψ k,i forthis state. Figure 5,c shows mean phase velocities for each oscillator calculated as ω i = 2 πM i /∆T , i = 1 , . . . , N , where M i is the number of complete rotations aroundthe origin performed by the i th unit during the time interval ∆T = 10 . Almostall values of ω i lie on a continuous curve. It means that all the oscillators makestheir spikes with the same periodicity. Small deviations from the constant value are ill be inserted by the editor 9 observed in the minima and maxima of the spatial profile, which indicates the possibledynamics with weak chaos in these regions. This feature will be discussed a while later.It may seem that each oscillator from the incoherent domain belong to top orbottom parts of spatial profile all the time. However, this is not true. At long timeintervals oscillators inside incoherent parts can change their belonging. The oscillatorslocated near the boundaries between incoherent and coherent regions are especiallysusceptible to this. It is very hard to see it on the space-time plots. Since the incoherentdomains are now stationary, we can use the cross-correlation function [33], whichshows temporal correlation between two oscillators (Fig. 5,d): Ψ k,i = (cid:104) ˜ x k ( t )˜ x i ( t ) (cid:105) (cid:112) (cid:104) ˜ x k ( t ) (cid:105)(cid:104) ˜ x i ( t ) (cid:105) , (5)where ˜ x i ( t ) = x ( t ) − (cid:104) ˜ x i ( t ) (cid:105) is a deviation from the mean value. The brackets (cid:104)· · ·(cid:105) denote the time averaging. This characteristic is equal to 1 and -1 for in-phase andanti-phase oscillations respectively and is not equal to it for non-synchronized dy-namics.The calculations of the cross-correlation function in Fig. 5,d is prepared in rela-tion to the oscillator k = 175 belonging to the top coherent region (Fig. 5,b). Thecross-correlation function is close to the value 1 only when comparing oscillators fromone coherent part. Oscillators from the opposite coherent part are not in anti-phasewith them, because the corresponding values of Ψ are not equal to -1. This can beexplained by the fact that each oscillator continues to demonstrate spiking behaviour,in which it is impossible to say about the presence of phase or antiphase synchro-nization, they are just shifted in time. However, Fig. 5,d clearly shows that the crosscorrelation decreases near the boundaries between coherent and incoherent domains.This confirms that some of the oscillators there may be thrown between what a partof the profile they belong to.Another interesting feature, which is clearly indicated by the cross-correlationfunction, is that there is some desynchronization in the middle of one of the coherentdomain. There is a local minima at i ∈ (200; 300) in Fig. 5,d. This is also one ofthe features of chimeras arising in ensembles of chaotic oscillators, where this type ofchimeras is called amplitude ones. As usual, they have a finite life time in the centerof coherent areas.The chimera state shown in Fig. 5 is obtained for the same parameter values asCR chimera but when the parameter φ is in a quite large range around the value5 .
4. At this value φ the coupling matrix has values: b uu = b vv ≈ . b uv = − b vu ≈ − . φ parameter on the newchimera state. Fig. 6 shows the areas of existence of the new chimera on the parameterplane ( φ, D ) for two values of the parameter a = 1 .
001 (near an Andronov–Hopfbifurcation) and a = 1 .
01 (far from it). The panel (a) shows the chimera at a =1 . t = 10 with the integration step of h = 0 . φ, D ) at a = 1 .
01. As for panel (a), several random initial conditions are considered. Thechimera coexists with the equilibrium state during the whole chimera region. Theclear area of stability, as it was for a = 1 . Fig. 6.
New chimera on the parameter plane ( φ, D ) shown by gray filled area. The hatchedregion corresponds to stable chimera state obtained from each random initial condition.Parameters: φ = 5 . ε = 0 . r = 0 . σ = 0 .
4. Initial conditions are random distributed inthe ring of radius 2.
In both cases the most often regime at zero value of noise intensity D = 0 consistsin weak oscillations near the equilibrium point. However there are rare cases whenfrom random initial conditions the spatial profile with the breaks appears (as inFig. 5, but without incoherent domains). At D = 0 this regime has a finite lifetime.The chimera can also be obtained, but only from specially prepared initial conditions,and it has a finite lifetime too.Another difference between the effect of noise at a=1.001 and a=1.01 is the fol-lowing. The increase of the noise near the bifurcation (Fig. 6,b at a = 1 . a = 1 .
01, the opposite ishappening. Noise leads to the decrease and change of the area of existence.The number of incoherent domains in chimera state can be increased by reducingthe coupling radius. In the case of chimera shown in Fig. 5 it is accompanied by theincrease of the spatial wave number. Figure 7 shows main regimes obtained on theplane of coupling parameters ( σ, r ). It can be seen from the figure that at a largecoupling strength σ and a large coupling radius r , weak oscillations near the stateof equilibrium (yellow region) are observed in the system. The same effect has beenfound out for CR chimera [22]. At small coupling impact the spatial incoherenceand various unstable modes are observed in the ensemble (remains white in Fig. 7).Chimera states are realized between these two main areas (purple area in Fig. 7).Inside this area different wave numbers can be obtained. Here only three are indicated: K = 1 (dotted), K = 2 (diagonal hatching) and K = 3 (vertical hatching). The mapof regimes is made at φ = 5 . D = 0 . The chimera states arising in the ring of nonlocally coupled excitable FitzHugh–Nagumo systems are considered in this article. Such system can demonstrate coherence-resonance chimera for some values of parameters [22]. Previously, it was assumed thatthe chimera was caused only by noise. However, no answer was found to the questionabout periodic switching of the incoherent cluster position. It is shown here that the ill be inserted by the editor 11
Fig. 7.
New chimera on the parameter plane ( σ, r ) shown by purple area. There are threeindicated spatial wave numbers: K = 3 (vertically hatched), K = 2 (diagonally hatched) and K = 1 (dotted). The yellow region correspond to the weak oscillations near the equilibriumpoint.Parameters: φ = 5 . ε = 0 . a = 1 . D = 0 . parameter φ , which is responsible for the influence of cross-linking and self-linking in-side the coupling terms, leads to a special type of spatial and temporal dynamics evenif the noise is switched off. In this case, a special spatial wave is formed. This mode isthe predecessor of the CR chimera. It contains several neighbouring oscillators, whichbegin to make a spiking event first. The other oscillators joins them later becauseof the coupling. The position of the first spiking oscillators changes periodically intime. Adding noise to the system causes these oscillators to spike less regularly. Thiscreates an incoherent cluster, and this wave regime is transforming into CR chimera.In addition, it is found that if the impact of self-linking and cross-linking is compa-rable ( φ ≈ . References
1. G. Hu, T. Ditzinger, C. Z. Ning, and H. Haken, Phys. Rev. Lett. , 807 (1993).2. A. Pikovsky and J. Kurths, Phys. Rev. Lett. , 775 (1997).3. A. B. Neiman, P. I. Saparin, and L. Stone, Phys. Rev. E , 270 (1997).4. B. Lindner, J. Garc´ıa-Ojalvo, A. B. Neiman, and L. Schimansky-Geier, Phys. Rep. , 321 (2004).5. O. V. Ushakov, H. J. W¨unsche, F. Henneberger, I. A. Khovanov, L. Schimansky-Geier, and M. A. Zaks, Phys. Rev. Lett. , 123903 (2005).6. A. Zakharova, T. Vadivasova, V. S. Anishchenko, A. Koseska, and J. Kurths,Phys. Rev. E , 011106 (2010).7. A. Zakharova, A. Feoktistov, T. Vadivasova, and E. Sch¨oll, Eur. Phys. J. Spec.Top. , 2481 (2013).8. R. Benzi, A. Sutera, and A. Vulpiani, J. Phys. A , L453 (1981).
9. M. Masoliver, N. Malik, E. Sch¨oll, and A. Zakharova, Chaos , 101102 (2017).10. M. E. Yamakou and J. Jost, Nonlinear Dynamics , 2121 (2018).11. A. V. Bukh, A. V. Slepnev, V. S. Anishchenko, and T. E. Vadivasova, Regularand Chaotic Dynamics , 325 (2018).12. S. Majhi, B. K. Bera, D. Ghosh, and M. Perc, Physics of Life Reviews , 100(2019).13. Y. Kuramoto and D. Battogtokh, Nonlin. Phen. in Complex Sys. , 380 (2002).14. D. M. Abrams and S. H. Strogatz, Phys. Rev. Lett. , 174102 (2004).15. D. M. Abrams, R. E. Mirollo, S. H. Strogatz, and D. A. Wiley, Phys. Rev. Lett. , 084103 (pages 4) (2008).16. M. J. Panaggio and D. M. Abrams, Nonlinearity , R67 (2015).17. Y. Maistrenko, O. Sudakov, O. Osiv, and V. Maistrenko, New J. Phys. , 073037(2015).18. B. K. Bera, S. Majhi, D. Ghosh, and M. Perc, Europhys. Lett. , 10001 (2017).19. S. Majhi, M. Perc, and D. Ghosh, Chaos: An Interdisciplinary Journal of Non-linear Science , 073109 (2017).20. I. Omelchenko, T. H¨ulser, A. Zakharova, and E. Sch¨oll, Frontiers in AppliedMathematics and Statistics , 67 (2019).21. I. Omelchenko, O. E. Omel’chenko, P. H¨ovel, and E. Sch¨oll, Phys. Rev. Lett. ,224101 (2013).22. N. Semenova, A. Zakharova, V. S. Anishchenko, and E. Sch¨oll, Phys. Rev. Lett. , 014102 (2016).23. R. Kozma, Phys. Lett. A , 85 (1998).24. M. M. Van Hulle, IEEE Trans. Neural Netw. , 850 (2004).25. J. A. Henderson and P. A. Robinson, Phys. Rev. Lett. , 018102 (2011).26. M. Wolfrum, O. E. Omel’chenko, S. Yanchuk, and Y. Maistrenko, Chaos ,013112 (2011).27. A. Zakharova, N. Semenova, V. S. Anishchenko, and E. Sch¨oll, in: Patterns of Dy-namics (eds. P. Gurevich, J. Hell, B. Sandstede, A. Scheel), Springer Proceedingsin Mathematics and Statistics (2017).28. C. Meunier and A. D. Verga, Journal of Statistical Physics , 345 (1988).29. E. M. Izhikevich, Dynamical Systems in Neuroscience (MIT Press, Cambridge,MA, 2007).30. V. Semenov, A. Feoktistov, T. Vadivasova, E. Sch¨oll, and A. Zakharova, Chaos , 033111 (2015).31. I. A. Shepelev, G. I. Strelkova, and V. S. Anishchenko, Chaos: An InterdisciplinaryJournal of Nonlinear Science , 063119 (2018).32. P. Chandran, R. Gopal, V. K. Chandrasekar, and N. Athavan, Chaos: An Inter-disciplinary Journal of Nonlinear Science , 053125 (2019).33. T. E. Vadivasova, G. I. Strelkova, S. A. Bogomolov, and V. S. Anishchenko,Chaos: An Interdisciplinary Journal of Nonlinear Science26