The mathematics of asymptotic stability in the Kuramoto model

1. Introduction

(a) The Kuramoto model of coupled oscillators

The Kuramoto model is the archetype of collective systems composed of heterogeneous individuals that are influenced by attractive pairwise interactions. Originally designed to mimic chemical instabilities [1,2], it has since become a standard of the transition to synchrony in agent-based systems, and has been applied to various examples in disciplines such as Condensed Matter, Neuroscience and Humanities [3,4].

In its simplest form, this model considers a collection of $N\in \mathbb{N}$ oscillators, represented by their phase θ_i, a variable in the unit circle ${\mathbb{T}}^1=\mathbb{R}/2\pi \mathbb{Z}$. The population dynamics is governed by the following set of globally coupled first-order ODEs

The time-independent frequencies ${\omega }_i\in \mathbb{R}$ are randomly drawn in order to account for individual heterogeneities. The parameter $K\in {\mathbb{R}}^{+}$ measures the interaction strength. While, up to time rescaling, K could be absorbed in the frequencies ω_i↦ω_i/K, it is more convenient to investigate the dependence of the dynamics upon this parameter, for a given frequency distribution.

Clever intuition, elaborate analytic considerations and extensive numerics have provided comprehensive insights into the Kuramoto phenomenology (e.g. [5–12]). Nonetheless, due to heterogeneities, statements about the full nonlinear dynamics, certified by complete mathematical proofs, are rather scarce. They can be summarized as follows.

For weak interactions, KAM theory for dissipative systems [13, thm 6.1] or [14, thm 3.1] asserts that, for frequencies in a Lebesgue positive set in ${\mathbb{R}}^N$, the dynamics for K small is conjugated to the system at K = 0.¹ This conjugacy implies infinite returns to arbitrary small neighbourhoods of the initial condition in ${\mathbb{T}}^N$.²

Results for strong interactions contrast with weak coupling recurrence, see [17–20] and also [18,21,22] for additional interesting statements. For $K>\underset{i,j}{max}|{\omega }_i-{\omega }_j|$ and provided that initial phase spreading is limited enough, full locking asymptotically takes place; i.e. the limit

exists for every pair (i, j). In the extreme case of homogeneous populations (i.e. ω_i independent of i), complete synchrony holds for every K > 0, viz. we have for almost every initial phases. However, for exceptional initial conditions, the population can cluster into two synchronized groups [17,23]. This clustering is not limited to homogeneous populations and may hold for symmetric frequency distributions [24].

No rigorous results exist about (1.1) in regimes when interactions and heterogeneities effects balance. However, insights can be obtained by considering the continuum limit approximation.

(b) The Kuramoto PDE: basic features

(i) Kuramoto dynamics at the continuum limit

The continuum approximation assumes that populations at the thermodynamic limit N →  + ∞ are described by absolutely continuous distributions f on the cylinder ${\mathbb{T}}^1\times \mathbb{R}$, more precisely, by their densities. Under this assumption, time evolution is governed by the following PDE [25,26]:

where For any trajectory of (1.1), the empirical measure ${\mu }_N(t)=(\frac{1}{N})\sum _{i=1}^N{\delta }_{{\theta }_i(t),{\omega }_i}$ (here δ_· ,· stands for the Dirac distribution) is a weak solution of (1.2). The continuum approximation is justified by the following basic features [27], which are common to mean-field models in classical mechanics, especially the Vlasov equation [28–30].

—

The Cauchy problem is globally well posed for (1.2), viz. for every initial probability measure f(0) on the cylinder, there exists a unique solution t↦f(t) defined for all t≥0. If f(0) is absolutely continuous, then so is f(t) for every t > 0.

—

The solution continuously depends on the initial condition, in the weak topology. More precisely, if dBL( · , · ) denotes the bounded Lipschitz distance of measures, then there exists C > 0 such that for every pair of solution t↦fi(t), i = 1, 2, we have $$d_{\text{B}L}(\hspace{0.17em}{f}_1(t),\hspace{0.17em}{f}_2(t))\le d_{\text{B}L}(\hspace{0.17em}{f}_1(0),\hspace{0.17em}{f}_2(0))\hspace{0.17em}{\mathrm{e}}^{Ct},\mathrm{\forall }t>0.$$

For N sufficiently large, μ_N(0) can be chosen close to an absolutely continuous distribution f(0), i.e. such that d_BL(μ_N(0), f(0)) is small. The inequality above then implies the continuum approximation on finite time interval, i.e. d_BL(μ_N(t), f(t)) remains small for t small enough.

In addition, the Kuramoto PDE has the following specific features.

—	Galilean invariance: if t↦f(t) is a solution, then t↦RΘ+Ωt,Ωf(t) is a solution for every $(\mathit{\Theta },\mathit{\Omega })\in {\mathbb{T}}^1\times \mathbb{R}$, where RΘ,Ω is the representation on measures of the map (θ, ω)↦(θ + Θ, ω + Ω). In particular, (1.2) is equivariant with respect to the rigid rotation RΘ: = RΘ,0.
—	For every solution t↦f(t), the frequency marginal ${\int }_{{\mathbb{T}}^1}f(t,\mathrm{d}\theta ,\hspace{0.17em}\mathrm{d}\omega )$ does not depend on t, and thus can be regarded as an input parameter on initial conditions.

(ii) Basic phenomenology

The degree of synchrony in Kuramoto dynamics can be characterized by the order parameter

In particular, stationary states can be classified accordingly: the case r = 0 corresponds to the homogeneous stationary state $f_{\mathrm{hom}}(\mathrm{d}\theta ,\hspace{0.17em}\mathrm{d}\omega )=\left(\frac{g(\omega )}{2\pi }\right)\hspace{0.17em}\mathrm{d}\theta \hspace{0.17em}\mathrm{d}\omega$ for which θ is uniformly distributed independently of ω, while r≠ 0 deals with partially locked states (PLS), for which angles such that |θ| ≤ |r|ω are uniquely attached to a single value of ω. (An expression of PLS is given below.)

The Kuramoto PDE displays a large phenomenology depending on the interaction strength and the (absolutely continuous) frequency marginal. The simplest case is when the corresponding density g is unimodal and symmetric. (Thanks to Galilean invariance, its maximum can be set at the origin 0.) Then, the phenomenology can be summarized as follows (figure 1a) [31]:

—	for K < Kc: = 2/πg(0), fhom is asymptotically stable. Hence, for every trajectory the order parameter asymptotically vanishes, limt → +∞r(t) = 0. This convergence is the analogue of the Landau damping phenomenon in the Vlasov equation [32]. (It is incompatible with KAM induced recurrence behaviours in (1.1). Hence, the continuum approximation mentioned above cannot hold for all t > 0.)
—	at K = Kc, fhom becomes unstable and a circle of stable stationary PLS emerges for K > Kc (together with a continuum of unstable PLS). In this regime, we have limt → +∞|r(t)| = |rpls|≠ 0 provided that f(0) is not in the stable manifold of fhom.

• Download figure
• Open in new tab
• Download PowerPoint

Illustrations of the dynamics in the finite-dimensional model, both for K < K_c and K > K_c, are provided in the electronic supplementary material, Movies. More elaborate bifurcation schemes occur for other marginals [33,34], especially for the bi-Cauchy distribution, see figure 1b and details in §5, and also for extensions of the model [35–39]. Even for asymmetric unimodal marginals, the phenomenology can be involved.

(iii) Proving the phenomenology: state-of-art and technical considerations

While the phenomenology above had been identified in early studies, full rigorous confirmation has remained elusive until recently. Mathematical studies have long been limited to linearized dynamics in strong topology. They have provided both stability criteria [26] and evidence that the relaxation rate, either algebraic or exponential, depends on g's regularity [11]. One should also mention that, for special frequency marginals, a rather impenetrable proof of f_hom stability is exposed in [40]. In addition, complete synchronization has been proved to hold when g is the Dirac distribution [41]. Besides, solid arguments have been provided for convergence of the order parameter dynamics to the corresponding one in the so-called Ott–Antonsen (OA) manifold [42,43]. There, the dynamics is governed by a finite-dimensional system when g is meromorphic with finitely many poles in the lower half-plane [9]; hence a standard analysis of stability and bifurcations can be developed in this case [34] (see also §5.6.2 in [44]).

The major obstacle to including nonlinearities in proofs is that, due to the free transport term ω∂_θf in (1.2), in strong topology, the linearized dynamics has continuous spectrum on the imaginary axis [45]. In fact, stationary states have all been shown to be nonlinearly unstable in the L²-norm [46,47]. Therefore, any proof of asymptotic stability must consider weaker topology.

For suitable norms in weak topology and analytic frequency marginals, the linearized dynamics essential spectrum is located to the left of the imaginary axis in the complex plane [48]. Provided that the remaining discrete spectrum is under control—hence the stability conditions—a standard strategy for asymptotic stability can be considered: since the linearized dynamics decays exponentially fast, it can dominate nonlinear instabilities for small enough perturbations. In practice, the proof is not so straightforward and needs adjustments, especially because angular derivatives in (1.2) imply that nonlinearities can be large even for small perturbations.

When the frequency marginal has only algebraic regularity, this strategy no longer applies because no spectral gap is at hand. Instead, the specific structure of linearized perturbation dynamics, which takes the form of a Volterra equation, needs to be exploited in order to prove algebraic damping via advanced bootstrap arguments [49].

Besides, PLS stability deals with circles of stationary states, by equivariance with respect to rotations R_Θ. The perturbation dynamics then must be neutral with respect to tangential perturbations. Asymptotic stability is to be proved for the relative equilibrium of the dynamics in the radial variable [50].

(c) Organization of the rest of the paper

This paper aims to review stability results and their proofs for stationary solutions of (1.2) that have been obtained in [44,46,48,49,51]. In few words, these results claim asymptotic convergence in the weak sense to either f_hom or to some PLS, depending on a corresponding stability condition. In addition, control of the order parameter relaxation speed will be given, which depends on the regularity of the initial condition (including the frequency distribution). Of note, thanks to Galilean invariance, all results immediately extend to globally rotating solutions.

The results are presented in §2. Section 3 provides insights into the main arguments of proofs, especially those that are likely to be of interest to readers not familiar with the analysis of PDEs. This includes linear stability analysis via considerations on Volterra equations and control of nonlinear terms by means of a Gearhart–Prüss-like argument. The key point is to obtain weighted L² estimates on the solution's Fourier transform. This is not only critical for the proofs, but it also implies both convergence to the centre manifold and to the OA manifold mentioned above (§4).

Stability conditions in the statements will be expressed in terms of K and g. Consistency considerations on these conditions are evaluated in §5, where bifurcation diagrams are also provided for various examples of frequency distributions, including original ones.

Finally, §6 mentions some extensions of the Kuramoto model for which the approaches presented here apply to yield rigorous results on asymptotic stability. Limitations and open questions are also briefly discussed.

2. Main results: asymptotic stability in the Kuramoto PDE

This section describes the behaviour of solutions t↦f(t) of (1.2), for a given absolutely continuous frequency marginal ${\int }_{{\mathbb{T}}^1}f(t,\mathrm{d}\theta ,\mathrm{d}\omega )=g(\omega )\hspace{0.17em}\mathrm{d}\omega$. More precisely, existence and stability conditions are given for the stationary states, together with the corresponding local basins of attraction.

The order parameter relaxation speed depends on the regularity of the frequency marginal and of the initial perturbation: more regularity implies faster decay. Regularity is usually quantified by Fourier transform decay. Given a function u on $\mathbb{R}$ and a measure v on ${\mathbb{T}}^1\times \mathbb{R}$, their Fourier transforms are defined by

Various constraints on Fourier transforms have been suggested, see end of §2a below. For simplicity, we shall express constraints in terms of weighted norms. Given a weight function $\varphi :{\mathbb{R}}^{+}\to {\mathbb{R}}^{+}$ and a sequence of functions $u=\{u_{\ell }(\tau )\}\in {\mathbb{C}}^{\mathbb{N}\times {\mathbb{R}}^{+}}$, let and let also Typical weights are ϕ(τ) = e^aτ (a > 0) and ϕ(τ) = (1 + τ)^b (b > 1). In ${\mathbb{C}}^{\mathbb{N}\times \mathbb{R}}$, the norm $\parallel \cdot {\parallel }_{{\mathcal{H}}_{\varphi }^1(\mathbb{N}\times \mathbb{R})}$ and space ${\mathcal{H}}_{\varphi }^1(\mathbb{N}\times \mathbb{R})$ are defined similarly. These norms are designed to accommodate the Fourier transforms of the appearing singular measures, such as PLS. Moreover, we shall need the following norms on frequency marginals Of note, together with $\hat{g}(-\tau )=\overline{\hat{g}(\tau )}$, the condition $\parallel \hat{g}{\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1({\mathbb{R}}^{+})}<+\mathrm{\infty }$ implies, via the Paley–Wiener theorem, that g must be analytic in a strip around the horizontal axis in $\mathbb{C}$.

While focus is on asymptotic stability of certain solutions, the developed exponential stability analysis also fits the setting of the theory of centre manifolds in infinite dimension [52]. In particular, for Banach spaces defined using weighted L^∞-norms for Fourier transforms, a centre-unstable manifold has been proved to exist for K∼K_c, which attracts all trajectories of (1.2) in a sufficiently small neighbourhood of f_hom (see theorem 7 in [46] and also [40]).

(a) Asymptotic relaxation to the homogeneous state

A unique homogeneous stationary state f_hom(dθ, dω) = (g(ω)/2π) dθ dω exists for every g and K, and its order parameter vanishes r_hom = 0. In addition to regularity requirements on g, the stability of this state relies on the following condition [46,51]:

which involves the Laplace transform of $\hat{g}$. When g is symmetric and unimodal, this requirement is equivalent to the inequality K < 2/πg(0) mentioned in the Introduction. In the general case, the argument principle and the Plemelj formula can be employed to show that (2.1) holds under the following analogue of the Penrose criterion in the Vlasov literature: This criterion is however not necessary for stability; the tri-Cauchy distribution at the end of §5 provides a counter-example.

Theorem 2.1.

Assume that$g\in C^2(\mathbb{R})$is such that$\parallel \hat{g}{\parallel }_{L_{{(1+\tau )}^b}^1({\mathbb{R}}^{+})}<+\mathrm{\infty }$for someb > 1 and letKbe so that (2.1) holds. Then, there existsϵ > 0 such that for every$f(0)\in C^2({\mathbb{T}}^1\times \mathbb{R})$with marginal densitygand satisfying$\parallel \widehat{f(0)}{\parallel }_{{\mathcal{H}}_{{(1+\tau )}^b}^1(\mathbb{N}\times {\mathbb{R}}^{+})}<ϵ$, we have, in the weak sense,

$$\underset{t\to +\mathrm{\infty }}{lim}f(t)=f_{\mathrm{hom}}.$$

This statement, which by definition of weak convergence, implies lim_{t → +∞}r(t) = 0 for the order parameter associated with f(t), is an immediate consequence of the combination of theorems 5 and 39 in [46].

The stability condition (2.1) is optimal, as least as far as linear stability is concerned. This means that, if there exists $z_0\in \mathbb{C}$ with Re(z₀) > 0 such that $(K/2){\int }_{{\mathbb{R}}^{+}}\hat{g}(\tau )\hspace{0.17em}{\mathrm{e}}^{-z_0\tau }\hspace{0.17em}\mathrm{d}\tau =1$, then the linearized Kuramoto equation around f_hom has a solution with exponentially growing order parameter [26].

The constraint $\parallel \widehat{f(0)}{\parallel }_{{\mathcal{H}}_{{(1+\tau )}^b}^1(\mathbb{N}\times {\mathbb{R}}^{+})}<ϵ$ actually impacts the perturbation f(0) − f_hom and is justified by the possible existence of PLS while (2.1) holds (e.g. [34,39,48]). However, such coexistence can only happen for relatively strong interaction and the conclusion of theorem 2.1 can be asserted for any C⁴ initial perturbation of finite weighted Sobolev norm ${\mathcal{H}}^4$, provided that K is small enough (proposition 3.2 in [51] combined with theorems 39 in [46]). Moreover, when focus is made on the observable r(t), the uniform constraint

ensures that lim_{t → +∞}r(t) = 0 holds for the trajectory of every initial perturbation in ${\mathcal{H}}_{{(1+\tau )}^b}^1(\mathbb{N}\times {\mathbb{R}}^{+})$ (see section 4.4 in [44]).

Asymptotic convergence of f(t) in theorem 2.1 is proved using accurate control of the relaxation rate of its Fourier transform. As anticipated in [11], the order parameter relaxation rate can be estimated based on the initial perturbation regularity.

Proposition 2.2.

(i)	Under the conditions of theorem 2.1, we have $$r(t)=O(t^{-b}).$$
(ii)	If, in addition, $\parallel \hat{g}{\parallel }_{L_{{\mathrm{e}}^{a\tau }}^1({\mathbb{R}}^{+})}<+\mathrm{\infty }$for some a > 0, then, there existϵ, a′ > 0 such that, for every$f(0)\in C^2({\mathbb{T}}^1\times \mathbb{R})$with$\parallel \widehat{f(0)}{\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1(\mathbb{N}\times \mathbb{R})}<ϵ$, we have $$r(t)=O({\mathrm{e}}^{-a^{\prime }t}).$$

Statement (i) (resp. (ii)) is a consequence of theorem 5 (resp. 4) in [46]. Under the conditions of (ii), the essential spectrum of the linearized dynamics operator is contained in the half-plane Re(z) ≤  − a. In the complement half-space, the spectrum consists of finitely many eigenvalues, all of them have negative real part. The rate a′ corresponds to these eigenvalue largest real part.

Finally, one can mention that the conclusions of (i) and (ii) apply under weaker assumptions on g and f(0) [46]. These conditions express as pointwise constraints on $\hat{g}$ and $\widehat{f(0)}$ and imply similar pointwise decay estimates for $\widehat{f(t)}$. Besides, inspired by a similar analysis for the Vlasov-HMF equation [53], Fernandez et al. [51] considers perturbations in the original space of measure densities via the weighted Sobolev norm defined by

In this setting, polynomial Landau damping is obtained as in (i) above, as well as algebraic convergence in a twisted reference frame, of the Sobolev norm of the density of f(t).

(b) Asymptotic relaxation to PLS

Recall that R_Θ denotes the representation on measures of the rigid rotation on the cylinder. In full generality, PLS can be defined as solutions of (1.2) of the form t↦R_Ωtf_pls, for some global frequency $\mathit{\Omega }\in \mathbb{R}$ and reference measure f_pls with non-vanishing order parameter

Galilean invariance implies that PLS are indifferent to the action of R_θ; hence we may assume that the state f_pls has real and positive order parameter r_pls and consider the circle ${\{R_{\mathit{\Theta }}{f}_{\mathrm{pls}}\}}_{\mathit{\Theta }\in {\mathbb{T}}^1}$. For the same reason, we may assume that Ω = 0, up to a translation of g. In this case, every R_Θf_pls is a stationary state and one can show that the corresponding expression of f_pls is [3,31,45] which confirms the singular nature of PLS [31]. Here the measurable function α:[ − Kr_pls, Kr_pls] → [0, 1] quantifies the relative contribution of the two equilibria $\arcsin (\omega /K{r}_{\mathrm{pls}})$ and $\pi -\arcsin (\omega /K{r}_{\mathrm{pls}})$ of the equation of characteristics The PLS with α = 1 a.e. is denoted by f_s, and its order parameter by r_s. The equilibrium $\arcsin (\omega /K{r}_{\mathrm{pls}})$ is stable for the one-dimensional dynamics (2.2), while the other one is unstable. This suggests that only f_s can be stable among possible f_pls [45]. This argument can be formally justified using the norm above. Indeed, provided that $\parallel \hat{g}{\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1(\mathbb{R})}<+\mathrm{\infty }$, the only f_pls whose Fourier transform lies in ${\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1(\mathbb{N}\times \mathbb{R})$ turns out to be f_s (proposition A.2 in [48]).

The explicit expression of f_pls yields an existence condition, which materializes as a self-consistency condition on r_pls [31,45,48]. For f_s, this condition writes

Of note, if g is symmetric around 0, then the imaginary part of the l.h.s. here automatically vanishes and the existence condition becomes [31,38,39,45] If also g∈C⁰, an analysis of this condition shows that a PLS f_s exists for every K > K_c [31]. Moreover, if g is unimodal, f_s is unique and does not exist for K ≤ K_c. For more general frequency marginals, PLS existence is not so simple but can be granted, once allowing for Ω≠0, under the condition that the homogeneous state is unstable, see §5.

For the stability condition, the following notations are needed: given $z\in \mathbb{C}$ with Re(z)≥ 0 and $r\in {\mathbb{R}}^{+}$, let M(z, r) be the 2 × 2 matrix defined by

(where the quantities J_k are defined by continuity for Re(z) = 0).

Theorem 2.3.

Given$b>\frac{3}{2}$, assume that$\parallel \hat{g}{\parallel }_{{\mathcal{H}}_{{(1+\tau )}^{b_g}}^1(\mathbb{R})}<+\mathrm{\infty }$for someb_g > b + 3 and letKbe such that a stationary PLSf_swith marginal densitygand order parameter$r_{\text{s}}\in {\mathbb{R}}^{+}$exists and satisfies

$$\begin{array}{rl}& det(\mathrm{Id}-\frac{K}{2}M(z,r_{\text{s}}))\ne 0,\hspace{0.17em}\mathrm{\forall }z\ne 0\hspace{0.17em}\mathrm{with}\hspace{0.17em}\mathrm{Re}(z)\ge 0,\\ \mathrm{and}& z=0\hspace{0.17em}\text{is a simple zero of the function}\hspace{0.17em}z\mapsto det(\mathrm{Id}-\frac{K}{2}M(z,r_{\text{s}})).\end{array}\}$$2.4

Then, there existsϵ > 0 such that for everyf(0) with marginal densitygand so that

$$\parallel \widehat{f(0)}-\widehat{R_{\mathit{\Theta }}{f}_{\text{s}}}{\parallel }_{{\mathcal{H}}_{{(1+\tau )}^b}^1(\mathbb{N}\times {\mathbb{R}}^{+})}<ϵ\hspace{0.17em}\text{for some}\hspace{0.17em}\mathit{\Theta }\in {\mathbb{T}}^1$$

there exists${\mathit{\Theta }}_{\mathrm{\infty }}\in {\mathbb{T}}^1$such that we have

$$\underset{t\to +\mathrm{\infty }}{lim}\hspace{0.17em}f(t)=R_{{\mathit{\Theta }}_{\mathrm{\infty }}}\hspace{0.17em}{f}_{\mathrm{s}}\hspace{0.17em}\text{(weak sense)}\mathrm{and}|r(t)-r_{\mathrm{s}}\hspace{0.17em}{\mathrm{e}}^{\mathrm{i}{\mathit{\Theta }}_{\mathrm{\infty }}}|=O(t^{1/2-b}).$$

As for (2.1), the stability condition (2.4) can be shown to be optimal. Moreover, the statement above is a simplification of theorem 2 in [49], which includes broader regularity conditions and provides quantitative control of the convergence in Fourier space. As for f_hom, an analogous statement holds in the exponential setting.

Theorem 2.4.

Assume that$\parallel \hat{g}{\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1({\mathbb{R}}^{+})}<+\mathrm{\infty }$for somea > 0 and letKbe such that a stationary PLSf_swith marginal densitygand order parameter$r_{\mathrm{s}}\in {\mathbb{R}}^{+}$exists and satisfies (2.4). Then, there existϵ, a′ > 0 such that for everyf(0) with marginal densitygso that

$$\parallel \widehat{f(0)}-\widehat{R_{\mathit{\Theta }}{f}_{\text{s}}}{\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1(\mathbb{N}\times \mathbb{R})}<ϵ\hspace{0.17em}\text{for some}\hspace{0.17em}\mathit{\Theta }\in {\mathbb{T}}^1,$$

there exists${\mathit{\Theta }}_{\mathrm{\infty }}\in {\mathbb{T}}^1$so that we have (in addition to weak convergence of measures)

$$|r(t)-r_{\mathrm{s}}\hspace{0.17em}{\mathrm{e}}^{\mathrm{i}{\mathit{\Theta }}_{\mathrm{\infty }}}|=O({\mathrm{e}}^{-a^{\prime }t}).$$

This statement is a consequence of theorem 2.1 in [48], which claims the following convergence of Fourier transforms (and hence the conclusion on the order parameter)

Detailed considerations on existence and stability of PLS will be given in §5. When g is unimodal and symmetric, (2.4) turns out to coincide with the existence condition K > K_c [48]. In particular, theorems 2.3 and 2.4 complete the proof of the bifurcation diagram in figure 1a.

In addition, global stability can never hold for PLS because f_hom is a distinct stationary state, which satisfies $\parallel \widehat{f_{\mathrm{hom}}}{\parallel }_{{\mathcal{H}}_{{(1+\tau )}^b}^1(\mathbb{N}\times {\mathbb{R}}^{+})}<+\mathrm{\infty }$ under the conditions of theorem 2.3 (or 2.4). Notice finally that for r_s → 0, not only the PLS expression reduces to that of f_hom, but PLS existence and stability conditions converge as well. We kept the exposition of stationary states separated for historical and pedagogical reasons.

3. Main ingredients of proofs

The asymptotic stability of stationary states has been proved using the formulation of the dynamics in Fourier space. Instead of providing all details, we focus here on the decay of the order parameter under the linearized dynamics, which turns out to be governed by a Volterra equation of the second kind. Conditions (2.1) and (2.4), and asymptotic decay as given in proposition 2.2 and in theorems 2.3 and 2.4, then follow from the corresponding theory [54]. In addition, we will comment on how to deal with the nonlinear terms in the exponential case.

(a) Volterra equation for the order parameter

Let $u={\{u_{\ell }\}}_{\mathbb{N}}={\{u_{\ell }(\tau )\}}_{\mathbb{N}\times \mathbb{R}}$ be an initial perturbation with u₀(τ) = 0 (so that the frequency marginal is preserved). Inserting the expression ${\hat{f}}_{\mathrm{s}}+u$ in the Kuramoto dynamics in Fourier space yields the following evolutionary equation:

where and the operator Q collects the nonlinear terms Naturally, for r_s = 0 and ${\hat{f}}_{\mathrm{s}}=\widehat{f_{\mathrm{hom}}}$, these equations describe the perturbation dynamics around f_hom. In this case, we have the simplification (L₁u)_ℓ = ℓ∂_τu_ℓ and ${(L_{2}u)}_{\ell }=(K/2)u_1(0)\hat{g}{\delta }_{\ell ,1}$, while Q remains unchanged.

Prior to any other consideration, this perturbation dynamics needs to be granted well-posed in the weighted norm setting. In this respect, proposition 3.1 in [48] claims that the subset of measures with Fourier transforms in ${\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1(\mathbb{N}\times \mathbb{R})$ has well-defined Kuramoto dynamics and is invariant under the flow (see [49] for a similar well-posedness result in the algebraic setting).

Moreover, one needs to incorporate the fact that L₂ is only $\mathbb{R}$-linear and not $\mathbb{C}$-linear if ${\hat{f}}_{\mathrm{s}}\ne \widehat{f_{\mathrm{hom}}}$. One way to proceed is to treat the real and imaginary components separately [39,45,49]. Here, we adopt a different but equivalent approach that substitutes complex conjugates by an independent variable. Given $u={\{u_{\ell }(\tau )\}}_{\mathbb{N}\times \mathbb{R}}$ and $v={\{v_{\ell }(\tau )\}}_{\mathbb{N}\times \mathbb{R}}$ (which is a substitute for $\overline{u}$), let

and consider the $\mathbb{C}$-linear operators ${\mathcal{L}}_i$ (i = 1, 2) defined by using the notations ${(u_{\mathrm{s},-})}_{\ell }=\ell {(\hspace{0.17em}{\hat{f}}_{\mathrm{s}})}_{\ell -1}$ and ${(u_{\mathrm{s},+})}_{\ell }=\ell {(\hspace{0.17em}{\hat{f}}_{\mathrm{s}})}_{\ell +1}$. These extended operators are defined in such a way that when $v_{\ell }={\overline{u}}_{\ell }$, we have Moreover, it can be checked that this extension does not generate unstable spurious modes [48].

Considerations on its resolvent in ${\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1(\mathbb{N}\times \mathbb{R})$ imply that L₁ generates a C⁰-semigroup in this space [48]. In case of r_s = 0, the semigroup is the free transport, namely

In ${\mathcal{H}}_{(1+{\tau }^b)}^1(\mathbb{N}\times \mathbb{R})$, the semigroup admits weak solutions and this is enough for our purpose. The same properties directly follow for the semigroup of ${\mathcal{L}}_1$ in the corresponding product spaces.

Both u_s,− and u_s,+ belong to ${\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1(\mathbb{N}\times \mathbb{R})$ [48, proposition A.2]; hence the operator ${\mathcal{L}}_2$ must be bounded. Therefore, ${\mathcal{L}}_1+{\mathcal{L}}_2$ similarly generates a C⁰-semigroup. In the product space associated with ${\mathcal{H}}_{(1+{\tau }^b)}^1(\mathbb{N}\times \mathbb{R})$, this semigroup is also well defined [49].

When regarding $t\mapsto {\mathcal{L}}_{2}u(t)$ as a forcing term in the linearized PDE

Duhamel's principle implies that the solution writes Using the linearity of e^tL₁ and the expression of ${\mathcal{L}}_2$, a self-consistent equation results for the coordinate (ℓ, τ) = (1, 0) of the solution's component u(t) = {u_ℓ(t, τ)}. This equation is the following Volterra equation: where the two-dimensional convolution by the kernel $\mathcal{K}$ is defined by and where $I(t)=({\mathrm{e}}^{t{\mathcal{L}}_1}u{(0))}_1(0)$ is regarded as an input (which contains the initial perturbation). In particular, for inputs with conjugated initial components $\{v_{\ell }(0,\tau )\}=\{\overline{u_{\ell }(0,\tau )}\}$, we have $u_1(t,0)=\left(\genfrac{}{}{0ex}{}{\overline{r(t)}}{r(t)}\right)$ and (3.2) describes the linearized evolution of the perturbation order parameter r(t). In the case r_s = 0, we have the simplification ${(u_{\mathrm{s},-})}_{\ell }=\hat{g}{\delta }_{\ell ,0}$ and (u_s,+)_ℓ = 0, and the order parameter linearized trajectory is governed by the following one-dimensional Volterra equation where * now denotes the standard convolution of complex functions.

(b) Asymptotic decay of solutions and stability conditions

Volterra equations have unique and explicit solutions provided that their kernel and forcing are locally bounded (Sect. 3, ch. 2 in [54]). In the case of (3.2), these properties are granted by the fact that e^tL₁ is itself locally bounded either in ${\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1(\mathbb{N}\times \mathbb{R})$ or in ${\mathcal{H}}_{{(1+\tau )}^b}^1(\mathbb{N}\times {\mathbb{R}}^{+})$, together with properties of the states u_s,− and u_s,+. The solution then writes

and the definition of the resolvent ${\mathcal{R}}_{\mathcal{K}}$ relies on the induction ${\mathcal{K}}^{\ast (k+1)}=\mathcal{K}\ast {\mathcal{K}}^{\ast k}$ with ${\mathcal{K}}^{\ast 1}=\mathcal{K}$.

(i) Analysis of the one-dimensional Volterra equation

Let ${\mathcal{R}}_{(K/2)\hat{g}}$ be the resolvent of the convolution by $(K/2)\hat{g}$. Equation (3.3) has the solution

whose asymptotic properties are readily accessible, using the following weighted norm The weights ϕ(t) = e^at and ϕ(t) = (1 + t)^b are sub-multiplicative functions. Together with Young's inequality, this property implies the following inequality [46] Therefore, if we can ensure that $\parallel {\mathcal{R}}_{(K/2)\hat{g}}{\parallel }_{L_{\varphi }^1({\mathbb{R}}^{+})}<+\mathrm{\infty }$, then quantified order parameter decay $\parallel r{\parallel }_{L_{\varphi }^{\mathrm{\infty }}({\mathbb{R}}^{+})}<+\mathrm{\infty }$ will follow from a similar feature $\parallel u_1(0,\cdot ){\parallel }_{L_{\varphi }^{\mathrm{\infty }}({\mathbb{R}}^{+})}<+\mathrm{\infty }$ of the initial perturbation. (By Sobolev embedding, the latter holds provided that $u(0)\in {\mathcal{H}}_{\varphi }^1(\mathbb{N}\times {\mathbb{R}}^{+})$.) In particular, the conclusions in proposition 2.2 for the linear dynamics will be instances of this property, when applied to the exponential and polynomial weights, respectively.

It remains to connect the constraint $\parallel {\mathcal{R}}_{(K/2)\hat{g}}{\parallel }_{L_{\varphi }^1({\mathbb{R}}^{+})}<+\mathrm{\infty }$ to the condition (2.1) in each case. This equivalence is given by the half-line Gelfand theorem (theorem 4.3, ch. 4 in [54]). Indeed, since ϕ is sub-multiplicative and the measure $\hat{g}(t)\hspace{0.17em}\mathrm{d}t$ is absolutely continuous, this statement implies that the desired constraint holds under the conditions $\hat{g}\in L_{\varphi }^1({\mathbb{R}}^{+})$ and

The conditions of proposition 2.2 (i) then immediately follow. In the exponential case, the constraint (3.4) appears more stringent than (2.1) (because −lim_t →+∞(lnϕ(t)/t) < 0). To see that the conditions of proposition 2.2 (ii) suffice, observe that $\parallel \hat{g}{\parallel }_{L_{{\mathrm{e}}^{a\tau }}^1({\mathbb{R}}^{+})}<+\mathrm{\infty }$ implies that the Laplace transform is holomorphic in every half-plane Re(z) > − a and continuous up to the boundary Re(z) = − a. By the Riemann–Lebesgue lemma, this function must be uniformly small outside a sufficiently large rectangular region of the form In particular, it cannot reach the value 2/K outside this domain. By analyticity, it can only reach this value at finitely many points with Re(z) > − a. Assuming (2.1), each of these points must satisfy Re(z) < 0. Therefore, all these points must satisfy Re(z) ≤  − a′ for some a′∈(0, a). It follows that (2.1) implies (3.4) for ϕ(t) = e^a′t, as desired.

(ii) Analysis of the two-dimensional Volterra equation

That PLS come in circles of stationary solutions implies that the linearized dynamics at $\widehat{f_{\mathrm{s}}}$ should be neutral with respect to perturbations that are tangent to the circle [45]. In fact, we have [48]

As a consequence, the solution of (3.2) cannot be decaying when the initial input lies along the corresponding 0-eigenmode of ${\mathcal{L}}_1+{\mathcal{L}}_2$ in the product space. Asymptotic decay can only hold for solutions whose input is initially transversal to this direction. In order to integrate this constraint, one should require that the analogue condition to (3.4) excludes the eigenvalue 0, i.e. Using that the Laplace transform of a semigroup is the resolvent of its generator and proceeding with algebraic manipulations on this resolvent [48, lemma 4.4], yield the equality and the first constraint in (2.4) follows suit. Together with imposing that 0 is a simple eigenvalue (second constraint in (2.4)) and the conditions on $\hat{g}$ in theorem 2.3 (resp. 2.4), another result in the theory of Volterra equations (theorem 3.7, ch. 7 in [54] and subsequent comment) implies that the resolvent writes where C is the constant 2 × 2 matrix corresponding to the 0-eigenmode above and where for ϕ(t) = (1 + t)^b (resp. ϕ(t) = e^at). The desired damping for transversal perturbations then results from the following inequality (obtained using a similar reasoning as above) where I(t) now denotes the first component of the input $({\mathrm{e}}^{t{\mathcal{L}}_1}u{(0))}_1(0)$ when assuming conjugated components $\{v_{\ell }(0,\tau )\}=\{\overline{u_{\ell }(0,\tau )}\}$ in the initial perturbation u(0).

Notice finally that, unlike in the previous section, estimates on $\parallel I(t){\parallel }_{L_{\varphi }^{\mathrm{\infty }}({\mathbb{R}}^{+})}$ are not immediate here, even when $\parallel u_1(0,\cdot ){\parallel }_{L_{\varphi }^{\mathrm{\infty }}({\mathbb{R}}^{+})}<+\mathrm{\infty }$. In the exponential case, these estimates follow from the semigroup exponential stability in ${\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1(\mathbb{N}\times \mathbb{R})$ [48], namely

for some a′′∈(0, a); hence the rate a′ < a′′ in theorem 2.4, when combined with (2.4) and similar analytic arguments to those in the previous section. In the algebraic case, no such property can exist for e^tL₁ in ${\mathcal{H}}_{{(1+\tau )}^b}^1(\mathbb{N}\times {\mathbb{R}}^{+})$ (otherwise, we would have exponential decay in this space). Instead, an energy estimate yields the following inequality (see lemma 4 in [49]) and algebraic decay follows from a control of the integral using the Cauchy–Schwarz inequality.

(c) Control of nonlinear terms in the exponential case

With full understanding at the linear level, nonlinearities remain to be accounted for. Here, focus is made on PLS. Similar considerations apply for f_hom. The fact that PLS come in circles requires to get rid of the angular coordinate and to consider the radial dynamics only [50]. It can be shown that the Fourier transform $\hat{f}$ of any measure close enough to ${\{{\hat{R}}_{\mathit{\Theta }}{\hat{f}}_{\mathrm{s}}\}}_{\mathit{\Theta }\in {\mathbb{T}}^1}$ can be written

and P_s is an appropriate projection on the complement of Ker(L₁ + L₂) [48]. By inserting this expression in (3.1) and by applying P_s, the following nonlinear equation results where Q′ is an updated nonlinear term, independent of the angular variable Θ. In the case of f_hom, no projection is needed, and the considerations below apply mutatis mutandis to (3.1).

The previous section showed that (2.4) implies in the exponential case that the semigroup e^t(L₁+L₂) is exponentially stable, namely

In a standard proof of sink asymptotic stability, the nonlinear terms are assumed to be sufficiently regular, say C², so that for sufficiently small perturbations, they can be dominated by the linear exponential stability and exponential decay of the full system solution follows.

Unfortunately, the nonlinearity Q (and hence Q′) is not regular at all in ${\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1(\mathbb{N}\times \mathbb{R})$; in fact it does not even map this space into itself. Instead, we have

and the new norm is an extension of the one introduced in §2Nonetheless, the linear terms (3.5) have enough regularizing effect to dominate nonlinearities. In fact, when regarding the nonlinearity in (3.5) as a forcing term, the following adaptation of the Gearhart–Prüss Theorem shows asymptotic decay. Given an Hilbert space H with norm ∥ · ∥_H, a number $\gamma \in {\mathbb{R}}^{+}$, and a mapping $w:\mathbb{R}\to H$, consider the norm defined by

Lemma 3.1 ([48]).

LetX↪Ybe Hilbert spaces andAbe a densely defined linear operator that generates aC⁰-semigroup on bothXandY . Assume the existence of$\gamma \in {\mathbb{R}}^{+}$such that the resolvent ofA over both spaces contains the half-plane Re(λ)≥ − γand satisfies

$$\underset{y\in \mathbb{R}}{sup}\parallel {((-\gamma +\mathrm{i}y)\hspace{0.17em}\mathrm{Id}-A)}^{-1}{\parallel }_{Y\to X}<+\mathrm{\infty }.$$

Then the unique mild solution$w\in C({\mathbb{R}}^{+},Y)$of the initial value problem

$$\frac{\mathrm{d}w}{\mathrm{d}t}=Aw+G,$$

assuming ∥G∥_Y,γ < + ∞ and ∥w_in∥_X < + ∞ forw(0) = w_in, has the following properties

—	w(t)∈Xfor a.e.$t\in {\mathbb{R}}^{+}$
—	∥w∥X,γ ≤ C(∥win∥X + ∥G∥Y,γ) for some $C\in {\mathbb{R}}^{+}$.

In short terms, under a suitable property of the resolvent, asymptotic decay of the solution can be ensured in X, even though control of the forcing only holds in the larger space Y . Now, one checks that the forced linear equation associated with (3.5) satisfies the conditions of this lemma, in appropriate product spaces, with $A={\mathcal{L}}_1+{\mathcal{L}}_2$ [48]. Therefore, its solution satisfies

To conclude the proof, it remains to show, using a localization procedure, that this L²-control in time implies L^∞-control for the original nonlinear equation, see section 5.4 in [48].

4. Convergence to the Ott–Antonsen manifold

In addition to the basic characteristics listed in the Introduction, (1.2) has another remarkable feature. Its solutions asymptotically approach the so-called Ott–Antonsen (OA) manifold, namely the set of measures for which the Fourier transform ${\{{\hat{f}}_{\ell }(\tau )\}}_{(\ell ,\tau )\in \mathbb{N}\times \mathbb{R}}$ satisfies [9],

where g is the frequency marginal associated with f, $h:\mathbb{R}\to \mathbb{C}$ is arbitrary, and * now denotes the convolution on the whole line, i.e. and h^*(ℓ+1) = h*h^*ℓ for $\ell \in {\mathbb{Z}}^{+}$, where h^*0 is the Dirac distribution.

Several motivations for the OA manifold have been given in the literature. As mentioned in the Introduction, the dynamics in this set is a finite-dimensional system when g is meromorphic with finitely many poles in the lower half-plane [9,34] (see also §5.6.2 in [44]). Moreover, this set captures the order parameter dynamics [42,43]. In addition, it selects suitable candidates for PLS stability, namely f_s is the only PLS f_pls contained in this set. Finally, when evaluated in the OA manifold, the corresponding PLS stability condition results to be identical to (2.4) [39,48]. Here, we provide a full proof that the OA manifold is a global attractor for appropriate measures. In order to formulate the statement, we first observe that this set can be regarded as the set of measures for which all functions

identically vanish. Accordingly, the distance to this set will be evaluated using

Proposition 4.1.

Assume thatf(0) is such that$\parallel w(0){\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1({\mathbb{N}}^2\times \mathbb{R})}<+\mathrm{\infty }$for somea > 0 and that the corresponding global solutiont↦f(t) exists. Then$\parallel w(t){\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1({\mathbb{N}}^2\times \mathbb{R})}<+\mathrm{\infty }$for all$t\in {\mathbb{R}}^{+}$and

$$\parallel w(t){\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1({\mathbb{N}}^2\times \mathbb{R})}\le \parallel w(0){\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1({\mathbb{N}}^2\times \mathbb{R})}\hspace{0.17em}{\mathrm{e}}^{-at}.$$

In particular, the following limit holds

$$\underset{t\to +\mathrm{\infty }}{lim}{w}_{n,m}(t,\tau )=0,\mathrm{\forall }n,n\in \mathbb{N},\tau \in \mathbb{R}.$$

Proof.

Using the relations

$${\mathrm{\partial }}_{\tau }{w}_{n,m}={\mathrm{\partial }}_{\tau }{\hat{f}}_{n+m}\ast \hat{g}-{\mathrm{\partial }}_{\tau }{\hat{f}}_n\ast {\hat{f}}_m={\mathrm{\partial }}_{\tau }{\hat{f}}_{n+m}\ast \hat{g}-{\hat{f}}_n\ast {\mathrm{\partial }}_{\tau }{\hat{f}}_m$$

and the Kuramoto dynamics in Fourier space, the following evolutionary equation results:

$${\mathrm{\partial }}_t{w}_{n,m}=(n+m){\mathrm{\partial }}_{\tau }{w}_{n,m}+\frac{K}{2}(\overline{r(t)}(n{w}_{n-1,m}+m{w}_{n,m-1})-r(t)(n{w}_{n+1,m}+m{w}_{n,m+1})).$$

Together with the relations w_0,m≡w_n,0≡0, this equation yields after standard manipulations

$$\begin{array}{rl}\frac{\mathrm{d}}{\mathrm{d}t}\sum _{n,m\in \mathbb{N}}{\int }_{\mathbb{R}}\frac{{\mathrm{e}}^{2a\tau }}{nm}|w_{n,m}(\tau ){|}^2\hspace{0.17em}\mathrm{d}\tau & =2\sum _{n,m\in \mathbb{N}}{\int }_{\mathbb{R}}\frac{(n+m)\hspace{0.17em}{\mathrm{e}}^{2a\tau }}{nm}\mathrm{Re}(w_{n,m}(\tau ){\overline{w^{\prime }}}_{m,n}(\tau ))\hspace{0.17em}\mathrm{d}\tau \\ & =-2a\sum _{n,m\in \mathbb{N}}{\int }_{\mathbb{R}}\frac{(n+m)\hspace{0.17em}{\mathrm{e}}^{2a\tau }}{nm}|w_{n,m}(\tau ){|}^2\hspace{0.17em}\mathrm{d}\tau .\end{array}$$

The same computations hold with w′_n,m instead of w_n,m. Adding the two results, it follows that

$$\frac{\mathrm{d}}{\mathrm{d}t}\parallel w{\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1({\mathbb{N}}^2\times \mathbb{R})}^2\le -2a\parallel w{\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1({\mathbb{N}}^2\times \mathbb{R})}^2,$$

from where the first part of the proposition results. The second part is a direct consequence of the first result together with the Sobolev embedding ${\mathcal{H}}^1([-\tau ,\tau ])\hookrightarrow C_0([-\tau ,\tau ])$, for every $\tau \in {\mathbb{R}}^{+}$. ▪

We do not know if the conditions of the proposition hold for every f(0) such that $\widehat{f(0)}\in {\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1(\mathbb{N}\times \mathbb{R})$. Yet, the next statement provides sufficient conditions for the proposition to apply.

Lemma 4.2.

Suppose that$\parallel \hat{g}{\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a^{\prime }\tau }}^1({\mathbb{R}}^{+})}<+\mathrm{\infty }$and$\parallel \widehat{f(0)}{\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a^{\prime }\tau }}^1(\mathbb{N}\times \mathbb{R})}<+\mathrm{\infty }$for alla′∈[a − ϵ, a + ϵ] whereϵ > 0 is (arbitrarily) small. Then, f(t) satisfies the assumptions of proposition 4.1 for everyt > 0.

Proof.

By splitting the convolution integral into the sum of an integral over ${\mathbb{R}}^{-}$ and one over ${\mathbb{R}}^{+}$, one easily gets the following estimate given any two functions $u,v:\mathbb{R}\to \mathbb{C}$

$$\parallel u\ast v{\parallel }_{L_{{\mathrm{e}}^{a\tau }}^2(\mathbb{R})}^2\le \frac{1}{2(a-a_{-})}\parallel u{\parallel }_{L_{{\mathrm{e}}^{a_{-}\tau }}^2(\mathbb{R})}^2\parallel v{\parallel }_{L_{{\mathrm{e}}^{a_{-}\tau }}^2(\mathbb{R})}^2+\frac{1}{2(a_{+}-a)}\parallel u{\parallel }_{L_{{\mathrm{e}}^{a_{+}\tau }}^2(\mathbb{R})}^2\parallel v{\parallel }_{L_{{\mathrm{e}}^{a_{+}\tau }}^2(\mathbb{R})}^2$$

for every a₋< a < a₊. Using this inequality in straightforward computations based on the definition of w_n,m, together with the notation in (3.6) and the estimate

$$\sum _{n,m\in \mathbb{N}}\frac{u_{n+m}^2}{nm}=\sum _{\ell =2}^{+\mathrm{\infty }}{u}_{\ell }^2\sum _{n=1}^{\ell -1}\frac{1}{n(\ell -n)}\le 2\sum _{\ell =2}^{+\mathrm{\infty }}\frac{u_{\ell }^2(1+\log \ell )}{\ell }\le 2\sum _{\ell =2}^{+\mathrm{\infty }}{u}_{\ell }^2.$$

we obtain

$$\begin{array}{rl}\parallel w(t){\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a\tau }}^1({\mathbb{N}}^2\times \mathbb{R})}^2& \le \frac{2\parallel \hat{g}{\parallel }_{L_{{\mathrm{e}}^{a_{-}\tau }}^2(\mathbb{R})}^2+\parallel \widehat{f(t)}{\parallel }_{L_{{\mathrm{e}}^{a_{-}\tau }}^2(\mathbb{N}\times \mathbb{R})}^2}{a-a_{-}}\parallel \widehat{f(t)}{\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a_{-}\tau }}^1(\mathbb{N}\times \mathbb{R})}^2\\ & +\frac{2\parallel \hat{g}{\parallel }_{L_{{\mathrm{e}}^{a_{+}\tau }}^2(\mathbb{R})}^2+\parallel \widehat{f(t)}{\parallel }_{L_{{\mathrm{e}}^{a_{+}\tau }}^2(\mathbb{N}\times \mathbb{R})}^2}{a_{+}-a}\parallel \widehat{f(t)}{\parallel }_{{\mathcal{H}}_{{\mathrm{e}}^{a_{+}\tau }}^1(\mathbb{N}\times \mathbb{R}).}^2\end{array}$$

The assumptions of the lemma and the fact that the Cauchy problem is well posed in ${\mathcal{H}}_{{\mathrm{e}}^{a^{\prime }\tau }}^1(\mathbb{N}\times \mathbb{R})$ for a′∈{a₋, a₊}⊂[a − ϵ, a + ϵ] [48, proposition 3.1] imply that the second terms in the r.h.s. above are bounded for every t > 0. ▪

5. Existence, stability and bifurcations

This section investigates the connections between the various existence and stability conditions of §2 and discusses their concrete materialization in some examples.

For g symmetric and unimodal, instability of f_hom is equivalent to existence and stability of stationary PLS (figure 1a). In other cases, this connection is not so tight. In particular, stable PLS could exist while f_hom is stable. This happens for instance for the bi-Cauchy distribution

when $\mathit{\Omega }>\mathrm{\Delta }/\sqrt{3}$ ($\mathrm{\Delta }\in {\mathbb{R}}^{+}$) so that the distribution is bimodal. In this case, f_hom is stable for all K ≤ K_c. Yet, stable and unstable stationary PLS f_s co-appear at some K < K_c. Moreover, the unstable PLS branch merges with f_hom via sub-critical bifurcation at K_c (figure 1b and see [34,48] for details). While this example shows that f_hom instability is not necessary for PLS existence, the next statement shows that it is sufficient, provided that globally rotating PLS are allowed.

Proposition 5.1.

Assume thatgis Lipschitz continuous and such that$\parallel \hat{g}{\parallel }_{L_{1+\tau }^1({\mathbb{R}}^{+})}<+\mathrm{\infty }$and that (2.1) fails for somezwith Re(z) > 0. Then, a PLS exists for some frequency$\mathit{\Omega }\in \mathbb{R}$and profile of typef_s.

We can have Ω≠ 0 even though g is symmetric around 0, see the example below.

Proof.

When combined with a suitable Galilean transformation, condition (2.3) immediately yields the following existence condition for PLS with frequency Ω and profile f_s

$$F_r(\mathit{\Omega })=1\mathrm{where}{F}_r(\mathit{\Omega })=\frac{1}{r}{\int }_{\mathbb{R}}\beta \left(\frac{\omega +\mathit{\Omega }}{Kr}\right)g(\omega )\hspace{0.17em}\mathrm{d}\omega .$$

In order to prove the existence of a solution (r, Ω) when f_hom is unstable, notice that

—	F is continuous at every $(r,\mathit{\Omega })\in (0,1]\times \mathbb{R}$, as a consequence of |β( · )| ≤ 1 and Lebesgue dominated convergence.
—	$\underset{\mathit{\Omega }\to \pm \mathrm{\infty }}{lim}\underset{r\in (0,1]}{sup}|F_r(\mathit{\Omega })|=0$ as a consequence of $F_r(\mathit{\Omega })=K{\int }_{\mathbb{R}}\beta (\omega )g(Kr\omega -\mathit{\Omega })\hspace{0.17em}\mathrm{d}\omega$ and dominated convergence again.

Extending F_r by continuity to $\overline{\mathbb{R}}$, the expression ${\{F_r(\mathit{\Omega })\}}_{\mathit{\Omega }\in \overline{\mathbb{R}}}$ defines, for every r∈(0, 1], a closed path in the complex plane. As the next statement reveals, the limit r → 0 also defines a closed path via the quantity involved in (2.1). ▪

Lemma 5.2.

Ifgis Lipschitz continuous, then the limitF₀₊₀(Ω) exists for every$\mathit{\Omega }\in \mathbb{R}$and we have

$$F_{0+0}(\mathit{\Omega })=\frac{K}{2}{\int }_{{\mathbb{R}}^{+}}\hat{g}(\tau )\hspace{0.17em}{\mathrm{e}}^{\mathrm{i}\mathit{\Omega }\tau }\hspace{0.17em}\mathrm{d}\tau ,\hspace{0.17em}\mathrm{\forall }\mathit{\Omega }\in \mathbb{R}.$$

The proof is given below. As argued in [46,51], continuity in Ω and the Riemann–Lebesgue lemma ensure that ${\{F_{0+0}(\mathit{\Omega })\}}_{\mathit{\Omega }\in \overline{\mathbb{R}}}$ is also a closed path. Moreover, these references showed that, assuming $\parallel \hat{g}{\parallel }_{L_{1+\tau }^1({\mathbb{R}}^{+})}<+\mathrm{\infty }$, this path winding number around the point z = 1 is non-zero if (2.1) fails for some z with Re(z) > 0,

On the other hand, the definition of β implies that Re(β(ω)) ≤ 1 for all $\omega \in \mathbb{R}$, with strict inequality when ω≠0. It follows that

The limits F₁( ± ∞) = 0 then imply that ${\{F_1(\mathit{\Omega })\}}_{\mathit{\Omega }\in \overline{\mathbb{R}}}$ winding number around z = 1 must but 0. By the uniform decay above, there must exist r∈(0, 1) for which the path ${\{F_r(\mathit{\Omega })\}}_{\mathit{\Omega }\in \overline{\mathbb{R}}}$ contains z = 1; hence the rotating PLS. The proof of the proposition is complete.

Proof of the lemma.

We rely on the Plemelj formula

$${\int }_{{\mathbb{R}}^{+}}\hat{g}(\tau )\hspace{0.17em}{\mathrm{e}}^{\mathrm{i}\mathit{\Omega }\tau }\hspace{0.17em}\mathrm{d}\tau =\pi g(-\mathit{\Omega })+\mathrm{i}\hspace{0.17em}\mathrm{PV}{\int }_{\mathbb{R}}\frac{g(\omega -\mathit{\Omega })}{\omega }\hspace{0.17em}\mathrm{d}\omega ,$$

and we separate the integral in F_r into the domain |ω| < r^2/3 and |ω|≥r^2/3. In the second domain, we have for small r

$$\frac{1}{r}\beta \left(\frac{\omega }{Kr}\right)=\frac{\mathrm{i}\omega }{K{r}^2}(1-\sqrt{1-{\left(\frac{Kr}{\omega }\right)}^2})=\frac{\mathrm{i}K}{2\omega }+O\left({\left(\frac{Kr}{\omega }\right)}^2\right),$$

and then, using that $g\in L^1(\mathbb{R})$,

$$\underset{r\to 0}{lim}\frac{1}{r}{\int }_{|\omega |\ge r^{2/3}}\beta \left(\frac{\omega }{Kr}\right)g(\omega -\mathit{\Omega })\hspace{0.17em}\mathrm{d}\omega =\frac{\mathrm{i}K}{2}\underset{r\to 0}{lim}{\int }_{|\omega |\ge r^{2/3}}\frac{g(\omega -\mathit{\Omega })}{\omega }\hspace{0.17em}\mathrm{d}\omega =\frac{\mathrm{i}K}{2}\hspace{0.17em}\mathrm{PV}{\int }_{\mathbb{R}}\frac{g(\omega -\mathit{\Omega })}{\omega }\hspace{0.17em}\mathrm{d}\omega .$$

In the first domain, we rely on g being Lipschitz continuous to write g(ω − Ω) = g( − Ω) + ωh (ω) where h is bounded. Then, we have for r small enough

$$\frac{1}{r}{\int }_{|\omega |<r^{2/3}}\beta \left(\frac{\omega }{Kr}\right)g(\omega -\mathit{\Omega })\hspace{0.17em}\mathrm{d}\omega =\frac{g(-\mathit{\Omega })}{r}{\int }_{|\omega |<Kr}\sqrt{1-{\left(\frac{\omega }{Kr}\right)}^2}\hspace{0.17em}\mathrm{d}\omega +\frac{1}{r}{\int }_{|\omega |<r^{2/3}}\beta \left(\frac{\omega }{Kr}\right)\omega h\hspace{0.17em}(\omega )\hspace{0.17em}\mathrm{d}\omega .$$

That h is bounded implies that the second integral vanishes in the limit r → 0. The lemma then follows from the fact that ${\int }_{|x|<1}\sqrt{1-x^2}\hspace{0.17em}\mathrm{d}x=\pi /2$. □

To conclude this section, we provide an example of intricate bifurcation diagram (figure 2b) similar to those reported in the Kuramoto–Sakaguchi model [38,39], but obtained in (1.2) for the tri-Cauchy distribution (figure 2a)

where g_Δ,Ω is the bi-Cauchy distribution defined above. In particular, while the second part of the diagram is reminiscent of the saddle-node bifurcation associated with the bi-Cauchy distribution, the first destabilization scheme of f_hom at K≃1.61 is original and generates a branch of rotating PLS pairs, with frequency Ω and −Ω, respectively. Also, f_hom becomes stable again for K≃2.12 and then suffers a pitchfork bifurcation at K≃2.27, as for a unimodal distribution.

• Download figure
• Open in new tab
• Download PowerPoint

6. Final comments and open questions

Proving asymptotic decay in the Kuramoto PDE has essentially consisted in reducing to a Volterra equation which captures stabilization mechanisms of the linearized dynamics. This approach, and the control of the remaining nonlinear terms, is not limited to the basic model. It can be shown to extend to various extensions such as when f also depends on an additional connectivity parameter $k\in \mathcal{D}$ ($\mathcal{D}\subset {\mathbb{R}}^n$, compact) and the potential writes

where $\beta \in \mathbb{R}$ and $\alpha :\mathcal{D}\times \mathcal{D}\to {\mathbb{R}}^{+}$ is assumed to be Lipschitz continuous.

When the connectivity parameter is irrelevant (i.e. α≡1), the resulting PDE governs the dynamics of empirical measures of the so-called Kuramoto–Sakaguchi model [55]. The very same analysis as in §3 can be developed to obtain stability conditions [38,39] and prove, without any additional conceptual obstacle, asymptotic stability of stationary solutions.

Otherwise, when $\mathcal{D}$ is a finite set, the equation describes the continuum limit of interacting communities of coupled oscillators [56,57]. The analysis repeats for the measure vector ${\{f_k(\mathrm{d}\theta ,\mathrm{d}\omega )\}}_{k\in \mathcal{D}}$, without any difficulty other than having to deal with multi-dimensional Volterra equations for the evolution of the order parameter vector ${\{r_k\}}_{k\in \mathcal{D}}$ with components

More general cases, when $\mathcal{D}$ is infinite, include modelling of networks with random interactions. For arbitrary α, explicit existence and stability conditions might be out of reach, especially for PLS. However, when this function decomposes into a product over individual variables [58,59] as when preferential attachment takes place [60], then a self-consistent Volterra equation holds for the integrated order parameter and the analysis entirely repeats in this case.

Finally, here are two problems that remain unsolved, if not unaddressed:

—

Prove asymptotic stability of other remarkable solutions of the Kuramoto PDE (1.2), such as the standing waves discussed in [61].

—

Prove asymptotic stability of stationary states (and other remarkable states) in extensions of the Kuramoto model for which interactions include several Fourier modes, as in the so-called Daido model [62]. The proof in [51] (inspired from [53]) of Landau damping to fhom straightforwardly extends to this case. However, the problem of asymptotic stability of singular states, such as PLS, remains entirely open.

Data accessibility

This article has no additional data.

Authors' contributions

All authors have contributed to the paper.

Competing interests

We declare we have no competing interests.

Funding

H.D. is supported by Université Sorbonne Paris Cité, in the framework of the ‘Investissements d'Avenir’, convention ANR-11-IDEX-0005, and the People Programme (Marie Curie Actions) of the European Union's Seventh Framework Programme (FP7/2007-2013) under REA grant agreement n. PCOFUND-GA-2013-609102, through the PRESTIGE programme coordinated by Campus France.

Footnotes

Electronic supplementary material is available online at https://dx.doi.org/10.6084/m9.figshare.c.4310738.