The effect of time-correlated noise on the Kuramoto model studied via the unified colored noise approximation

TL;DR

This paper investigates how finite-time correlated noise influences synchronization in the Kuramoto model, developing an approximation to predict critical coupling changes and revealing that correlated noise promotes synchronization.

Contribution

The study introduces a unified colored noise approximation to analyze the impact of time-correlated noise on the Kuramoto model's synchronization threshold.

Findings

Time-correlated noise decreases the critical coupling $k_c$.

The approximation qualitatively matches numerical results near $k_c$.

Order parameter curves scale on a master curve regardless of $ au$.

Abstract

Many natural and social phenomena are characterized by synchronization. The Kuramoto model, taking into account the basic ingredients for observing synchronized states, allows to study mathematically synchronization in a simplified but nontrivial picture. Here we study how a noise that is correlated on a finite time-scale $\tau$ impacts the ability of the Kuramoto model to achieve synchronization. We develop an approximated theory that allows to compute the critical coupling constant $k_c$ as a function of the correlation time $\tau$. We obtain that that $k_c(\tau)$ decreases as $\tau$ increases indicating that time-correlated noise promotes synchronization. Moreover, we show that theory describes qualitatively well the degree of synchronization near $k_c$ obtained numerically. Finally, we show that, independently on the value of $\tau$, the curves of the order parameter versus $k$ scale on the same master curve even at values of $k$ very far from $k_c$.

Figures

Equations

$${\dot{\varphi }}_i={\omega }_i-\frac{k}{N}\sum _{j=1}^N\sin ({\varphi }_i-{\varphi }_j)+{\eta }_i$$

$${\dot{\eta }}_i=-{\tau }^{-1}{\eta }_i+D^{1/2}{\tau }^{-1}{\zeta }_i$$

$$\sigma \exp (i{\varphi }_0)=\frac{1}{N}\sum _{j=1}^N\exp (i{\varphi }_j)$$

$$\dot{\psi }=\omega -k\sigma \sin \psi +\eta$$

$$\ddot{\psi }+\dot{\psi }[{\tau }^{-1}-f^{\prime }(\psi )]-f(\psi )/\tau =D^{1/2}\zeta /\tau$$

$$\ddot{\psi }+\dot{\psi }[{\tau }^{-1/2}-{\tau }^{1/2}{f}^{\prime }(\psi )]-f(\psi )=D^{1/2}\zeta /{\tau }^{1/4}.$$

$$\dot{\psi }=[1-\tau f^{\prime }(\psi ){]}^{-1}f(\psi )+[1-\tau f^{\prime }(\psi ){]}^{-1}{D}^{1/2}\zeta .$$

$$\dot{n}=-{\partial }_{\psi }\left[f\gamma n-D\gamma {\partial }_{\psi }(\gamma n)\right]$$

$$f\gamma n-D\gamma {\partial }_{\psi }(\gamma n)=\mathrm{const}.$$

$$n(\psi ,\omega )$$

$$F(\psi ,\omega )$$

$$G(\psi )$$

$$\sigma ={\int }_{-\infty }^{\infty }d\omega g(\omega ){\int }_0^{2\pi }d\psi n(\psi ,\omega )\exp (i\psi ).$$

$$\sigma -k\sigma {\int }_{-\infty }^{\infty }d\omega Q(\omega )+k^3{\sigma }^3{\int }_{-\infty }^{\infty }d\omega P(\omega )=0$$

$$Q(\omega )=\frac{1}{2}g(\omega )\left(\frac{D}{D^2+{\omega }^2}+\tau \right)$$

$$P(\omega )=\frac{g(\omega )\left[2D^3-4D{\omega }^2+{\tau }^2\left(6D^5+6D^3{\omega }^2\right)+\tau \left(4D^4+5D^2{\omega }^2+{\omega }^4\right)\right]}{8{\left(D^2+{\omega }^2\right)}^2\left(4D^2+{\omega }^2\right)}.$$

$$\sigma ={\left[\frac{k\int d\omega Q(\omega )-1}{k^3\int d\omega P(\omega )}\right]}^{1/2}$$

$$k_c(\tau )=\frac{2}{\tau +D{\int }_{-\infty }^{\infty }d\omega \frac{g(\omega )}{D^2+{\omega }^2}}$$

$$\sigma =\mathcal{A}(\tau ){\left[\frac{\kappa }{(1+\kappa {)}^3}\right]}^{1/2}$$

$$\sigma \approx \mathcal{A}(\tau ){\kappa }^{1/2}$$

$$g(\omega )=\frac{1}{\pi }\frac{\lambda }{{\omega }^2+{\lambda }^2}$$

$$k_c(\tau )=\frac{2(D+\lambda )}{1+\tau (D+\lambda )}$$

$$k_c(\tau \to 0)=2(D+\lambda )$$

$$k_c(\tau \to \infty )=2\lambda$$

$$k_c(\tau )=k_c(\infty )+\frac{k_c(0)-k_c(\infty )}{1+\alpha \tau }$$

Full text

I Introduction

The firing of neurons in the visual cortex Eckhorn et al. (1988); Gray et al. (1989), the frequency locking in Josephson arrays Wiesenfeld et al. (1998), and the flashing in large groups of fireflies Buck and Buck (1968) are just a few examples, ranging from biology to physics, where synchronization plays a fundamental role Strogatz (2003). Winfree realized that nonlinear interactions are a key ingredient for synchronization phenomena that happen above a threshold value of the coupling constant Winfree (1967, 2001). The Kuramoto model (KM), originally introduced by Kuramoto in 1975, contains all the basic ingredients needed to observe synchronization Kuramoto (1975). In the KM a population of oscillators, each of them characterized by its natural frequency, are coupled globally through a non-linear interaction. The model is thus general enough to describe diverse situations, contains a small number of control parameters, and, because of the mean-field coupling, is analytically tractable.

The KM undergoes a continuous phase transition towards a partially synchronous state Acebrón et al. (2005) as the control parameter of the transition (the coupling constant $k$ between the oscillators) is increased above a critical value. When white noise is taken into account, this acts as a random perturbation that prevents the system from reaching a perfectly synchronized state and makes the critical value of the coupling constant noise-dependent, i. e., $k_{c}=k_{c}(T)$, with the “temperature” $T$ representing the strength of the noise. However, in many biological systems Nozaki et al. (1999), the external noise is not characterized by a flat spectrum. A noise source that is exponentially correlated on a finite time scale $\tau$, is a step towards a more realistic description of the system. In this paper, we study the KM in the presence of a Gaussian exponentially correlated noise Tönjes (2010); Bag et al. (2007); Moro and Sánchez (1998). Following the standard terminology in the field of noise-driven dynamical systems Hänggi (1995), we refer to this type of noise as Gaussian colored noise (the adjective “colored” refers to the fact the power spectrum of the noise is not flat/white as the noise spectrum of standard Brownian motion).

We show that, within the Unified Colored Noise Approximation (UCNA) Jung and Hänggi (1987); Hänggi (1995); Maggi et al. (2015), it is possible to solve the colored-noise driven KM for the critical coupling constant $k_{c}$ which becomes a decreasing function of $\tau$. This is in agreement with early numerical and theoretical works Bag et al. (2007); Moro and Sánchez (1998) and we further test it by computer simulations which, however, reveal significant quantitative deviations from the UCNA predictions at large values of $\tau$. Moreover, we compute analytically the critical behavior of the order parameter that is found to grow as $[k-k_{c}(\tau)]^{1/2}$, i.e. with the same mean-field exponent of the KM in presence of thermal noise. Finally, we show that, intriguingly, the synchronization order parameter of the KM displays a universal behavior that holds well even when $k$ is very far from $k_{c}$, at all values of $\tau$ simulated.

II Theory

We consider the KM defined by the set of equations Acebrón et al. (2005)

[TABLE]

where $\phi_{i}\in[0,2\pi]$ is the phase of the $i-$th oscillator ($i=1,...,N$), $k$ is the coupling strength and the $\omega_{i}$ is the (random) natural frequency distributed among the oscillators according to the probability function $g(\omega)$ which is assumed to be symmetric and with zero mean. The $\eta_{i}$ are a set colored noise sources evolving according to the stochastic equations:

[TABLE]

where $\tau$ is the correlation time of the noise, $D$ is the noise amplitude and the $\zeta_{i}$ are a set of uncorrelated standard white noise sources: ${\langle\zeta_{i}(0)\zeta_{j}(t)\rangle=2\delta_{ij}\,\delta(t)}$.

Following Ref. Sakaguchi (1988) we define the (complex) order parameter $\sigma\,\exp(i\,\phi_{0})$ as

[TABLE]

where the amplitude $\sigma$ of the order parameter can be interpreted as the degree of synchronization of the oscillators population ($0\leq\sigma\leq 1$). By using Eq. (3) we can rewrite Eq. (1) as

[TABLE]

where $\psi=\phi_{i}-\phi_{0}$. Since from this point we focus on one single oscillator (the $i$-th oscillator), we have dropped the index $i$ for simplicity (setting ${\eta_{i}=\eta}$ and ${\omega_{i}=\omega}$) as in Ref. Sakaguchi (1988).

It is well known that Hänggi (1995), for non-linear colored noise-driven systems as the one described by Eq. (4), further progress can be made only by using some approximation scheme. We now follow ref. Hänggi (1995) to illustrate how the UCNA applies to our problem. We first write Eqs. (4) and (2) as one single equation of higher order

[TABLE]

where $f(\psi)=\omega-k\,\sigma\,\sin\psi$ is the “mean-field force” appearing in Eq. (4) and the prime indicates differentiation with respect to $\psi$. Eq. (13) shows that the dynamics of the phase $\psi$ is governed by the field $f$, the white noise $\zeta$ and the non-homogeneous friction $[\tau^{-1}-f(\psi)]$. It is now convenient to introduce the time-scale $s=\tau^{-1/2}t$, so that Eq. (5) becomes

[TABLE]

The UCNA consists in eliminating adiabatically the term $\ddot{\psi}$ in Eq. (6) which is justified whenever $\tau$ is small. If $\tau$ is large the validity of the approximation is restricted to the case where $f^{\prime}(\psi)$ is negative and large, i.e. when the external potential $u(\psi)=-\int^{\psi}d\psi f(\psi)$ has high positive curvature. In the present case $u(\psi)$ has curvature $u^{\prime\prime}(\psi)=\sigma\,k\,\cos(\psi)$, which is negative for $\pi/2<\psi<3\pi/2$, so we expect that the UCNA is accurate only for small values of $\tau$ (as confirmed by numerical simulations in the next section). Having neglected the inertial term $\ddot{\psi}$ in Eq. (6), and going back from the time variable $s$ to $t$, we get the final UCNA stochastic evolution equation

[TABLE]

Crucially, from Eq. (7) we can write the Fokker-Planck equation Risken (1996) determining the evolution of the phase probability distribution $n=n(t,\psi,\omega)$

[TABLE]

where $\gamma=[1-\tau\partial_{\psi}f]^{-1}$. The stationary probability distribution $n(\psi,\omega)$ is found from Eq. (8) by setting the probability current to a constant:

[TABLE]

The solution of Eq. (9) can be obtained by imposing periodic the boundary condition $n(0,\omega)=n(2\pi,\omega)$. Following the method of Ref. Marconi et al. (2017) we thus have

[TABLE]

and the constant $\mathcal{N}$ is set by the normalization condition $\int d\psi\,n(\psi,\omega)=1$. The self-consistent equation for the order parameter $\sigma$ is obtained by integrating over the “disorder” $\omega$, i.e. by integrating over all oscillators each characterized by a natural frequency $\omega$

[TABLE]

To find the critical coupling strength $k_{c}$ and the behavior of $\sigma$ near $k_{c}$ we expand the r.h.s. of Eq. (11) up to $\sigma^{3}$ so that this equation becomes

[TABLE]

where

[TABLE]

Obviously, $\sigma=0$ is always a solution of Eq. (12). However, an additional real positive root of the form

[TABLE]

appears when ${k\int d\omega\,Q(\omega)>1}$ yielding the equation ${k_{c}=1/\int d\omega\,Q(\omega)}$, i.e.

[TABLE]

Note that Eq. (14) suggest that $k_{c}$ is a decreasing function of $\tau$ and reduces to the well known formula for the critical coupling of the KM in presence of white noise in the limit $\tau\rightarrow 0$ (see Ref.s Sakaguchi (1988) and Bag et al. (2007)).

It is also interesting to note that the equation for the order parameter derived within the UCNA (Eq. (13)) depends on $\tau$ only via the integrals of $P$ and $Q$ while keeping the $k$-dependence untouched. In practice $\sigma(k)$ at $\tau>0$ is a shifted and scaled version of the standard $\sigma(k)$ curve of the white-noise driven KM as sketched in Fig. 1.

This implies that $\sigma$ is a “universal” function of $k$ when properly scaled. To see this let us introduce the reduced coupling constant $\kappa=[k-k_{c}(\tau)]/k_{c}(\tau)$ so that Eq. (13) becomes

[TABLE]

where ${\mathcal{A}(\tau)=(\int d\omega\,Q)^{3/2}/(\int d\omega\,P)^{1/2}}$. It is evident from Eq. (15) that the scaled order parameter $\sigma/\mathcal{A}$ as a function of $\kappa$ is independent on $\tau$. We stress that the universal functional form of Eq. (15) is expected to be valid for small $\kappa$, however the validity of Eq. (15) is not restricted only to “critical” regimes where $k$ is extremely close to $k_{c}$ but it should hold also for larger values of $\kappa$ (see the inset of Fig. 1). In the near-critical regime, $\sigma$ can be further approximated by

[TABLE]

Summarizing our results, by using the UCNA, we have found that, while the critical coupling strength decreases with $\tau$ (Eq. (14)), the behavior of the order parameter near $k_{c}$ should be qualitatively unaffected by the memory of the noise (Eq. (15)) and also the mean-field critical exponent $1/2$ should remain the same as the one of the white-noise driven KM (Eq. (16)).

III Simulations

In order to test the results of the previous section, we integrate numerically the equations of motion (1) and (2) by using an Euler algorithm that, taking advantage of the peculiar form of the interaction, results in a computational cost $\mathcal{O}(N)$ Januszewski and Kostur (2010). Following previous works Bag et al. (2007); Tönjes (2010) we simulate $N=5000$ oscillators fixing $D=1$ scanning several values of $k$ and $\tau$. We make the usual choice of distributing the $\omega_{i}$ according to a Cauchy-Lorentz probability function

[TABLE]

setting $\lambda=1/2$. In all simulations the time-step is set to $10^{-3}$, each simulation runs at fixed $k$,$\tau$ and $D$ for $2\times 10^{6}$ steps and the order parameter $\sigma$ is computed as in Eq. (3) over the last $10^{6}$ steps averaging over all spins.

Simulation results are shown in Fig. 2 where it is clear that the order parameter goes to zero at $k$ values that become progressively lower as $\tau$ increases. In Fig. 2 we also compared the simulation data with the numerical solution of the theoretical Eq. (11) for $\tau=0$. To obtain this solution we evaluate the integrals in Eq. (11) by the trapezoidal rule and solve at each value of $k$ by using the Newton method. The close agreement between the simulations and Eq. (11) confirms that a system with $N=5000$ behaves essentially as in the thermodynamic limit (as noticed e.g. in Ref. Bag et al. (2007)) with some deviation observed for $k<k_{c}$.

From Fig. 2 follows that $k_{c}$ is a decreasing function of $\tau$ as predicted by the theory (see Eq. (14)). However, a more accurate comparison shows clear quantitative deviations. To see this, we calculate $k_{c}(\tau)$ from Eq. (14) with the choice (17)

[TABLE]

and compare this formula with the values of $k_{c}$ obtained numerically in Fig. 3. We see that Eq. (18) is close to the numerical $k_{c}$ for small $\tau$ yielding the exact result $k_{c}(0)=2(D+\lambda)$ at $\tau=0$, however it quickly deviates from the numerical data as $\tau$ increases.

In particular, while Eq. (18) predicts that $\lim_{\tau\rightarrow\infty}k_{c}(\tau)=0$, the numerical value of $k_{c}$ approaches an asymptotic value $k_{c}(\infty)>0$. This behavior is expected if we consider that the colored noise “intensity” $\langle\eta^{2}\rangle=D/\tau$ (from Eq. (2)) goes to zero in the $\tau\rightarrow\infty$ limit at fixed $D$. Therefore we expect that, in this limit, the critical coupling coincides with the $k_{c}=2\lambda$ of the KM in absence of noise Bag et al. (2007). Moreover, the breakdown of the UCNA for large $\tau$ values is expected from the discussion of the previous section: when the inertial term in the Langevin equation (6) cannot be neglected the mapping of the original colored noise problem into a white noise-driven dynamical system fails at large $\tau$.

To underline even more how the UCNA quantitatively deviates from the simulations we fit the $k_{c}$ data with a simple empirical formula (see Fig. 3). Assuming that the actual $k_{c}(\tau)$ has the form of a hyperbola and requiring that this coincides with the known limiting values

[TABLE]

we fit the simulation data points with

[TABLE]

where $\alpha$ is an adjustable parameter. In Fig. 3 we show that such a simple equation interpolates well the data (with $\alpha\approx 0.5$) signaling that, while the UCNA leads to a quantitatively wrong result, the functional form of $k_{c}(\tau)$ is still well captured by a hyperbola.

Finally, we discuss the scaling of the order parameter $\sigma$ near the transition. As shown in Eq. (15) when $k-k_{c}(\tau)$ is small the critical behavior of $\sigma$ should be invariant with respect to the correlation time of the noise $\tau$. Therefore all data-points, collected varying $\tau$, should collapse in the same power-law when properly scaled. To test this in Fig. 4 we plot the scaled $\sigma$ as a function of the scaled coupling $[k-k_{c}(\tau)]/k_{c}(\tau)$ (using the $k_{c}$ obtained numerically). To scale the $y$-value in this plot we minimize numerically the difference between the simulation values of $\sigma$ and the theoretical curved $\sigma(k-k_{c})$ obtained by solving numerically Eq. (11) for $\tau=0$ (also shown in Fig. 2) in the near-critical region $(k-k_{c})/k_{c}<0.5$. Quite surprisingly we find that a very good data collapse is observed even far from the critical region (i.e. up to values of $(k-k_{c})\approx 2k_{c}$). To visualize this we also plot Eq. (15) in Fig. 4 and we see that a good data collapse is found even where Eq. (15) deviates substantially from the data points. Justifying this observation from a theoretical point of view is, however, very challenging (even within the simplification of the UCNA) since one in principle should solve analytically the transcendental self-consistency equation defining the order parameter (Eq. (11))

In conclusions simulations show that the prediction of the UCNA is qualitatively correct but quantitatively inaccurate for locating the critical coupling strength in presence of colored noise. The scaling of the order parameter predicted by the theory is confirmed numerically and it is found to extend well beyond the expected regime of validity suggesting that the functional form of the order parameter $\sigma(k)$ is nearly invariant with respect to the correlation time of the noise for all values of $\tau$.

IV Conclusions

In this paper, we have studied the KM driven by exponentially correlated noise. We have shown that the UCNA allows to obtain analytical predictions that are in qualitatively agreement with numerical simulations. In particular, we computed the critical coupling $k_{c}(\tau)$ that results to be a decreasing function of the correlation time of the noise $\tau$. This finding is confirmed by numerical data that clearly show the increasing of magnitude of the order parameter $\sigma$ as $\tau$ increases. The observed phenomenology is in agreement with early works on similar models Bag et al. (2007); Tönjes (2010); Moro and Sánchez (1998). Moreover, the theory predicts that colored noise does not change the critical behavior of the order parameter $\sigma$, as confirmed by numerical simulations that provide evidence for an even more general universal scaling of $\sigma$ as a function scaled coupling $(k-k_{c})/k_{c}$.

What we learn in general from this work is that a persistent noise helps the system to reach a more synchronized state (without affecting the universality class of the model). More specifically we find that the UCNA allows to derive analytically this general result. A similar scenario has been also observed in a closely related lattice-model driven out-of-equilibrium by correlated noise, i.e. the $xy$-model in the presence of a Gaussian colored noise Paoluzzi et al. (2018). The picture is also similar to the one of a colored-noise driven $\phi^{4}$ theory where, however, it was found that the critical coupling does not change monotonically with $\tau$ Paoluzzi et al. (2016). Finally, we want to stress that the fact that these models remain qualitatively identical to the corresponding white-noise driven systems is not a general feature of all colored-noise driven systems. This picture differs substantially, for example, from the one found in active repulsive particle systems where the color of noise completely changes the nature of the effective interactions Farage et al. (2015); Marconi et al. (2016) and can induce phase transitions that are not present at $\tau=0$ Cates and Tailleur (2015). In this context, the UCNA appears as an invaluable tool to understand qualitatively the impact of the noise memory on the phase behavior.

Acknowledgments

MP acknowledges funding from Regione Lazio, Grant Prot. n. 85-2017-15257 (”Progetti di Gruppi di Ricerca - Legge 13/2008 - art. 4”). This work was also supported by the Joint Laboratory on “Advanced and Innovative Materials”, ADINMAT, WIS-Sapienza (MP).