Robustness of chimera states for coupled FitzHugh-Nagumo oscillators

Robustness of chimera states for coupled FitzHugh-Nagumo oscillators

Iryna Omelchenko Institut für Theoretische Physik, Technische Universität Berlin, Hardenbergstraße 36, 10623 Berlin, Germany    Astero Provata    Johanne Hizanidis Institute of Nanoscience and Nanotechnology, National Center for Scientific Research “Demokritos”, 15310 Athens, Greece    Eckehard Schöll Institut für Theoretische Physik, Technische Universität Berlin, Hardenbergstraße 36, 10623 Berlin, Germany    Philipp Hövel Institut für Theoretische Physik, Technische Universität Berlin, Hardenbergstraße 36, 10623 Berlin, Germany Bernstein Center for Computational Neuroscience, Philippstraße 13, 10115 Berlin, Germany

Chimera states are complex spatio-temporal patterns that consist of coexisting domains of spatially coherent and incoherent dynamics. This counterintuitive phenomenon was first observed in systems of identical oscillators with symmetric coupling topology. Can one overcome these limitations? To address this question, we discuss the robustness of chimera states in networks of FitzHugh-Nagumo oscillators. Considering networks of inhomogeneous elements with regular coupling topology, and networks of identical elements with irregular coupling topologies, we demonstrate that chimera states are robust with respect to these perturbations, and analyze their properties as the inhomogeneities increase. We find that modifications of coupling topologies cause qualitative changes of chimera states: additional random links induce a shift of the stability regions in the system parameter plane, gaps in the connectivity matrix result in a change of the multiplicity of incoherent regions of the chimera state, and hierarchical geometry in the connectivity matrix induces nested coherent and incoherent regions.

nonlinear systems, dynamical networks, coherence, spatial chaos
05.45.Xt, 87.18.Sn, 89.75.-k

I Introduction

Synchronization in systems of coupled oscillators has been investigated in various fields such as nonlinear dynamics, network science, and statistical physics, and found its applications in physics, biology, and technology. Fascinating dynamical scenarios, which emerge for coupled systems, have attracted much attention in the scientific community PIK01 (); BOC06a (). For instance, a very peculiar type of dynamics called chimera states, when a network exhibits a hybrid state combining both coherent and incoherent parts, was first reported for coupled phase oscillators KUR02a (); ABR04 (). The surprising aspect of this phenomenon is that these states were detected in systems of identical oscillators coupled in a symmetric ring topology with a symmetric interaction function, and coexist with a completely synchronized state. The last decade has seen an increasing interest in chimera states LAI09 (); MOT10 (); MAR10 (); OME10a (); OME12a (); MAR10b (); WOL11a (); BOU14 (). It was shown that they are not limited to phase oscillators, but can be found in a large variety of different systems including time-discrete maps OME11 () and time-continuous chaotic models OME12 (), and neural systems OME13 (); HIZ13 (). Moreover, chimera states were found in systems with higher spatial dimensions OME12a (); SHI04 (); PAN13 (); PAN14 (). Together with initially reported chimera states, which consist of one coherent and one incoherent domain, new types of these peculiar states having multiple incoherent regions SET08 (); OME13 (), as well as amplitude-mediated SET13 (); SET14 (), and pure amplitude chimera and chimera death states ZAK14 () were discovered.

Potential applications of chimera states in nature include the phenomenon of unihemispheric sleep in birds and dolphins RAT00 (), bump states in neural systems LAI01 (); SAK06a (), in power grids FIL08 (), or in social systems GON14 (). Many works considering chimera states were mostly based on numerical results. A deeper bifurcation analysis OME13a () and even a possibility to control chimera states SIE14c () were obtained only recently.

Ten years after the first numerical observation, the experimental verification of chimera states was demonstrated in chemical TIN12 (); NKO13 () and optical HAG12 () systems. Further experiments involved mechanical MAR13 (), electronic LAR13 (), and electrochemical SCH14a (); WIC13 () oscillator systems.

Identical elements and symmetric coupling topology were widely assumed to be the necessary ingredients for chimera states. Recent studies show that one can overcome these limitations and chimera-like states can be found also when the elements of the system are nonidentical LAI10 (), or when the topology is not regular KO08 (); SHA10 (); LAI12 (); YAO13 (); ZHU14 () or even global SCH14a (); YEL14 (). Then, the spatial order in the system is lost and the assignment to coherent or incoherent regions is based, for instance, on the dynamical variables. Hence, the intriguing phenomena of chimera states continues to show new aspects and is far from being completely explored.

In the present study, we discuss the robustness of chimera states in systems of nonlocally coupled FitzHugh-Nagumo (FHN) oscillators FIT61 (); NAG62 () studied in Ref. OME13 (). We consider two types of inhomogeneities. The first type is realized by inhomogeneous parameters of the local elements in a regular ring topology. As a second type, we study systems with identical units, but irregular coupling topologies. We find that for small inhomogeneities, chimera states are robust under these perturbations, and we discuss qualitative changes of the chimera states caused by different types of irregularities.

The motivation for studying irregular coupling topologies comes from recent results in the area of neuroscience. Diffusion Tensor Magnetic Resonance Imaging (DT-MRI) studies revealed an intricate architecture in the neuron interconnectivity of the human and mammalian brain, which has already been used in simulations  VUK14 (). The analysis of DT-MRI images (resolution of the order of mm) has shown that the connectivity of the neuron axons network represents a hierarchical geometry with fractal dimensions varying between 2.3 and 2.8, depending on the local properties, on the subject, and on the noise reduction threshold KAT09 (); EXP11 (); PRO12 (). Based on these findings, we study the development of chimera states in coupled neurons operating in networks which involve connectivity gaps as well as in topologies with hierarchical connectivity. In both cases, we assess the influence of the irregular connectivity on the properties of the chimera state.

The rest of the paper is organized as follows: In the next section we introduce the model, a set of nonlocally coupled FHN oscillators, and we briefly review the conditions for the appearance of chimera states. In Sec. III, we introduce an inhomogeneity in the system parameters via randomly chosen frequencies of the local units, and discuss the robustness of chimera and multi-chimera states. In Sec. IV, we include additional random links on top of the regular nonlocal coupling scheme and obtain stability regions for chimera states depending on the probability of the introduction of new links. In Sec. V, we consider symmetric and asymmetric gaps in the nonlocal coupling between each node and its neighbors. Numerical results show that when gaps are introduced, the multiplicity of the chimera states may change depending on the position of the gaps. Hierarchical geometry in the connectivity between each node and all other nodes of the network is discussed in Sec. VI. Multi-chimeras with nested shapes appear as a result of hierarchical connectivity and appropriate choice of parameters. The main conclusions are summarized in Sec. VII and open problems are discussed.

Ii The model

We consider a 1-dimensional ring of nonlocally coupled FHN oscillators, where each element is coupled to neighbors on either side OME13 ():


where and are the activator and inhibitor variables, respectively FIT61 (); NAG62 (), and is a small parameter characterizing a timescale separation, which we fix at throughout the paper. denotes the coupling strength. It is convenient to consider the ratio , called coupling radius, which ranges from (nearest-neighbor coupling) to 0.5 (all-to-all coupling). All indices are modulo . Depending upon the threshold parameter , , each individual FHN unit exhibits either oscillatory () or excitable () behavior. In this study, we assume that the elements are in the oscillatory regime, i.e., .

The important feature of Eqs. (1a,b) is that they contain not only direct - and - couplings, but also cross-couplings between activator () and inhibitor () variables, which we model by a rotational coupling matrix:


Therefore, the matrix is determined by the coupling phase .

The existence of chimera states for the system (1a,b) with identical elements was reported in Ref. OME13 (). There, applying a phase-reduction technique IZH00b (), we identified nonlocal off-diagonal coupling to be a crucial ingredient for the occurrence of chimera states in the system. In addition, multi-chimera states consisting of multiple domains of incoherence were observed.

Figure 1: (Color online) Snapshots of the variables and mean phase velocities for inhomogeneous oscillators. (a),(b),(c) ; (d),(e),(f) ; (g),(h),(i) ; values of shown for each panel above the mean phase velocities plots. System parameters: , , , and .

Iii Inhomogeneous elements

In this section we consider system  (1a,b) in the case when the local FHN units are nonidentical and have different frequencies. This feature can be achieved by random choice of the parameters , because the parameter determines the period of the individual element BRA09 (). Our attention is focused on the case when the threshold parameters are drawn from a normal (Gaussian) distribution with mean value and variance . In our numerical simulations, we have also considered the case when the threshold parameters were drawn from a uniform random distribution, and obtained qualitatively similar results.

Together with the coupling radius , coupling strength , and phase of the rotational local interaction matrix, the parameters and now control the system dynamics. In the following, we choose , and vary the variance . This quantity will define the amount of inhomogeneity of the elements in our system. To compare our results with the system of identical elements, we will follow the choice of system parameters as considered in Ref. OME13 ().

A significant feature of chimera states is the difference of average frequencies of the oscillators that belong to the coherent or incoherent parts, respectively. This feature can be visualized by the mean phase velocities for each oscillator calculated as , , where denotes the number of oscillations and is the simulation time. The profile of the mean phase velocities is usually characterized by a plateau of equal frequencies that correspond to the coherent domain of the chimera state, and an arc-like part corresponding to the incoherent part.

Figure 1 shows typical snapshots and corresponding profiles of the mean phase velocity for chimera states with one, two, and three incoherent regions, as the inhomogeneity in the system increases. The number of incoherent regions depends on the coupling strength and coupling range OME13 (). Starting the simulation for each panel with initial conditions that are randomly distributed on the circle , we obtain chimera states with well pronounced mean phase velocity profiles that demonstrate a clear difference between coherent and incoherent domains. When the system elements have very close frequencies and the variance is small [Fig. 1(a),(d),(g)], chimera states with one and multiple incoherent regions are robust. The difference between the maximum mean phase velocity of the incoherent domain and the -value of the coherent domain is larger in the single chimera state than in the multi-chimera states. For a small increase of the inhomogeneity, e.g. [Fig. 1(b),(e),(h)], the chimera state with one incoherent domain remains visible, but the profiles of the mean phase velocity for chimera states with two and three incoherent domains become more noisy around the incoherent parts. For even larger inhomogeneities such as [Fig. 1(c),(f),(i)], the multiple incoherent domains merge into one larger region. However, the corresponding mean phase velocity profile does not resemble that of the original single chimera state [Fig. 1(c)], which maintains the arc-shaped form of the homogeneous system. For the chimera state with one incoherent domain we still clearly observe the difference of the mean frequencies for incoherent and coherent domain. These observations can be explained as follows: In multi-chimera states frequencies are less distinguishable (smaller ). Therefore, introducing inhomogeneities leads to the collapse of the profile to a noisy configuration, where incoherent oscillators take over the coherent ones.

In the numerical simulations of the system (1a,b) for fixed values of and , we use different realizations of the threshold parameters Moreover, as we consider randomly chosen initial conditions, it is important to understand how their choice influences the obtained system solutions.

Figure 2: (Color online) Averaged mean phase velocities for chimera states with one and multiple incoherent domains with increasing system inhomogeneity. Black denotes phase velocities for the system of identical oscillators , red (gray)- mean phase velocities averaged over realizations, gray (light gray)- variations of phase velocities in realizations. (a),(b),(c) ; (d),(e),(f) ; (g),(h),(i) ; values of are shown inside each panel. Other system parameters as in Fig. 1.

To address this issue, Fig. 2 depicts the comparison of mean phase velocities for chimera states with one, two, and three incoherent domains with increasing of the system inhomogeneity. Red (gray) lines in Fig. 2 show the averaged mean phase velocities over realizations, and gray (light gray) error bars denote the variations of the obtained values, which is given by the minimum and maximum. Black lines show the values of mean phase velocities for system with identical elements for comparison.

For a small value of [Fig. 2(a),(d),(g)], we obtain mean phase velocity profiles with pronounced differences between the coherent and incoherent domains and the variations of the velocities are in general small. As inhomogeneity increases, [Fig. 2(b),(e),(h)], the variations for mean phase velocities of multi-chimera states become more pronounced, although the average value still shows maxima for incoherent domains. For larger inhomogeneities, [Fig. 2(c),(f),(i)], the variations are strong and appear even in the coherent domain and multiple incoherent regions collapse into one. The chimera state with one incoherent domain, however, appears to be more robust than in the case of multi-chimera.

The oscillation period of each FHN system strongly depends on the value of its threshold parameter  BRA09 (). With a larger range for the values of , we increase the spectrum of their frequencies, which are then additionally influenced by the coupling term. This results in the changes of the individual frequencies of the oscillators in the inhomogeneous system. Therefore, the mean phase velocities for the oscillators that belong to the coherent domains of multi-chimera states are smaller and then become larger compared to the case of identical elements when the value of increases. See black and gray (red) curves in Fig. 2.

To summarize this section, inhomogeneous system parameters in rings of nonlocally coupled FHN oscillators still allow us to observe chimera states. Classical chimera states with one incoherent domain are more robust under this perturbation than multi-chimera states with multiple incoherent domains. Chimera states can be observed only for small inhomogeneities. Otherwise, their incoherent parts merge and they form a chimera state with a single, large incoherent part.

Iv Irregular topology: additional random links

Studies of the connectivity in neural networks show that local couplings often coexist with long-range connections between the neurons SCA95 (); HEL00 (); HOL03 (). In such coupling structures, the interplay between the local interactions and long-range shortcuts can give birth to nontrivial dynamical effects. Starting from this idea, we consider in this section a perturbation of the regular ring coupling of system (1a,b). Keeping the original coupling radius, we add new links, i.e. shortcuts, with probability between the elements. No multiple edges between the same nodes are allowed. The transformed system can be written as


where again denotes strength of the coupling. All indices are modulo and the individual FHN elements are in the oscillatory regime and identical, i.e., . is the adjacency matrix with elements or , denotes the number of nearest neighbors, and is the number of new links added to the -th element, which is . The number of connections for each element increases for larger probability . Thus, we rescale the coupling strengths by the updated number of neighbors such that every connection has the same strength, for comparison with the regular ring topology of system (1a,b) OME13 ().

The introduction of shortcuts shares similarities to a well-known small-world network WAT98 () with one difference. Here, we do not rewire already existing links when adding new ones, and the regular nonlocal structure is kept. For large networks, both strategies lead to similar topologies NEW99b ().

The dynamics of the system (3a,b) is defined by five control parameters: the threshold parameter of the individual uncoupled FHN unit, the strength of the coupling , the phase of the rotational matrix defining local interaction scheme, the regular coupling radius , and the probability of adding new links. The last two parameters give us an approximate number of interactions for each system element, that is, an effective coupling radius , which can be calculated as

Figure 3: (Color online) Snapshots of the variables (left panels) and mean phase velocities (right panels) for the system Eqs. (3a,b). Probability of adding new links: (a) , (b) , (c) , (d) . Parameters: coupling radius , coupling strength , , and other parameters as in Fig. 1.

Figure 3 shows examples of chimera states obtained in the system (3a,b) for intermediate regular coupling radius and small coupling strength , with increasing probability for new links. This parameter choice corresponds to chimera states in the regularly coupled system with no additional links shown in Fig. 1(a). Left and right panels in Fig. 3 depict snapshots and mean phase velocity profiles, respectively. Looking at the snapshots, chimera states appear to be robust in system (3a,b). Increasing the number of new links results in the following changes in the phase velocity profiles: The characteristic arc-like shape of the incoherent region is still visible, but the border to the coherent domain becomes less sharp. Furthermore, we find that the difference between the maximum and minimum of the average phase velocities, , decreases. This effect can be explained by considering the additional input by the oscillators of the coherent part applied to the oscillators of the incoherent part. Elements, which had a larger difference in the phase velocities in the unperturbed case (Sec. III), are now in contact. This leads to homogenization of their frequencies and thus to a smaller .

Chimera states with multiple incoherent regions are much more sensitive with respect to this topological perturbation of the regular ring topology. Even a small number of shortcuts results in the loss of smaller coherent domains in between the incoherent ones (cf. Figs. 1 and 2). This leads to the formation of a chimera state with one larger incoherent part (results not shown).

Figure 4: (Color online) (a) Stability regions for chimera states with one incoherent domain in the plane of coupling radius  and coupling strength . Gray region corresponds to the system with regular nonlocal coupling topology (), and lines denote the borders of stability region for . The inset shows the same plots rescaled to the effective coupling radius . (b) Fraction of successful realizations of chimera state starting from random initial conditions in numerical experiments for each parameter set. The parameter values denoted by black dots in panel (a) are: , , and . Other parameters as in Fig. 1 with . (c) Fraction of successful realizations of chimera state in the network with rewired links, all parameters as in panel (b).

Figure 4(a) depicts the stability regions for chimera states with increasing probability for random shortcuts. The gray region shows the stability region of chimera state with one incoherent domain in the original system (1a,b) without additional random links, reproduced from Ref. OME13 () as a reference case. The increase in the number of random links in the network results in a transformation of the stability region: the regions change their form and their left border moves in the direction of smaller coupling radii. Hence, additional random shortcuts allow for the existence of chimera states for smaller regular coupling radius, where they were not possible in the system with regular coupling topology.

The left border of the stability regions for different values of , which can be conceived as the critical coupling radius, is given by [see Fig. 4(a)]. According to Eq. (4), the corresponding values of the critical effective coupling radius are , , , and . They appear to be close to the critical coupling radius of the regular system: . In addition, we observe a data collapse by rescaling to the effective coupling radius [see inset of Fig. 4(a)].

Due to the presence of the randomness in the system, we can statistically analyze the persistence of chimeras. More precisely, we investigate the probability to obtain a chimera state starting from random initial conditions with random realization of the system topology for a fixed value of probability . Let us fix three points inside the stability region for the unperturbed regular system , , and . Figure 4(b) shows the fraction of successful realizations of chimera states, which we obtained starting from random initial conditions and random setup of shortcuts over 500 realizations for each parameter set. For a small number of shortcuts, we always obtain a chimera state. This also confirms the fact that points A, B, and C remain inside the stability regions shown in Fig. 4(a). With increasing and correspondingly a larger number of shortcuts, the fraction of realizations with chimera states decreases. The behavior of the system becomes highly multistable and the basin of attraction for chimera states shrinks. Surprisingly, even for very large number of shortcuts, there is still a small probability to obtain a chimera state in the system, although the coupling is close to all-to-all.

For comparison, let us consider now a classical small-world network WAT98 (), where the existing links are rewired with probability to a random node. In this case, the effective coupling radius will coincide with the coupling radius of the original regular network. The probability of adding random links serves again as a measure of the network irregularity. Figure 4(c) depicts the fraction of successful realizations of chimera states over 500 numerical realizations in a similar way as in Fig. 4(b), and for the same system parameters. Even for a small amount of random links, only a part of the realizations lead to chimera states, in contrast to the previous case. This can be explained by the fact that rewiring of the existing links destroys the strong effect of the nonlocal coupling kernel, still present in networks in which the regular ring is maintained as a backbone. Increasing the probability , the irregularity in the coupling topology becomes stronger, and the basin of attraction for chimera states decreases.

We conclude that adding of random links to the system of nonlocally coupled FHN oscillators still allows us to observe chimera states. For small probabilities, we are able to obtain stability regions for chimera states in the case when we do not rewire already existing links, and the effect of nonlocal coupling kernel in the system predominates. For larger numbers of additional links, the basin of attraction of chimera states decreases rapidly, although chimera states still can be observed even for a large amount of added random links, in both cases when the existing links are rewired or kept.

The numerical simulations presented in Fig. 4(b),(c) were obtained for systems of elements. We performed analogous simulations for other systems sizes, and obtained similar results for larger systems (). For smaller system sizes () we found smaller fractions of successful realizations of chimera states. This can be explained by the fact that for smaller systems chimera states perform a random spatial motion of their coherent/incoherent regions (drift) OME10a (). Then, the averaged mean phase velocity profiles do not show the characteristic arc shape, and chimera states are barely detectable.

V Homogeneous and Inhomogeneous Gap Distributions

After the investigation of shortcuts added to the regular ring, we will address another modification of the topology Eqs. (3a,b) in this section that gives rise to a different form of inhomogeneities. These are associated with the presence of gaps in the connectivity matrix. This is inspired by current research in neuroscience, and in particular from functional magnetic resonance imaging (fMRI) task experiments used to identify connectivity between different parts of the brain. Early fMRI investigations have shown that for simple tasks, such as parroting or finger tapping, only some parts of the brain show coherent electromagnetic activity, while the rest shows normal disordered activity. This observation indicates that the overall connectivity is fragmented with gaps in the connectivity matrix. In other words, each node has a preferable attachment to some groups of nodes and a weak attachment to others.

Before introducing inhomogeneities in the connectivity matrix, we consider the question if synchronization is related to the symmetry in the link arrangement around each node. To answer this question we create an asymmetric, directed network, where each node is linked only to its right neighbors. We compare this symmetry-broken configuration with the case of a symmetric connectivity network, keeping all other parameters fixed. For notational convenience, the node , around which the connectivity is described, will be termed reference node. Next, we consider the case where the reference node is not linked to its direct neighbors, but to neighbors positions away. For one-sided connectivity, node is linked to nodes . For comparison with previous studies OME13 (), we consider the following set of parameters: , , , , , and , which has been studied earlier OME13 ().

Figure 5: (Color online) Snapshots and corresponding mean phase velocity profiles for different connectivity distributions. (a) Symmetric connectivity distribution , , (b) asymmetric connectivity distribution , , and (c) displaced asymmetric connectivity distribution, , , . Other parameters are: , , , , . All simulations start from the same initial conditions.

Figure 5 shows the snapshots and mean phase velocity profiles for the following three cases: (a) the links are symmetrically placed around each node to the left and to the right, (b) all links are placed to the right of each node and , (c) all links are placed to the right of each node [as in case (b)] but there is a displacement in the connectivity of each node. Note that the overall number of neighbors remains constant at for each node.

From Fig. 5 it is obvious that the multiplicity of the chimera state does not change qualitatively. This result suggests that symmetry around the reference node is not essential for the formation of chimeras and has been verified for other coupling radii (results not shown). These observations do not depend on the specific arrangement of the links as long as they are densely packed (without gaps). Effects of gaps are discussed next.

Gaps may be introduced in the connectivity matrix as schematically depicted in Fig. 6. Each oscillator is coupled with a finite number of oscillators, to the left and to the right. The gaps are arranged symmetrically or asymmetrically on the left and the right of each node. If () is the size of a gap to the left (right) of a node, it follows that .

In the symmetric case (), the adjacency matrix , , is given by:


where all indices are taken modulo .

Each node will have a number of connections symmetrically placed to both its left and right with a number of gaps equal to also symmetrically placed around it [see Fig. 6(a)]. Hence, the total number of connections remains . In principle, one can introduce many symmetric gaps on the left and on the right of each element in a similar way, but this is beyond the scope of this study.

Figure 6: Schematic distribution of gaps around a reference node in a ring geometry: (a) symmetrically distributed gaps with ; (b) asymmetrically distributed gaps.

In the general case of inhomogeneous gap distributions, the size of the gaps as well as the size of the linked regions can be different on the left and the right of each element [see Fig. 6(b)], i.e.


Now, and denote the number of links to the left (right) of the th element before and after the gap. Thus, one might have . In the same way, the sizes of the linked regions , , , and need not be related. Multiple inhomogeneous gaps may be introduced in a similar way.

Numerical simulations indicate that the presence of gaps tend to stabilize existing chimeras and to create multi-chimera states. As in the previous section, we always start from initial values randomly placed on the circle of radius , that is . Each oscillator is coupled with others, distributed on equal numbers to the left and to the right, i.e. .

Figure 7: (Color online) Snapshots of the variables and corresponding mean phase velocities for different gap positions. Each element is coupled with others elements symmetrically distributed around each node. The connectivity matrix contains two gaps of size symmetrically placed around each node while the values of and are indicated on each panel. Other parameters: , and . All simulations start from the same initial conditions.

Figure 7 depicts on the left panels (a)-(f) typical snapshots obtained in the case of symmetrically distributed gaps, with different values of and keeping constant the gap sizes . This way only the topology of the connectivity changes and not the actual size of the gaps. The right panels depict the corresponding mean phase velocities and the values of are indicated. Starting with the top panels (shapshot (a) and corresponding phase velocities), where the case of , is depicted, we recover the case of symmetrically distributed links without gaps and a 2-chimera state is observed, as in Ref. OME13 (). As the gaps are introduced and moving closer and closer to the reference node, the mean phase velocities of the (in)coherent regions increase [Fig. 7(a) and (b)]. Moving the gaps close to the center, the number of coherent and incoherent regions changes from 2 to 3. It is interesting to note that the difference in the mean phase velocity between coherent and incoherent regions becomes maximum when the gaps (on the left and right) split the coupling regions approximately into equal parts. This will be discussed in greater detail below.

Overall, Fig. 7 indicates that the position of the gaps may change the multiplicity of the chimera as well as the relative mean phase velocity of coherent and incoherent regions. This effect depends on the complex interplay between the specific pattern of the links and the initial conditions.

Figure 8: (Color online) Measures of coherence-incoherence of a chimera state as a function of the size of the inner coupling regions . Each element is coupled with others elements symmetrically distributed around each node. The connectivity matrix contains two gaps of size symmetrically placed around each node. Other parameters as in Fig. 7. Averages are taken over 10 runs.

For a better understanding of the role of the coupling topology for the development of a chimera state, we introduce and calculate two measures that account for the relative size of the coherent versus incoherent part of the chimera. Let us denote by the mean phase velocity of element and by the mean phase velocity of the coherent parts. Then, the relative size of the incoherent parts of the chimera state is calculated as:


where is the step function which takes the value 1 when its argument takes positive values and zero otherwise. is a small tolerance, which in this case we set to 0.05.

We also define the extensive cumulative size of the incoherent parts as:


This is an extensive measure which represents the area below the arcs in the mean phase velocity profiles. demonstrates the degree of incoherence of the chimera state and is zero for purely coherent states.

Figure 8 displays four measures of coherence as a function of the inner coupling regions , for the above case of a symmetric connectivity matrix with elements and two gaps of equal size symmetrically placed around the reference node: (a) , the mean phase velocity of the coherent region, (b) , the maximum difference between incoherent and coherent values, (c) , the relative size of incoherent regions, and (d) , which provides information about the amount of influence of the incoherent regions. We calculate these measures as a function of , changing the position of the gaps symmetrically around the ring. Averages over 10 runs are shown, starting from different initial conditions randomly chosen on the circle . Averages are necessary to account for stochastic deviations as well as multistability. Error bars in the figure indicate the standard deviations.

Figure 8(a) shows that the average mean phase velocity of the coherent parts increases as the two gaps move symmetrically away from the reference node and assumes a maximum around . At this value the coupling regions before and after the gap are of equal size. As the gaps depart further from the reference node, that is, for further increase of , decreases again to attain the value of the regular ring without gaps for and . When becomes maximum, the difference of the mean phase velocity profiles becomes largest [see Fig. 8(b)] and so does the relative number of incoherent oscillators [see Fig. 8(c)]. At the same value, also attains its highest value.

From this discussion it becomes evident that the gaps may change the position and/or multiplicity of the chimera state. Moreover, the gap position along the ring determines the relative “weight” of the coherent versus incoherent part. Therefore, using appropriate gap positions, one may influence the number of elements belonging to the incoherent parts and also the difference in the mean phase velocities .

Figure 9: (Color online) Snapshots and mean phase velocity of the variable for one-sided connectivity. Each element is coupled with elements to its right, i.e. . The connectivity matrix is asymmetric and the coupling pattern contains one gap of size . The values of and are indicated in the panels. Other parameters as in Fig. 7. All realizations start from the same initial conditions.

Next, we study the question of inhomogeneous gap distribution [Eq. (6)]. In order to avoid introducing too many parameters in the system, we keep , thus allowing links only on the right of each element. We consider one gap of size and various values of and under the condition . This keeps the same number of links as above, for comparison with the previous case of the symmetrically placed gaps. Due to the asymmetry in the coupling, the network becomes directed: if a node receives a signal from a node , the opposite does not necessarily hold as in the symmetric case. In addition, the reference node is not always centrally located as in the symmetric case.

Figure 9 depicts typical snapshots and corresponding mean phase velocity profiles of asymmetric connectivity matrices with one gap for different values of the parameters and , while the sum is kept constant, . Phenomena of splitting and merging of the incoherent parts are typically observed for different spatial arrangements of the gap. By shifting the gap towards the reference node, the difference in the mean phase velocity profile attains maximum values when the gap separates the coupling regions into equal parts. As the gap moves further away from the reference node, the mean phase velocity decreases again. This behavior is verified for different values of and for .

Figure 10: (Color online) Measures of coherence-incoherence of a chimera state. Each element is coupled with elements to its right. The connectivity matrix is asymmetric and contains one gap of size . All other parameters as in Fig. 9. Averages are taken over 10 runs.

To investigate further the properties of coherent-incoherent regions as a function of the position of the gap with respect to the reference node, we calculate the four measures that account for the presence of chimera states, in analogy to the case of symmetrically placed gaps shown in Fig. 8. In the simulations, all other parameters are kept constant, as in Fig. 9, while varying the size of the inner and outer coupling radii and , retaining their sum constant, . For each parameter pair averages were taken over 10 realizations starting from different initial conditions, to take into account stochastic effects and multistability. The calculations of the various measures of coherence are shown in Fig. 10 as a function of .

It is interesting to note that as in the case of symmetric connectivity (Fig. 8), both the mean phase velocity of the coherent regions and the maximum is achieved when the gaps are placed centrally, separating the coupling regions in approximately equal parts [Fig. 10(a) and (b)]. The same observation holds for the other two measures, and . Therefore, the introduction of gaps can be used as a control parameter to modulate the synchronization and the magnitude of phase velocities in the coherent and incoherent parts of the chimera states. The above results are indicative of the complexity induced in the chimera states by the spatial arrangement of the connections even if all the other system parameters are kept fixed.

Vi Hierarchical Gap Distribution and Chimeras

In the previous section, we have presented first evidence that the presence and the position of gaps may modify the coherence properties of the chimera state and its multiplicity. In the simulations, the position and the size of the gaps did not follow any precise distribution, but was specifically designed to exemplarily demonstrate the changes in the chimera properties with the gap position and symmetry. We now proceed to consider distributions of the gaps with specific hierarchical architecture. The study of different hierarchical architectures in the neuron connectivity is motivated by MRI results of the brain structure which show that the neuron axons network spans the brain area fractally and not homogeneously KAT09 (); EXP11 (); KAT12 (). In the remain of this section the word ’fractal’ will be employed to denote mainly hierarchical structures of finite orders , since the human brain has finite size and does not cover all orders, (as in the exact definition of a fractal set).

The fractal connectivity dictates a hierarchical ordering in the distribution of neurons which is essential for the fast and optimal handling of information in the brain. When we study phenomena, which involve coupled neurons, the need to take into account the fractal architecture of the neuron network becomes apparent. This is a problem that has emerged after the recent development of DT-MRI techniques allowing for the detection of the neuron network architecture with high levels of accuracy, while traditionally the synchronization properties of coupled neurons are studied in ring architectures with local, nonlocal, or global couplings.

Simple hierarchical structures can be constructed using the classic Cantor fractal construction process MAN83 (); FED88 (). Using the iterative bottom-up procedure to construct the Cantor set, we create a symbolic sequence consisting of two symbols hierarchically nested in one another. Starting with a base containing symbols (0’s or 1’s) we iterate it a number of times and thus obtain systems of size . By closing this string of symbols in a ring we construct a hierarchical connectivity matrix considering that the symbol 1 denotes the existence of a link, while the symbol “0” the absence of a link, namely,


In this way a connectivity matrix of size is constructed, which contains a hierarchical distribution of gaps with a variety of sizes. The number of times the symbol “1” appears in the base, denoted by , defines formally the fractal dimension of the structure, as . This measure describes perfectly the fractal structure when the number of iterations , while for the structure is called hierarchical, because it has been constructed based on a hierarchical algorithm. In Fig. 11 we present schematically the construction of a ring connectivity using the triadic (101) Cantor set FED88 ().

Figure 11: Schematic distributions of hierarchical gaps in a ring geometry based on the triadic (101) Cantor set.

In the current study we use a construction base . Thus, the size of the network is given in powers of . For comparison with previous results we use the iteration of the fractal giving a total size of the network nodes. We change the value of (number of times the symbol “1” is encountered within the base ) to vary the number of hierarchically coupled elements in the structure and thus to vary the fractal dimension of the structure. Fixing still leaves one degree of freedom, that is, the position of the symbols “1” in the base pattern. This issue will be addressed later. All other parameters, apart from the number and position of the coupled elements, are kept constant as in the previous section and Figs. 7-10.

Figure 12: (Color online) Snapshots of the variable and corresponding mean phase velocities for different hierarchical connectivity matrices. The fractal dimensions of the hierarchical structures are: (a)  (base pattern: ), (b)  (base pattern: ) and (c)  (base pattern: ). Parameters: and other parameters as in Fig. 7.

In Fig. 12(a) the calculations are performed with . The connectivity matrix is very sparse and has a fractal dimension . In the iteration it contains only links of each element with the others, each element failing to link with the remaining elements. The links are mostly isolated while the gaps cover most of the structure. A multi-chimera with multiplicity and clearly identified coherent/incoherent parts is observed. This result agrees with previously published works indicating that the multiplicity of a chimera state is high when the number of links in the ring network is small VUE14a ().

In Fig. 12(b) the number of links was increased, using . The connectivity matrix has a fractal dimension and each element is connected to others. Here, the chimera represents a nested structure, containing 10 coherent/incoherent regions clustered into two parts. The two incoherent parts show a substructure consisting of five closely packed incoherent regions. The two clusters are separated by large coherent regions. This phenomenon will be further explored below. To continue with the dependence of the qualitative and quantitative features of the chimera state on the hierarchical architecture, in Fig. 12(c) we consider . The number of links increases further and now each element is coupled with elements. The incoherent parts seem to merge into a -chimera, but the calculation of the phase velocity demonstrates that this single incoherent region has a substructure with three maxima.

The current results support and extend those of Ref. OME13 (), where it was demonstrated that the increase in the number of links leads to decrease in the chimera multiplicity. In addition, here we give first evidence that if the links are distributed in a hierarchical manner, the corresponding chimera state shows nested complex incoherent patterns. Similar structure of chimera states were found in systems of phase oscillators with repulsive coupling MAI14 ().

Figure 13: (Color online) Snapshots of the the variable (on the left) corresponding mean phase velocities (middle) and details of mean phase velocities (on the right) for different Cantor patterns. The fractal dimensions in all hierarchical structures is constant . The initiation strings are : (a) 110011 and (b) 111010. The left panels are zooms of the corresponding mean phase velocities (middle panels) to demonstrate structural ramifications. Parameters as in Fig. 12.

Next, we investigate the influence of the local structure of the connectivity matrix without changing the fractal dimension. In other words, we only change the local configuration of the fractal. As an example, we use two different fractal realizations that have the same dimension , but start from different initialization strings (base patterns): (a) the string 110011 (symmetric) and (b) the string 111010 (asymmetric). The results are presented in Fig. 13.

In particular, Fig. 13(a) depicts the chimera state when the connectivity matrix is initiated by the string 110011. The corresponding mean phase velocity profile (middle panel) shows six incoherent regions. The right panel of Fig. 13(a) depicts an enlargement for the first 150 oscillators. From this panel, it is evident that the structure is ramified and that within each incoherent regions there are nested coherent and incoherent ones. This is certainly reminiscent of the hierarchical structure of the connectivity matrix which includes nested gaps between linked elements.

We now turn to the second example, the initiation string 111010. The fourth iteration of this string results in producing a connectivity matrix with multiple sizes of gaps and solid interacting regions around each reference element. The resulting picture (snapshot) of the amplitude of the oscillators, Figs. 13(b), demonstrates the presence of a chimera state with fewer coherent/incoherent regions. Zooming into the incoherent regions [Fig. 13(b), right panel] does not demonstrate any trace of ramifications. This picture suggests that (i) the fractal dimension is not the only parameter which determines the multiplicity or the pattern of the chimera state, and (ii) different initiation strings, producing the same fractal dimension, give rise to different chimera patterns.

These observations allow us to propose that the presence or absence of ramifications in Fig. 13 may be attributed to the presence of symmetries in the pattern of the original string. If the size of the system is increased by using higher iteration levels in the connectivity matrix, then the hierarchical structure will be more evident and then nested structures may arise even for original strings which do not present obvious symmetries. This hypothesis needs to be tested with extended simulations in future publications.

Finally, we investigate the change of coherence in the chimera state when we change the coupling strength . We keep the fractal dimension fixed at and consider four iterations. Thus, the coupling architecture contains hierarchically nested links. Figure 14 shows how the chimera multiplicity changes with . In fact, multiplicity does not change, but it becomes ramified. Starting from , two well-defined coherent (and two incoherent) regions develop. The same picture continues for , but the coherent regions slowly acquire a finer structure inside, which can be better quantified by the mean phase velocity profiles (see Fig. 14, right panels). The fine structures become best developed for , while they slowly decay for larger values. Note that even for we may observe a wavy structure in the mean phase velocities around the maxima of the incoherent regions, which is not observed in the well formed arc shape of the classical chimeras (see Figs. 2 and 3). Therefore, the coupling strength acts as a control parameter for turning the ramifications on and off.

These findings constitute some first observations of ramified structures emerging when the link architectures are not solid but contain nested, hierarchical gaps. It is important to add here that for a given system size it is not possible to predict the size and number of coherent and incoherent gaps, solely from the fractal dimension and/or the number of links in the system. The interplay between the hierarchical geometry of the links and the randomness of the initial conditions may lead to different chimera profiles as was seen in Fig. 13 where in cases (a) and (b) all parameter were identical including identical initial conditions and number of links, the only difference being in the spatial/hierarchical arrangement of the links. For this reason it is not helpful to study the chimera statistics, but rather one needs to address the basin of attraction that leads to specific chimera patterns.

The general conclusion of this section is that connectivity matrices, which contain hierarchically nested gaps around the reference element, commonly induce nested coherent and incoherent regions. These nested structures are best demonstrated by the mean phase velocity diagrams. In particular, for the same fractal dimension the structure of the ramifications in the incoherent regions is tuned by the value of , and by the pattern of the initial string used to construct the connectivity matrix in connection with the initial conditions.

Figure 14: (Color online) Snapshots of the the variable and the corresponding mean phase velocities for different coupling strengths. The fractal dimensions of the hierarchical structure is constant . The panels correspond to: (a) , (b) , (c) (d) [same as Fig. 12(b)], (e) , (f) . Initiation string is and other parameters are as in Fig. 12.

Vii Conclusions

In the current study, we have addressed the issue of robustness of chimera states in systems of nonlocally coupled units with respect to perturbations of the frequencies of individual elements, and structural transformations of the network topology. We demonstrated that in the system of nonidentical FitzHugh-Nagumo oscillators with regular nonlocal coupling topology, chimera states are robust for small inhomogeneity. As the inhomogeneity increases, chimera states with multiple incoherent regions transform into chimera states with one incoherent region. This finding could be important from the point of view of experiments and applications. In real-world situations, one hardly finds a system with absolutely identical elements.

Random structural perturbations of the network topology do not destroy the chimera states immediately, and for relatively small numbers of random links, chimera states, as well as their stability regions, can be determined. We have compared the two cases when already existing links in the nonlocally coupled network are either rewired or kept, and provide a statistical analysis of chimera occurrences as the irregularity of the network topology increases. Chimera states might still be realized even for a large number of random links.

The introduction of connectivity gaps has allowed us to analyze the symmetric and asymmetric network topologies, when gaps are on both or only on one side of the reference node, respectively. We have demonstrated numerically that the position of the gaps influences the multiplicity and the relative weight of the incoherent and coherent regions of chimera states. Moreover, when the gaps separate the connected areas into equal parts, the difference between the mean phase velocities of the coherent and incoherent domains attains a maximum.

In the case of multiple connectivity gaps, we have shown that when the gap distribution is hierarchical, the resulting chimera pattern exhibits nested incoherent parts. The specific structure depends on the complex interplay between the hierarchical link geometry and the initial conditions.

Our findings corroborate the universal existence of chimera states and extend the variety of systems where chimera states can be found. These results can be useful from the point of view of applications, where inhomogeneity of the elements, and more complex coupling topologies are common.

This work was supported by the German Academic Exchange Service DAAD and the Greek State Scholarship Foundation IKY within the PPP-IKYDA framework. IO and PH acknowledge support by BMBF (grant no. 01Q1001B) in the framework of BCCN Berlin (Project A13). ES, IO, and PH acknowledge support by Deutsche Forschungsgemeinschaft in the framework of SFB 910. The project was also funded by the John S. Latsis Public Benefit Foundation.


  • (1) A. Pikovsky, M. G. Rosenblum, and J. Kurths, Synchronization, A Universal Concept in Nonlinear Sciences (Cambridge University Press, Cambridge, 2001).
  • (2) S. Boccaletti, V. Latora, Y. Moreno, M. Chavez, and D. U. Hwang, Phys. Rep. 424, 175 (2006).
  • (3) Y. Kuramoto and D. Battogtokh, Nonlin. Phen. in Complex Sys. 5, 380 (2002).
  • (4) D. M. Abrams and S. H. Strogatz, Phys. Rev. Lett. 93, 174102 (2004).
  • (5) C. R. Laing, Physica D 238, 1569 (2009).
  • (6) A. E. Motter, Nature Physics 6, 164 (2010).
  • (7) E. A. Martens, C. R. Laing, and S. H. Strogatz, Phys. Rev. Lett. 104, 044101 (2010).
  • (8) O. E. Omel’chenko, M. Wolfrum, and Y. L. Maistrenko, Phys. Rev. E 81, 065201(R) (2010).
  • (9) O. E. Omel’chenko, M. Wolfrum, S. Yanchuk, Y. L. Maistrenko, and O. Sudakov, Phys. Rev. E 85, 036210 (2012).
  • (10) E. A. Martens, Chaos 20, 043122 (2010).
  • (11) M. Wolfrum, O. E. Omel’chenko, S. Yanchuk, and Y. L. Maistrenko, Chaos 21, 013112 (2011).
  • (12) T. Bountis, V. Kanas, J. Hizanidis, and A. Bezerianos, Eur. Phys. J. Special Topic 223, 721 (2014).
  • (13) I. Omelchenko, Y. L. Maistrenko, P. Hövel, and E. Schöll, Phys. Rev. Lett. 106, 234102 (2011).
  • (14) I. Omelchenko, B. Riemenschneider, P. Hövel, Y. L. Maistrenko, and E. Schöll, Phys. Rev. E 85, 026212 (2012).
  • (15) I. Omelchenko, O. E. Omel’chenko, P. Hövel, and E. Schöll, Phys. Rev. Lett. 110, 224101 (2013).
  • (16) J. Hizanidis, V. Kanas, A. Bezerianos, and T. Bountis, Int. J. Bifurcat. Chaos 24, 1450030 (2014).
  • (17) S.-I. Shima and Y. Kuramoto, Phys. Rev. E 69, 036213 (2004).
  • (18) M. J. Panaggio and D. M. Abrams, Phys. Rev. Lett. 110, 094102 (2013).
  • (19) M. J. Panaggio and D. M. Abrams, arXiv:1403.6204 (2014).
  • (20) G. C. Sethia, A. Sen, and F. M. Atay, Phys. Rev. Lett. 100, 144102 (2008).
  • (21) G. C. Sethia, A. Sen, and G. L. Johnston, Phys. Rev. E 88, 042917 (2013).
  • (22) G. C. Sethia and A. Sen, Phys. Rev. Lett. 112, 144101 (2014).
  • (23) A. Zakharova, M. Kapeller, and E. Schöll, Phys. Rev. Lett. 112, 154101 (2014).
  • (24) N. C. Rattenborg, C. J. Amlaner, and S. L. Lima, Neurosci. Biobehav. Rev. 24, 817 (2000).
  • (25) C. R. Laing and C. C. Chow, Neural Computation 13, 1473 (2001).
  • (26) H. Sakaguchi, Phys. Rev. E 73, 031907 (2006).
  • (27) A. E. Filatova, A. E. Hramov, A. A. Koronovskii, and S. Boccaletti, Chaos 18, 023133 (2008).
  • (28) J. C. Gonzalez-Avella, M. G. Cosenza, and M. S. Miguel, Physica A 399, 24 (2014).
  • (29) O. E. Omel’chenko, Nonlinearity 26, 2469 (2013).
  • (30) J. Sieber, O. E. Omel’chenko, and M. Wolfrum, Phys. Rev. Lett. 112, 054102 (2014).
  • (31) M. R. Tinsley, S. Nkomo, and K. Showalter, Nature Physics 8, 662 (2012).
  • (32) S. Nkomo, M. R. Tinsley, and K. Showalter, Phys. Rev. Lett. 110, 244102 (2013).
  • (33) A. M. Hagerstrom, T. E. Murphy, R. Roy, P. Hövel, I. Omelchenko, and E. Schöll, Nature Physics 8, 658 (2012).
  • (34) E. A. Martens, S. Thutupalli, A. Fourrière, and O. Hallatschek, Proc. Nat. Acad. Sciences 110, 10563 (2013).
  • (35) L. Larger, B. Penkovsky, and Y. L. Maistrenko, Phys. Rev. Lett. 111, 054103 (2013).
  • (36) L. Schmidt, K. Schönleber, K. Krischer, and V. Garcia-Morales, Chaos 24, 013102 (2014).
  • (37) M. Wickramasinghe and I. Z. Kiss, PLoS ONE 8, e80586 (2013).
  • (38) C. R. Laing, Phys. Rev. E 81, 066221 (2010).
  • (39) T. W. Ko and G. B. Ermentrout, Phys. Rev. E 78, 016203 (2008).
  • (40) M. Shanahan, Chaos: An Interdisciplinary Journal of Nonlinear Science 20, 013108 (2010).
  • (41) C. R. Laing, K. Rajendran, and I. G. Kevrekidis, Chaos 22, 013132 (2012).
  • (42) N. Yao, Z.-G. Huang, Y.-C. Lai, and Z. Zheng, Scientific Reports 3, 3522 (2013).
  • (43) Y. Zhu, Z. Zheng, and J. Yang, Phys. Rev. E 89, 022914 (2014).
  • (44) A. Yeldesbay, A. Pikovsky, and M. Rosenblum, Phys. Rev. Lett. 112, 144103 (2014).
  • (45) R. FitzHugh, Biophys. J. 1, 445 (1961).
  • (46) J. Nagumo, S. Arimoto, and S. Yoshizawa, Proc. IRE 50, 2061 (1962).
  • (47) V. Vuksanović and P. Hövel, NeuroImage 97, 1 (2014).
  • (48) P. Katsaloulis, D. A. Verganelakis, and A. Provata, Fractals 17, 181 (2009).
  • (49) P. Expert, T. S. Evans, V. D. Blondel, and R. Lambiotte, PNAS 108, 7663 (2011).
  • (50) A. Provata, P. Katsaloulis, and D. A. Verganelakis, Chaos, Solitons & Fractals 45, 174 (2012).
  • (51) E. M. Izhikevich, SIAM J. Appl. Math. 60, 1789 (2000).
  • (52) S. A. Brandstetter, M. A. Dahlem, and E. Schöll, Phil. Trans. R. Soc. A 368, 391 (2010).
  • (53) J. W. Scannell, C. Blakemore, and M. P. Young, J. Neurosci. 15(2), 1463 (1995).
  • (54) B. Hellwig, Biol. Cybern. 82, 111 (2000).
  • (55) C. Holmgren, T. Harkany, B. Svennenfors, and Y. Zilberter, J. Physiol. 551.1, 139 (2003).
  • (56) D. J. Watts and S. H. Strogatz, Nature 393, 440 (1998).
  • (57) M. E. J. Newman and D. J. Watts, Phys. Lett. A 263, 341 (1999).
  • (58) P. Katsaloulis, A. Ghosh, A. C. Philippe, A. Provata, and R. Deriche, Eur. Phys. J. B 85, 1 (2012).
  • (59) B. B. Mandelbrot, The fractal geometry of nature, 3 ed. (W. H. Freeman and Comp., New York, 1983).
  • (60) J. Feder, Fractals (Plenum Press, New York, 1988).
  • (61) A. Vüllings, J. Hizanidis, I. Omelchenko, and P. Hövel, New J. Phys. (in print, 2014), ArXiv:1407.5304.
  • (62) Yu. Maistrenko, A. Vasylenko, O. Sudakov, R. Levchenko, V. Maistrenko, Int. J. Bif. Chaos 24 (8), 1440014 (2014).
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description