Heteroclinic switching between chimeras

TL;DR

This paper analyzes how higher-order interactions in phase oscillator networks lead to metastable chimera states connected by heteroclinic trajectories, explaining switching dynamics in neuronal-like systems.

Contribution

It provides exact results on the emergence of metastable chimeras due to higher-order interactions and elucidates the mechanisms behind switching dynamics in phase oscillator networks.

Findings

Higher-order interactions induce metastable chimeras.

Heteroclinic connections enable switching between states.

Mechanisms are applicable to neuronal network models.

Abstract

Functional oscillator networks, such as neuronal networks in the brain, exhibit switching between metastable states involving many oscillators. We give exact results how such global dynamics can arise in paradigmatic phase oscillator networks: higher-order network interaction gives rise to metastable chimeras - localized frequency synchrony patterns - which are joined by heteroclinic connections. Moreover, we illuminate the mechanisms that underly the switching dynamics in these experimentally accessible networks.

Figures

Equations

$${\theta }_1\cdots {\theta }_{\sigma -1}\text{S}{\theta }_{\sigma +1}\cdots {\theta }_M$$

$${\theta }_1\cdots {\theta }_{\sigma -1}\text{D}{\theta }_{\sigma +1}\cdots {\theta }_M$$

$$g(\vartheta )=\sin (\vartheta +\alpha )+r\sin (2(\vartheta +\alpha ))$$

$$\begin{array}{rl}\mathrm{snp}(\varphi ,\vartheta ;{\theta }_{\tau })& =\cos ({\theta }_{\tau ,1}-{\theta }_{\tau ,2}+\varphi -\vartheta +\alpha )\\ & +\cos ({\theta }_{\tau ,2}-{\theta }_{\tau ,1}+\varphi -\vartheta +\alpha ).\end{array}$$

$$\begin{array}{rl}{\dot{\theta }}_{\sigma ,1}& =\omega +g({\theta }_{\sigma ,2}-{\theta }_{\sigma ,1})-K\mathrm{snp}({\theta }_{\sigma ,2},{\theta }_{\sigma ,1};{\theta }_{\sigma -1})\\ & +K\mathrm{snp}({\theta }_{\sigma ,2},{\theta }_{\sigma ,1};{\theta }_{\sigma +1})=:X_{\sigma ,1}(\theta ),\end{array}$$

$$\begin{array}{rl}{\dot{\theta }}_{\sigma ,2}& =\omega +g({\theta }_{\sigma ,1}-{\theta }_{\sigma ,2})-K\mathrm{snp}({\theta }_{\sigma ,1},{\theta }_{\sigma ,2};{\theta }_{\sigma -1})\\ & +K\mathrm{snp}({\theta }_{\sigma ,1},{\theta }_{\sigma ,2};{\theta }_{\sigma +1})=:X_{\sigma ,2}(\theta ),\end{array}$$

$${\dot{\theta }}_{\sigma ,k}=X_{\sigma ,k}(\theta )+\delta S_{\sigma ,k}(\theta )+\eta W_{\sigma ,k}$$

$$g(\vartheta )=\sin (\vartheta +\alpha )+r\sin (a(\vartheta +\alpha ))$$

$${\dot{\theta }}_{\sigma ,k}=\omega +\frac{1}{N}\sum _{j\ne k}g({\theta }_{\sigma ,j}-{\theta }_{\sigma ,k}+\Delta {\alpha }_{\sigma })$$

$$\Delta {\alpha }_{\sigma }=K((1-R_{\sigma -1}^2)-(1-R_{\sigma +1}^2)),$$

$$\begin{array}{rl}g(\vartheta +\Delta {\alpha }_{\sigma })& =g(\vartheta )+K(R_{\sigma +1}^2-R_{\sigma -1}^2)\cos (\vartheta +\alpha )\\ & +O\left(K^2\right)+O\left(Kr\right).\end{array}$$

$$\mathrm{snps}(\varphi ,\vartheta ;{\theta }_{\tau })=\frac{1}{N^2}\sum _{p,q=1}^N\cos ({\theta }_{\tau ,p}-{\theta }_{\tau ,q}+\varphi -\vartheta +\alpha ).$$

$$R_{\tau }^2\cos \left({\theta }_{\sigma ,j}-{\theta }_{\sigma ,k}+\alpha \right)$$

$${\dot{\theta }}_{\sigma ,k}$$

$$\displaystyle\qquad\quad+K\operatorname{snps}(\theta_{\sigma,j},\theta_{\sigma,k};\theta_{\sigma+1})\Big{)}=:X_{\sigma,k}(\theta)$$

$$\Delta {\alpha }_1$$

$$\Delta {\alpha }_2$$

$$\Delta {\alpha }_3$$

$$\Delta {\alpha }_4$$

$$\dot{\theta}_{\sigma,k}=\omega+g\big{(}\theta_{\sigma,3-k}-\theta_{\sigma,k}+\Delta\alpha_{\sigma}\big{)}+\eta W_{\sigma,k}$$

Full text

Networks of (almost) identical nonlinear oscillators give rise to fascinating collective dynamics where populations of localized oscillators exhibit distinct frequencies and levels of phase synchronization Litwin-Kumar and Doiron (2012); Martens et al. (2013). In neuronal networks, the location of such localized frequency synchrony patterns can encode information Hubel and Wiesel (1959); Tognoli and Kelso (2014); Bick and Martens (2015). Thus, sequential switching between distinct localized dynamics has been associated with neural computation Ashwin and Timme (2005); Rabinovich et al. (2008); Britz et al. (2010); Ashwin et al. (2016); sequential dynamics in the hippocampus in the absence of external input Wilson and McNaughton (1994) are a striking example. Most efforts to understand switching dynamics between localized frequency synchrony patterns rely on averaged models which neglect the contributions of individual oscillators to the network dynamics Kiebel et al. (2009); Bick and Rabinovich (2009); Rabinovich et al. (2015); Horchler et al. (2015); Schaub et al. (2015) or are statistical Wildie and Shanahan (2012). For finite networks, however, the dynamics of individual oscillators cannot be neglected.

In this article we give explicit results for the emergence of switching between synchrony patterns that are characterized by localized frequency synchrony—commonly known as (weak) chimeras Panaggio and Abrams (2015); Schöll (2016)—in phase oscillator networks with higher-order interactions. More precisely, we prove the existence of saddle weak chimeras which are joined by heteroclinic connections; nearby trajectories exhibit sequential switching of localized frequency synchrony. Our results directly relate two distinct dynamic phenomena, heteroclinic switching and chimeras, and thus give a number of insights into the global dynamics of oscillator networks. First, they elucidate how network topology and the functional form of the oscillator coupling facilitate switching dynamics: the heteroclinic structures arise through an interplay of higher harmonics in the phase coupling function and interaction terms which depend on the phase differences of more than two oscillators (nonpairwise interaction). Although such generalized forms of network coupling arise naturally in phase reductions of generically coupled limit cycle oscillator networks Ashwin and Rodrigues (2016), they are neglected in classical Kuramoto-type networks Acebrón et al. (2005); Rodrigues et al. (2016). Hence, our results emphasize how higher-order interaction terms can shape the phase dynamics of many physical systems, from oscillator networks Rosenblum and Pikovsky (2007); Bick et al. (2016); Tanaka and Aoyagi (2011) to ecological systems Levine et al. (2017). Second, switching between metastable chimeras is an explicit dynamical mechanism how networks of neural oscillators may encode sequential information and give rise to dynamics similar to hippocampal replay. Third, we provide a theoretical foundation to understand self-organized switching between chimeras that was recently observed in numerical simulations Haugland et al. (2015); Maistrenko et al. (2017). Finally, relating heteroclinic switching and chimeras opens up a range of questions; for example, whether any given itinerary can be realized as a heteroclinic structure between chimeras.

In the following, we consider networks of $M$ populations of $N$ phase oscillators. Let $\theta_{\sigma,k}\in\mathbf{T}:=\mathbb{R}/2\pi\mathbb{Z}$ denote the phase of oscillator $k$ in population $\sigma$. Write $\theta=(\theta_{1},\dotsc,\theta_{M})\in\mathbf{T}^{MN}$ where $\theta_{\sigma}=(\theta_{\sigma,1},\dotsc,\theta_{\sigma,N})\in\mathbf{T}^{N}$ is the state of population $\sigma$. The set \textrm{S}:=\big{\{}(\phi_{1},\dotsc,\phi_{N})\in\mathbf{T}^{N}\,\big{|}\;\phi_{k}=\phi_{k+1}\big{\}} corresponds to phase synchrony and \textrm{D}:=\big{\{}(\phi_{1},\dotsc,\phi_{N})\in\mathbf{T}^{N}\,\big{|}\;\phi_{k+1}=\phi_{k}+\frac{2\pi}{N}\big{\}} denotes the splay phase where phases are distributed uniformly on the circle. Following Martens (2010) we use the shorthand notation

[TABLE]

to indicate that population $\sigma$ is phase synchronized or in splay phase. Hence, $\textrm{S}\dotsb\textrm{S}$ ($M$ times) is the set of cluster states and $\textrm{D}\dotsb\textrm{D}$ is the set where all populations are in splay phase. Given a dynamical system on $\mathbf{T}^{MN}$ and a trajectory $\theta(t)$ with initial condition $\theta(0)=\theta^{0}$, define the asymptotic average angular frequency $\Omega_{\sigma,k}(\theta^{0}):=\lim_{t\to\infty}\frac{1}{t}\theta_{\sigma,k}(t)$. The characterizing feature of a weak chimera as an invariant set $A\subset\mathbf{T}^{MN}$ is localized frequency synchrony: for all $\theta^{0}\in A$ we have oscillators $(\sigma,k)$, $(\tau,j)$, $(\rho,\ell)$ such that $\Omega_{\sigma,k}(\theta^{0})=\Omega_{\tau,j}(\theta^{0})\neq\Omega_{\rho,\ell}(\theta^{0})$; see also Ashwin and Burylko (2015); Bick and Ashwin (2016); Bick (2017).

*Heteroclinic cycles in small networks.—*Consider a network of $M=3$ populations of $N=2$ identical phase oscillators where the interaction within populations is pairwise and determined by the coupling function

[TABLE]

parametrized by $\alpha,r\in\mathbb{R}$, whereas different populations interact at coupling strength $K$ through the sinusoidal nonpairwise interaction function

[TABLE]

More specifically, the dynamics of population $\sigma\in\left\{1,2,3\right\}$ is given by

[TABLE]

where $\omega$ is the oscillators’ intrinsic frequency 111Note that we can set $\omega$ to any value without loss of generality by going in a suitable co-rotating reference frame. In the figures we set $\omega=-\sum_{j=1}^{N-1}g(2\pi j/N)$ so that for $K=0$ the splay configuration appears stationary. and indices are taken modulo $M$.

The coupling induces symmetries of the oscillator network. For each of the $M$ populations, let $\mathbf{T}$ act by shifting all phases of that population by a common constant and let the symmetric group $\mathbf{S}_{N}$ permute its $N$ oscillators. Suppose that $\mathbb{Z}_{M}:=\mathbb{Z}/M\mathbb{Z}$ permutes populations cyclically. The equations of motion (4) are invariant under the group of transformations $(\mathbf{S}_{N}\times\mathbf{T})^{M}\rtimes\mathbb{Z}_{M}$ of $\mathbf{T}^{MN}$. The semidirect product “$\rtimes$” indicates that actions do not necessarily commute Ashwin and Swift (1992). These symmetries induce invariant subspaces Golubitsky and Stewart (2002): in particular SSS, DDD as well as DSS, DDS and their images under permutations of populations are dynamically invariant.

We can now give conditions for (4) to have the heteroclinic cycle depicted in Fig. 1(a) between saddle weak chimeras DSS, DDS and their symmetric counterparts. Because of symmetry, it suffices to consider DSS, DDS. We proceed in three steps. First, we want DSS, DDS to be weak chimeras. Second, we give conditions for the invariant sets to be saddles. Third, we show that they are connected by heteroclinic orbits. Here we focus on the case of $\alpha=\frac{\pi}{2}$ and refer to Bick and et al (2017) for more generality and a proof that there is in fact an open set of parameters $(\alpha,K,r)$ for which this heteroclinic cycle between weak chimeras exists.

First, for DSS, DDS to be weak chimeras, we calculate the frequencies $\Omega_{\sigma,k}$ for (4). For $K=0$ we have $\Omega_{1,k}(\theta^{0})=\omega+1$ for $\theta^{0}\in\textrm{S}\theta_{2}\theta_{3}$ and $\Omega_{2,k}(\theta^{0})=\omega-1$ for $\theta^{0}\in\theta_{1}\textrm{D}\theta_{3}$. In other words, without coupling between populations, the frequency difference between a synchronized and an anti-phase population is $\left|(\omega+1)-(\omega-1)\right|=2$. With coupling, $K>0$, the maximal change in frequency difference is proportional to $K$. Specifically, using the triangle inequality in (4) yields that ${\Omega_{1,k}(\theta^{0})\neq\Omega_{2,k}(\theta^{0})}$ for $\theta^{0}\in\textrm{S}\textrm{D}\theta_{3}$ if $2-8K>0$. At the same time, ${\Omega_{\sigma,k}(\theta^{0})=\Omega_{\sigma,j}(\theta^{0})}$ for all $\theta^{0}\in\mathbf{T}^{MN}$ with $\theta^{0}_{\sigma}\in\textrm{S},\textrm{D}$. Hence, DSS, DDS are weak chimeras for (4) on $\mathbf{T}^{MN}$ if $2K<\frac{1}{2}$.

Second, we need ${\textrm{D}\textrm{S}\textrm{S}},{\textrm{D}\textrm{D}\textrm{S}}$ to be saddle invariant sets. Reduce the phase-shift symmetries by rewriting (4) in terms of phase differences $\psi_{\sigma,k}:=\theta_{\sigma,k+1}-\theta_{\sigma,1}$, $k=1,\dotsc,N-1$. (Consequently, we may replace all $\theta$ by the phase differences $\psi$ in (1).) Since $N=2$ here, $\psi_{\sigma}=\psi_{\sigma,1}$ determines the state of population $\sigma$ and the effective dynamics of (4) are three-dimensional. In the reduced system ${\textrm{D}\textrm{S}\textrm{S}}=(\pi,0,0)$, ${\textrm{D}\textrm{D}\textrm{S}}=(\pi,\pi,0)$ are equilibria. Linearizing at DSS yields eigenvalues $\lambda^{\textrm{D}\textrm{S}\textrm{S}}_{1}=4r$, $\lambda^{\textrm{D}\textrm{S}\textrm{S}}_{2}=8K+4r$, $\lambda^{\textrm{D}\textrm{S}\textrm{S}}_{3}=-8K+4r$ that correspond to linear stability of the first, second, and third population, respectively. Similarly, for DDS we obtain the eigenvalues $\lambda^{\textrm{D}\textrm{D}\textrm{S}}_{1}=8K+4r$, $\lambda^{\textrm{D}\textrm{D}\textrm{S}}_{2}=-8K+4r$, $\lambda^{\textrm{D}\textrm{D}\textrm{S}}_{3}=4r$. Observe that if $0<-r<2K$ we have $\lambda^{\textrm{D}\textrm{S}\textrm{S}}_{1}=\lambda^{\textrm{D}\textrm{D}\textrm{S}}_{3}<0$, $\lambda^{\textrm{D}\textrm{S}\textrm{S}}_{2}=\lambda^{\textrm{D}\textrm{D}\textrm{S}}_{1}>0$, $\lambda^{\textrm{D}\textrm{S}\textrm{S}}_{3}=\lambda^{\textrm{D}\textrm{D}\textrm{S}}_{2}<0$ and thus ${\textrm{D}\textrm{S}\textrm{S}},{\textrm{D}\textrm{D}\textrm{S}}$ are saddle invariant sets with two-dimensional stable and one-dimensional unstable manifolds.

Third, we obtain conditions for heteroclinic connections between ${\textrm{D}\textrm{S}\textrm{S}},{\textrm{D}\textrm{D}\textrm{S}}$ given their stability above. Observe that $\lambda^{\textrm{D}\textrm{S}\textrm{S}}_{2}>0$, $\lambda^{\textrm{D}\textrm{D}\textrm{S}}_{2}<0$ implies that the unstable manifold of DSS and the stable manifold of DDS both intersect the invariant subspace $\textrm{D}\psi_{2}\textrm{S}$ on which the dynamics reduce to $\dot{\psi}_{2}=\sin(\psi_{2})(8K-4r\cos(\psi_{2}))$. Thus, if $-r<2K$ there are no equilibria other than $\psi_{2}\in\left\{0,\pi\right\}$ (these are DSS and DDS) in $\textrm{D}\psi_{2}\textrm{S}$ and we have a heteroclinic connection. Indeed, we get the same condition for there to be no additional equilibria in $\psi_{1}\textrm{D}\textrm{S}$. To summarize, for $\alpha=\frac{\pi}{2}$ the heteroclinic cycle sketched in Fig. 1(a) exists if $0<-r<2K<\frac{1}{2}$. Moreover, one can show by evaluating the saddle values that for $K<-r$ the cycle is expected to attract nearby initial conditions Bick and et al (2017).

The switching dynamics between weak chimeras persists when the particular nonpairwise coupling scheme of (4) is broken. With noise given by a Wiener process $W_{\sigma,k}$ (Brownian motion) and a symmetry breaking coupling term $S_{\sigma,k}(\theta)=\Delta\omega_{\sigma,k}+\frac{1}{MN}\sum_{\tau=1}^{M}\sum_{j=1}^{N}\sin(\theta_{\tau,j}-\theta_{\sigma,k})$ with normally distributed frequency deviations $\Delta\omega_{\sigma,k}$ (mean zero and variance one), we integrated the system

[TABLE]

numerically in XPP Ermentrout (2002) where $X_{\sigma,k}$ as in (4). For $\eta>0$, $\delta=0$ we obtain heteroclinic switching where transition times scale with the noise amplitude $\eta$ as expected Stone and Holmes (1990); cf. Fig. 1(b). Setting $\delta>0$ breaks all symmetries to a single phase-shift symmetry acting as a common phase shift for all oscillators. Although this breaks the invariant subspaces containing the heteroclinic connections, we still obtain sequential dynamics prescribed by the heteroclinic network as shown in Fig. 1(c).

*Order parameter dependent coupling induces switching.—*The dynamical mechanism which leads to heteroclinic cycles in (4) can be best understood if the oscillator network is seen as individual populations coupled through their mean fields. Let $i=\sqrt{-1}$. The absolute value of the Kuramoto order parameter R_{\sigma}:=R(\theta_{\sigma})=\big{|}\frac{1}{N}\sum_{j=1}^{N}\exp(i\theta_{\sigma,j})\big{|} gives information about synchronization: $\theta_{\sigma}\in\textrm{S}$ iff $R(\theta_{\sigma})=1$ and $\theta_{\sigma}\in\textrm{D}$ implies $R(\theta_{\sigma})=0$. For $a\in\mathbb{N}$ let

[TABLE]

generalize the coupling function (2). Now consider a system of $M$ populations of $N$ phase oscillators each where the dynamics of oscillator $k$ in population $\sigma$ are given by

[TABLE]

and $\Delta\alpha_{\sigma}$ modulates the phase-shift $\alpha$ of the coupling function (6). If $r=0$ then either full synchrony S or the phase configurations with $R_{\sigma}=0$ are globally attracting for (7) depending on the value of $\alpha+\Delta\alpha_{\sigma}$ 222Except for some set of initial conditions of zero Lebesgue measure.. In particular, the global attractors swap stability at $\alpha+\Delta\alpha_{\sigma}=\pm{\frac{\pi}{2}}$. Hence, for $r=0$ and $\alpha\approx\frac{\pi}{2}$ the order parameter-dependent modulation of $\Delta\alpha_{\sigma}$ by

[TABLE]

$0<K\lessapprox\frac{\pi}{2}$, yields a mechanism for sequential synchronization: If population $\sigma-1$ is synchronized ($R_{\sigma-1}=1$) and population $\sigma+1$ is in splay phase ($R_{\sigma+1}=0$) then S is asymptotically stable for population $\sigma$. Conversely, if $R_{\sigma+1}=1$ and $R_{\sigma-1}=0$ then $R_{\sigma}=0$ is asymptotically stable for population $\sigma$. Whereas the system is degenerate for $R_{\sigma-1}=R_{\sigma+1}$ if $\alpha=\frac{\pi}{2}$ and $r=0$, an appropriate choice of $a$ and $r\neq 0$ to induce bistability of S and D will resolve the degeneracy below.

A network with nonpairwise coupling approximates the system (7) with state-dependent phase shift (8). We have

[TABLE]

Generalizing (3), define the sinusoidal nonpairwise scaled interaction function

[TABLE]

Note that $R_{\tau}^{2}=\frac{1}{N^{2}}\sum_{p,q=1}^{N}\cos(\theta_{\tau,p}-\theta_{\tau,q})$ which implies

[TABLE]

Substituting (9) and (10) into (7) and dropping the O\big{(}K^{2}\big{)}, $O\!\left(Kr\right)$ terms yields the phase dynamics

[TABLE]

as an approximation of (7). Note that for $M=3$, $N=2$, the system (4) with coupling function (2) is—up to rescaling of $K$ and time—exactly this approximation (Heteroclinic switching between chimeras) with (6) and harmonic $a=2$ that yields hyperbolic saddles.

*Switching dynamics for larger networks.—*The derivation of the nonpairwise coupling suggests a general mechanism to obtain switching dynamics in systems with population sizes $N>2$. Indeed, we obtain sequential switching dynamics for example for $M=3$, $N=11$: integrating (5) with $X_{\sigma,k}$ as in (Heteroclinic switching between chimeras) yields sequential switching even when the system symmetries are broken, $\delta>0$; cf. Fig. 2. Note that the transitions now take place along high-dimensional invariant subspaces.

*From heteroclinic cycles to networks.—*Generalizing the order parameter-dependent coupling (8) for the dynamics (7) leads to switching similar to those observed for the Kirk–Silber heteroclinic network Kirk and Silber (1994) which contains more than one cycle; cf. Fig. 3(a). Similar to (8), set

[TABLE]

Consider $M=4$ populations of $N=2$ oscillators where oscillator $(\sigma,k)$ evolves according to

[TABLE]

with coupling function $g$ as in (2). The $\Delta\alpha_{\sigma}$ given by (12) are now chosen to allow for switching from SDSS to either SSDS or SSSD: if population 2 is desynchronized, $R_{2}=0$, and all other populations are synchronized, $R_{\sigma}=1,\sigma\neq 2$ then D will be attracting for both populations 3 and 4 (in the limiting case $r=0$). Fig. 3(b) shows noise-induced switching in (13). A full analysis of this system (and its nonpairwise approximation) is beyond the scope of this article.

*Discussion—*Phase oscillator networks with nonpairwise coupling have surprisingly rich dynamics Rosenblum and Pikovsky (2007); Tanaka and Aoyagi (2011); Ashwin and Rodrigues (2016); Bick et al. (2016); here, nonpairwise interaction allows to show the existence of heteroclinic connections between weak chimeras. Here nonpairwise coupling arises through a bifurcation parameter that depends on local order parameters of different populations. By contrast, the dynamics of a network with a bifurcation parameter depending on the global order parameter has been studied in their own right Burylko and Pikovsky (2011) and exploited for applications Sieber et al. (2014). In contrast to sequential switching of phase synchrony for nonidentical oscillators Komarov and Pikovsky (2011), here we observe switching of localized frequency synchrony in a network of indistinguishable phase oscillators (the symmetry action is transitive). Moreover, since the system is close to bifurcation for small $K$, small perturbations to the vector field allow for going from one switching sequence to another.

Our results open up a range of questions relating both chimeras and heteroclinic networks. Are there heteroclinic cycles between saddle weak chimeras with chaotic dynamics Bick and Ashwin (2016)? Is it possible to realize any heteroclinic network in a phase oscillator network where the saddles are weak chimeras, see also Ashwin and Postlethwaite (2013); Field (2015)? How do the dynamics of (13) relate to results obtained for the Kirk–Silber network Castro and Lohse (2016)?

Heteroclinic switching between localized frequency synchrony patterns is of direct relevance for real-world systems. On the one hand, note that the small networks considered here are accessible for experimental realizations: weak chimeras have recently been observed in electrochemical systems Bick et al. (2017) with linear and quadratic interactions interactions Kori et al. (2008). Thus, we are interested in whether switching of localized frequency synchrony is observed these experimental setups. On the other hand, sequential switching of localized frequency synchrony may be an important aspect of functional dynamics in networks of neurons. Our results elucidate the features of network interaction (e.g., symmetries and nonpairwise interactions) and the dynamical mechanisms that facilitate switching dynamics. Thus, our insights may open up ways to restore and control functional dynamics, for example, if the network becomes pathologically synchronized.

*Acknowledgements—*The author would like to thank M. Field, E. A. Martens, O. Omel’chenko, T. Pereira, M. Rabinovich, M. Wolfrum, and in particular P. Ashwin for many helpful discussions. This work has received funding from the People Programme (Marie Curie Actions) of the European Union’s Seventh Framework Programme (FP7/2007–2013) under REA grant agreement no. 626111.