# Robustness of chimera states for coupled FitzHugh-Nagumo oscillators

## Abstract

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.

###### pacs:

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 (1); (2). 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 (3); (4). 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 (5); (6); (7); (8); (9); (10); (11); (12). 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 (13) and time-continuous chaotic models (14), and neural systems (15); (16). Moreover, chimera states were found in systems with higher spatial dimensions (9); (17); (18); (19). 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 (20); (15), as well as amplitude-mediated (21); (22), and pure amplitude chimera and chimera death states (23) were discovered.

Potential applications of chimera states in nature include the phenomenon of unihemispheric sleep in birds and dolphins (24), bump states in neural systems (25); (26), in power grids (27), or in social systems (28). Many works considering chimera states were mostly based on numerical results. A deeper bifurcation analysis (29) and even a possibility to control chimera states (30) were obtained only recently.

Ten years after the first numerical observation, the experimental verification of chimera states was demonstrated in chemical (31); (32) and optical (33) systems. Further experiments involved mechanical (34), electronic (35), and electrochemical (36); (37) 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 (38), or when the topology is not regular (39); (40); (41); (42); (43) or even global (36); (44). 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 (45); (46) studied in Ref. (15). 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 (47). 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 (48); (49); (50). 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 (15):

(1a) | ||||

(1b) |

where and are the activator and inhibitor variables, respectively (45); (46), 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:

(2) |

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. (15). There, applying a phase-reduction technique (51), 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.

## 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 (52). 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. (15).

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 (15). 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.

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 (52). 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 (53); (54); (55). 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

(3a) | ||||

(3b) |

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) (15).

The introduction of shortcuts shares similarities to a well-known small-world network (56) 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 (57).

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

(4) |

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(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. (15) 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 (56), 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) (8). 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 (15), we consider the following set of parameters: , , , , , and , which has been studied earlier (15).

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:

(5) |

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.

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.

(6) |

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 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. (15). 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.

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:

(7) |

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:

(8) |

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 .

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 .

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 (48); (49); (58). 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 (59); (60). 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,

(9) |

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 (60).

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.

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 (61).

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. (15), 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 (62).

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.

## 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.

###### Acknowledgements.

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.### References

- A. Pikovsky, M. G. Rosenblum, and J. Kurths, Synchronization, A Universal Concept in Nonlinear Sciences (Cambridge University Press, Cambridge, 2001).
- S. Boccaletti, V. Latora, Y. Moreno, M. Chavez, and D. U. Hwang, Phys. Rep. 424, 175 (2006).
- Y. Kuramoto and D. Battogtokh, Nonlin. Phen. in Complex Sys. 5, 380 (2002).
- D. M. Abrams and S. H. Strogatz, Phys. Rev. Lett. 93, 174102 (2004).
- C. R. Laing, Physica D 238, 1569 (2009).
- A. E. Motter, Nature Physics 6, 164 (2010).
- E. A. Martens, C. R. Laing, and S. H. Strogatz, Phys. Rev. Lett. 104, 044101 (2010).
- O. E. Omel’chenko, M. Wolfrum, and Y. L. Maistrenko, Phys. Rev. E 81, 065201(R) (2010).
- O. E. Omel’chenko, M. Wolfrum, S. Yanchuk, Y. L. Maistrenko, and O. Sudakov, Phys. Rev. E 85, 036210 (2012).
- E. A. Martens, Chaos 20, 043122 (2010).
- M. Wolfrum, O. E. Omel’chenko, S. Yanchuk, and Y. L. Maistrenko, Chaos 21, 013112 (2011).
- T. Bountis, V. Kanas, J. Hizanidis, and A. Bezerianos, Eur. Phys. J. Special Topic 223, 721 (2014).
- I. Omelchenko, Y. L. Maistrenko, P. Hövel, and E. Schöll, Phys. Rev. Lett. 106, 234102 (2011).
- I. Omelchenko, B. Riemenschneider, P. Hövel, Y. L. Maistrenko, and E. Schöll, Phys. Rev. E 85, 026212 (2012).
- I. Omelchenko, O. E. Omel’chenko, P. Hövel, and E. Schöll, Phys. Rev. Lett. 110, 224101 (2013).
- J. Hizanidis, V. Kanas, A. Bezerianos, and T. Bountis, Int. J. Bifurcat. Chaos 24, 1450030 (2014).
- S.-I. Shima and Y. Kuramoto, Phys. Rev. E 69, 036213 (2004).
- M. J. Panaggio and D. M. Abrams, Phys. Rev. Lett. 110, 094102 (2013).
- M. J. Panaggio and D. M. Abrams, arXiv:1403.6204 (2014).
- G. C. Sethia, A. Sen, and F. M. Atay, Phys. Rev. Lett. 100, 144102 (2008).
- G. C. Sethia, A. Sen, and G. L. Johnston, Phys. Rev. E 88, 042917 (2013).
- G. C. Sethia and A. Sen, Phys. Rev. Lett. 112, 144101 (2014).
- A. Zakharova, M. Kapeller, and E. Schöll, Phys. Rev. Lett. 112, 154101 (2014).
- N. C. Rattenborg, C. J. Amlaner, and S. L. Lima, Neurosci. Biobehav. Rev. 24, 817 (2000).
- C. R. Laing and C. C. Chow, Neural Computation 13, 1473 (2001).
- H. Sakaguchi, Phys. Rev. E 73, 031907 (2006).
- A. E. Filatova, A. E. Hramov, A. A. Koronovskii, and S. Boccaletti, Chaos 18, 023133 (2008).
- J. C. Gonzalez-Avella, M. G. Cosenza, and M. S. Miguel, Physica A 399, 24 (2014).
- O. E. Omel’chenko, Nonlinearity 26, 2469 (2013).
- J. Sieber, O. E. Omel’chenko, and M. Wolfrum, Phys. Rev. Lett. 112, 054102 (2014).
- M. R. Tinsley, S. Nkomo, and K. Showalter, Nature Physics 8, 662 (2012).
- S. Nkomo, M. R. Tinsley, and K. Showalter, Phys. Rev. Lett. 110, 244102 (2013).
- A. M. Hagerstrom, T. E. Murphy, R. Roy, P. Hövel, I. Omelchenko, and E. Schöll, Nature Physics 8, 658 (2012).
- E. A. Martens, S. Thutupalli, A. Fourrière, and O. Hallatschek, Proc. Nat. Acad. Sciences 110, 10563 (2013).
- L. Larger, B. Penkovsky, and Y. L. Maistrenko, Phys. Rev. Lett. 111, 054103 (2013).
- L. Schmidt, K. Schönleber, K. Krischer, and V. Garcia-Morales, Chaos 24, 013102 (2014).
- M. Wickramasinghe and I. Z. Kiss, PLoS ONE 8, e80586 (2013).
- C. R. Laing, Phys. Rev. E 81, 066221 (2010).
- T. W. Ko and G. B. Ermentrout, Phys. Rev. E 78, 016203 (2008).
- M. Shanahan, Chaos: An Interdisciplinary Journal of Nonlinear Science 20, 013108 (2010).
- C. R. Laing, K. Rajendran, and I. G. Kevrekidis, Chaos 22, 013132 (2012).
- N. Yao, Z.-G. Huang, Y.-C. Lai, and Z. Zheng, Scientific Reports 3, 3522 (2013).
- Y. Zhu, Z. Zheng, and J. Yang, Phys. Rev. E 89, 022914 (2014).
- A. Yeldesbay, A. Pikovsky, and M. Rosenblum, Phys. Rev. Lett. 112, 144103 (2014).
- R. FitzHugh, Biophys. J. 1, 445 (1961).
- J. Nagumo, S. Arimoto, and S. Yoshizawa, Proc. IRE 50, 2061 (1962).
- V. Vuksanović and P. Hövel, NeuroImage 97, 1 (2014).
- P. Katsaloulis, D. A. Verganelakis, and A. Provata, Fractals 17, 181 (2009).
- P. Expert, T. S. Evans, V. D. Blondel, and R. Lambiotte, PNAS 108, 7663 (2011).
- A. Provata, P. Katsaloulis, and D. A. Verganelakis, Chaos, Solitons & Fractals 45, 174 (2012).
- E. M. Izhikevich, SIAM J. Appl. Math. 60, 1789 (2000).
- S. A. Brandstetter, M. A. Dahlem, and E. Schöll, Phil. Trans. R. Soc. A 368, 391 (2010).
- J. W. Scannell, C. Blakemore, and M. P. Young, J. Neurosci. 15(2), 1463 (1995).
- B. Hellwig, Biol. Cybern. 82, 111 (2000).
- C. Holmgren, T. Harkany, B. Svennenfors, and Y. Zilberter, J. Physiol. 551.1, 139 (2003).
- D. J. Watts and S. H. Strogatz, Nature 393, 440 (1998).
- M. E. J. Newman and D. J. Watts, Phys. Lett. A 263, 341 (1999).
- P. Katsaloulis, A. Ghosh, A. C. Philippe, A. Provata, and R. Deriche, Eur. Phys. J. B 85, 1 (2012).
- B. B. Mandelbrot, The fractal geometry of nature, 3 ed. (W. H. Freeman and Comp., New York, 1983).
- J. Feder, Fractals (Plenum Press, New York, 1988).
- A. Vüllings, J. Hizanidis, I. Omelchenko, and P. Hövel, New J. Phys. (in print, 2014), ArXiv:1407.5304.
- Yu. Maistrenko, A. Vasylenko, O. Sudakov, R. Levchenko, V. Maistrenko, Int. J. Bif. Chaos 24 (8), 1440014 (2014).