Linear stability of the incoherent solution and thetransition formula for the Kuramoto-Daido model
By

Hayato Chiba∗
October, 19, 2009

Abstract
The Kuramoto-Daido model, which describes synchronization phenomena, is a systemof ordinary diﬀerential equations on N -torus deﬁned as coupled harmonic oscillators, whose
natural frequencies are drawn from some distribution function. In this paper, the continuous
model for the Kuramoto-Daido model is introduced and the linear stability of its trivial solution(incoherent solution) is investigated. Kuramoto’s transition point Kc , at which the incoherent
solution changes the stability, is derived for an arbitrary distribution function for natural
frequencies. It is proved that if the coupling strength K is smaller than Kc , the incoherent
solution is asymptotically stable, while if K is larger than Kc , it is unstable.

§ 1.

Introduction

Collective synchronization phenomena are observed in a variety of areas such as
chemical reactions, engineering circuits and biological populations [16]. In order toinvestigate such a phenomenon, Kuramoto [9] proposed a system of ordinary diﬀerential
equations
(1.1)

N
K
dθi
sin(θj − θi ), i = 1, · · · , N,
= ωi +dt
N j=1

where θi ∈ [0, 2π) denotes the phase of an i-th oscillator on a circle, ωi ∈ R denotes
its natural frequency, K > 0 is the coupling strength, and where N is the numberof oscillators. Eq.(1.1) is derived by means of the averaging method from coupled
dynamical systems having limit cycles, and now it is called the Kuramoto model.
It is obvious that when K = 0, θi (t) and θj (t) rotate on a circle at diﬀerent
velocities unless ωi is equal to ωj , and it is true for suﬃciently small K > 0. On the∗ Faculty of Mathematics, Kyushu University, Fukuoka, 819-0395, Japan.

e-mail: chiba@math.kyushu-u.ac.jp

2

Hayato Chibarsynchronization

de-synchronization

Figure 1. The order parameter for the Kuramoto model.

other hand, if K is suﬃciently large, it is numerically observed that some of oscillatorsor all of them tend to rotate at the same velocity on average, which is called thesynchronization [16, 18, 14]. If N is small, such a transition from de-synchronization
to synchronization may be well revealed by means of the bifurcation theory [3, 11, 12].
However, if N is large, it is diﬃcult to investigate the transition from the view point of
the bifurcation theory and it is still far from understood.
In order to evaluate whether synchronization occurs or not, Kuramoto introduced√
the order parameter r(t)e −1ψ(t) by
√

(1.2)

r(t)e −1ψ(t) :=1  √−1θj (t)
e
,
N j=1
N

which gives the centroid of oscillators, where r, ψ ∈ R. It seems that if synchronous
state is formed, r(t) takes a positive number, while if de-synchronization is stable, r(t) iszero on time average (see Fig.1). Based on this observation and some formal calculation,
Kuramoto conjectured a bifurcation diagram of r(t) as follows:
Kuramoto’s conjecture
Suppose that N → ∞ and natural frequencies ωi ’s are distributed according toa probability density function g(ω). If g(ω) is an even and unimodal function, thenthe bifurcation diagram of r(t) is given as Fig.2 (a); that is, if the coupling strength
K is smaller than Kc := 2/(πg(0)), then r(t) ≡ 0 is asymptotically stable. On theother hand, if K is larger than Kc , there exists a positive constant rc such that r(t) =rc is asymptotically stable. Near the transition point Kc , the scaling law of rc is ofO((K − Kc )1/2 ).
Now the value Kc = 2/(πg(0)) is called the Kuramoto’s transition point. See [10]and [18] for the Kuramoto’s discussion.
Signiﬁcant papers of Strogatz et al. [19, 20, 15] partially conﬁrmed the Kuramoto’s
conjecture. Though their arguments are not rigorous from a mathematical view point,Linear stability of the incoherent solution for the Kuramoto-Daido model

3

r

r

Kc
(a)

K

K
(b)

Figure 2. Typical bifurcation diagrams of the order parameter for the cases that (a)g(ω) is even and unimodal (b) g(ω) is even and bimodal. Solid lines denote stable
solutions and dotted lines denote unstable solutions.

almost all of them are justiﬁed as will be done in this paper. In [20], they introduced
the continuous model for the Kuramoto model and investigated the linear stability of atrivial solution called the incoherent solution, which corresponds to the de-synchronous
state r ≡ 0. They derived the Kuramoto’s transition point Kc = 2/(πg(0)) and showedthat if K > Kc , the incoherent solution is unstable in the linear level (i.e. nonlinear
terms are neglected). When K < Kc , the linear operator T , which deﬁnes the linearized
equation of the continuous model around the incoherent solution, has no eigenvalues.
However, in [19], they found that an analytic continuation of the resolvent (λ−T )−1 may
have poles (resonance poles) on the left half plane, and they remarked a possibility that
resonance poles induce exponential decay of the order parameter. In [15], the stabilityof the partially locked state, which corresponds to a solution with positive constantr = rc , is investigated in the linear level.
Despite the active interest in the case that the distribution function g(ω) is evenand unimodal, bifurcation diagrams of r for g(ω) other than the even and unimodal
cases are not revealed well. Martens et al. [13] investigated the bifurcation diagram fora bimodal g(ω) which consists of two Lorentzian distributions. In particular, they foundthat stable partially locked states can coexist with stable incoherent solutions if K isslightly smaller than Kc (see Fig.2 (b)). Such a diagram seems to be common for any
bimodal distributions.
A simple extension of the Kuramoto model deﬁned to be
(1.3)

N
K
dθi
f (θj − θi ), i = 1, · · · , N,
= ωi +dt
N j=1

is called the Kuramoto-Daido model [4, 5, 6, 7], where the 2π-periodic function f :R → R is called the coupling function. Daido [7] investigated bifurcation diagrams of
the order parameter for the Kuramoto-Daido model with even and unimodal g(ω) by4Hayato Chiba

a similar argument to Kuramoto’s one. He found that if f (θ) = sin θ, partially locked
states may coexist with stable incoherent solutions even if g(ω) is even and unimodal.
All such studies by physicist are based on formal calculations and numerical simulations. The purpose of this paper is to justify and extend their results as mathematics
for the Kuramoto-Daido model with any distribution function g(ω). The continuous
model for the Kuramoto-Daido model is introduced and the linear stability of the incoherent solution is studied. In particular, the spectrum and the semigroup of a linear
operator T , which is obtained by linearizing the continuous model around the incoherent
solution, will be investigated in detail. At ﬁrst, a formula for obtaining the transition
point Kc for an arbitrary distribution g(ω) is derived. As a corollary, the Kuramoto’s
transition point Kc = 2/(πg(0)) is obtained if g(ω) is an even and unimodal function. If
K > Kc , it is proved that the incoherent solution is unstable because the operator T has
eigenvalues on the right half plane. It means that if the coupling strength K is large,
the de-synchronous state is unstable and thus synchronization may occur. On the other
hand, if 0 < K < Kc , it will be shown that the spectrum of the operator T consistsof the continuous spectrum and it lies on the imaginary axis. Thus the stability of theincoherent solution is nontrivial. Despite this fact, under appropriate assumptions forg(ω), the order parameter proves to decay exponentially because of existence of resonance poles on the left half plane as was expected by Strogatz et al. [19]. It suggeststhat in general, linear stability of a trivial solution of a linear equation on an inﬁnite
dimensional space is determined by not only the spectrum of the linear operator but
also its resonance poles.
§ 2.

Continuous model

In this section, we introduce a continuous model of the Kuramoto-Daido model andshow a few properties of it.
Let us consider the Kuramoto-Daido model (1.3). We suppose that the coupling
function f is a C 1 periodic function with the period 2π. It is expanded in a Fourier
series as∞√fl e −1lθ , fl ∈ C.
(2.1)
f (θ) =
l=−∞

We can suppose that f0 = 0 without loss of generality because f0 is renormalized into√
the constants ωi . For the Kuramoto model (f (θ) = sin θ), f±1 = ±1/(2 −1) andfl = 0 (l = ±1). Following Daido [7], we introduce the generalized order parameters Ẑk0by1  √−1kθj (t)e, k = 0, ±1, ±2, · · · .
N j=1
N

(2.2)

Ẑk0 (t) =Linear stability of the incoherent solution for the Kuramoto-Daido model5

In particular, Ẑ10 is the order parameter deﬁned in Section 1. By using them, Eq.(1.3)
is rewritten as
∞

(2.3)

√

dθi
fl Ẑl0 (t)e− −1lθi .
= ωi + K
dtl=−∞

Motivated by these equations, we introduce a continuous model of the Kuramoto-Daidomodel, which is an evolution equation of a probability density function ρt = ρt (θ, ω) on
S 1 × R parameterized by t ∈ R, as


∞


√


∂
∂ρ
t

0
− −1lθ
ρt = 0,
+
ω
+
K
f
Z
(t)e

l
l

∂θ

 ∂t
(2.4)

l=−∞
2π √
e −1kθ ρt (θ, ω)dθdω,
Zk0 (t) =



R 0



ρ0 (θ, ω) = h(θ, ω),

where h(θ, ω) is an initial density function. The Zk0 (t) is a continuous version of Ẑk0 (t),and we also call it the generalized order parameter. We can prove that Eq.(2.4) is properin the sense that Ẑk0 (t) → Zk0 (t) as N → ∞ under some assumptions, although the proof
is not given in this paper. If we regard
v := ω + K∞
√fl Zl0 (t)e− −1lθl=−∞

as a velocity ﬁeld, Eq.(2.4) provides an equation of continuity ∂ρt /∂t + ∂(ρt v)/∂θ = 0
known in ﬂuid dynamics. It is easy to prove the low of conservation of mass:
2π

(2.5)
R 0

2π

ρt (θ, ω)dθ =h(θ, ω)dθ =: g(ω).
R 0

A function g deﬁned as above gives a probability density function for natural frequencies
ω ∈ R such that R g(ω)dω = 1.
By using the characteristic curve method, Eq.(2.4) is formally integrated as follows:
Consider the equation
∞

(2.6)

√

dx
fl Zl0 (t)e− −1lx , x ∈ S 1 ,
=ω+K
dt
l=−∞

which deﬁnes a characteristic curve. Let x = x(t, s; θ, ω) be a solution of Eq.(2.6)
satisfying x(s, s; θ, ω) = θ. Then, ρt is given as(2.7)

ρt (θ, ω) = h(x(0, t; θ, ω), ω) exp K∞t 
0 l=−∞√√−1 l fl Zl0 (s)e− −1lx(s,t;θ,ω) ds .

6Hayato Chiba

By using Eq.(2.7), it is easy to show the equality
2π

(2.8)
R 0

2π

a(θ, ω)ρt (θ, ω)dθdω =a(x(t, 0; θ, ω), ω)h(θ, ω)dθdω,R 0

for any continuous function a(θ, ω). In particular, the generalized order parameters
Zk0 (t) are rewritten as
Zk0 (t) =

(2.9)

2π √
−1kx(t,0;θ,ω)eh(θ, ω)dθdω.

R 0

Substituting it into Eqs.(2.6), (2.7), we obtain(2.10)

2πdx(t, s; θ, ω) = ω + K
dt

R 0

f (x(t, 0; θ , ω  ) − x(t, s; θ, ω))h(θ , ω  )dθ dω  ,andρt (θ, ω) = h(x(0, t; θ, ω), ω) ×
2πt

(2.11)

exp K
0

ds ·

R 0


∂f 
x(s, 0; θ , ω  ) − x(s, t; θ, ω) h(θ , ω  )dθ dω  ,∂θrespectively. Even if h(θ, ω) is not diﬀerentiable, we consider Eq.(2.11) to be a weaksolution of Eq.(2.4). It is easy in usual way to prove that the integro-ODE (2.10) has aunique solution for any t > 0, and this proves that the continuous model Eq.(2.4) has aunique weak solution (2.11) for an arbitrary initial data h(θ, ω).
Throughout this paper, we suppose that the initial date h(θ, ω) is of the formh(θ, ω) = ĥ(θ)g(ω). This assumption corresponds to the assumption for the KuramotoN
Daido model (1.3) that initial values {θj (0)}N
j=1 and natural frequencies {ωj }j=1 are
independently distributed. This is a physically natural assumption used in many literatures. In this case, ρt (θ, ω) is written as ρt (θ, ω) = ρ̂t (θ, ω)g(ω), whereρ̂t (θ, ω) = ĥ(x(0, t; θ, ω)) ×
2π

t(2.12)

exp K
0

ds ·

R 0


∂f 
x(s, 0; θ , ω  ) − x(s, t; θ, ω) ĥ(θ )g(ω  )dθ dω  ,∂θand ρ̂t (θ, ω) satisﬁes the same equation as Eq.(2.4).
§ 3.

Linear stability of the incoherent solution

A trivial solution of the continuous model (2.4), which is independent of θ andt, is given by ρt (θ, ω) = g(ω)/(2π), or equivalently ρ̂t (θ, ω) = 1/(2π). It is called theincoherent solution, which corresponds to the de-synchronized state. Note that in thisLinear stability of the incoherent solution for the Kuramoto-Daido model

7

case r = 0. In this section, we investigate the stability of the incoherent solution andthe order parameter.
Let
2π √

(3.1)

Zj (t, ω) :=0

e −1jθ ρ̂t (θ, ω)dθ

be the Fourier coeﬃcients of ρ̂t (θ, ω). Then, Z0 (t, ω) = 1 and Zj , j = ±1, ±2, · · · satisfy
the diﬀerential equations
∞

√
dZj √
fl Zl0 (t)Zj−l= −1jωZj + −1jK
dt
−∞

√
√
√
= −1jωZj + −1jKfj Zj0 (t) + −1jK
fl Zl0 (t)Zj−l .
l=jThe incoherent solution corresponds to the zero solution Zj ≡ 0 for j = ±1, ±2, · · · .
Since |Zj (t, ω)| ≤ 1, Zj (t, ω) is in the Hilbert space L2 (R, g(ω)dω) for every t :||Zj ||2L2 (R,g(ω)dω) =R|Zj (t, ω)|2 g(ω)dω ≤ 1.

Thus we linearize the above equation as an evolution equation on L2 (R, g(ω)dω)
(3.2)

√

 √
dZj
= j −1M + j −1Kfj P Zj , j = ±1, ±2, · · · ,dtwhere M : q(ω) → ωq(ω) is the multiplication operator on L2 (R, g(ω)dω) and P is the
projection on L2 (R, g(ω)dω) deﬁned to be(3.3)

Pq(ω) =q(ω)g(ω)dω.
R

If we put P0 (ω) ≡ 1, P is also written as Pq(ω) = (q, P0 ), where ( , ) is the innerproduct on L2 (R, g(ω)dω):(3.4)(q1 , q2 ) :=R

q1 (ω)q2 (ω)g(ω)dω.

Note that the order parameter is given as Z10 = (Z1 , P0 ). To determine the stability of
the incoherent solution and the order parameter, we have to investigate the spectrum√√
and the semigroup of the operator Tj := j −1M + j −1Kfj P.
√§ 3.1. Analysis of the operator −1M√
If fj = 0, Tj = j −1M. It is known that the multiplication operator M on
L2 (R, g(ω)dω) is self-adjoint and its spectrum is given by σ(M) = supp(g) ⊂ R, where√supp(g) is a support of the density function g. Thus the spectrum of j −1M is
√
√
√
√
(3.5)
σ(j −1M) = j −1 · supp(g) = {j −1λ | λ ∈ supp(g)} ⊂ −1R.

8Hayato Chiba√√√√
The semi-group ej −1Mt generated by j −1M is given as ej −1Mt q(ω) = ej −1ωt q(ω).
In particular, we obtain√(ej −1Mt q1 , q2 ) =

(3.6)

√

R

ej −1ωt q1 (ω)q2 (ω)g(ω)dω

for any q1 , q2 ∈ L2 (R, g(ω)dω). This is the Fourier transform of the function q1 (ω)q2 (ω)g(ω).
Thus if q1 (ω)q2 (ω)g(ω) is real analytic on R and has an analytic continuation to a neigh√borhood of the real axis, then (ej −1Mt q1 , q2 ) decays exponentially as t → ∞, while if
q1 (ω)q2 (ω)g(ω) is C r , then it decays as O(1/tr ) (see Vilenkin [21]).
These facts are summarized as follows:

Proposition 3.1.
Suppose that fj = 0 and Eq.(3.2) is reduced to dZj /dt =j −1MZj . A solution of this equation with an initial value q(ω) ∈ L2 (R, g(ω)dω)√√
is given by Zj (t) = ej −1Mt q(ω) = ej −1ωt q(ω). In particular the linearized order√parameter Z10 (t) = (e −1Mt q, P0 ) decays exponentially as t → ∞ if g(ω) and q(ω) haveanalytic continuations to a neighborhood of the real axis.
√√
The resolvent (λ − j −1M)−1 of the operator j −1M is calculated as
√

(3.7)

√
((λ − j −1M)−1 q1 , q2 ) =Rλ − j1√−1ω

q1 (ω)q2 (ω)g(ω)dω.

We deﬁne the function D(λ) to be(3.8)√
D(λ) = ((λ − j −1M)−1 P0 , P0 ) =Rλ − j1√−1ω

g(ω)dω

√(recall that P0 (ω) ≡ 1). It is holomorphic in C\σ(j −1M) and will play an importantrole in the later calculation.
§ 3.2.

√√
Analysis of the operator Tj = j −1M + j −1Kfj P

In what follows, we suppose that fj = 0. The domain D(Tj ) of Tj is given by
D(M) ∩ D(P) = D(M). Since M is self-adjoint and since P is bounded, Tj is a closed
operator [8]. Let (Tj ) be the resolvent set of Tj and σ(Tj ) = C\ (Tj ) the spectrum.
Since Tj is closed, there is no residual spectrum. Let σp (Tj ) and σc (Tj ) be the point
spectrum (the set of eigenvalues) and the continuous spectrum of Tj , respectively.

Proposition 3.2.
(3.9)(i) Eigenvalues λ of Tj are given as roots of√1
D(λ) = √, λ ∈ C\σ(j −1M).
j −1KfjLinear stability of the incoherent solution for the Kuramoto-Daido model

9(ii) The continuous spectrum of Tj is given by√√(3.10)σc (Tj ) = σ(j −1M) = j −1 · supp(g).
√Proof.
(i) Suppose that λ ∈ σp (Tj )\σ(j −1M). Then, there exists x ∈L2 (R, g(ω)dω) such that√√
λx = (j −1M + j −1Kfj P)x, x = 0.
√√
Since λ ∈/ σ(j −1M), (λ − j −1M)−1 is deﬁned and the above is rewritten as√√
x = (λ − j −1M)−1 j −1Kfj Px√√= j −1Kfj (x, P0 )(λ − j −1M)−1 P0 (ω).
By taking the inner product with P0 (ω), we obtain√√√(3.11)
1 = j −1Kfj ((λ − j −1M)−1 P0 , P0 ) = j −1Kfj D(λ).
This proves that roots of Eq.(3.9) is in σp (Tj ). The corresponding eigenvector is given√√√
by x = (λ − j −1M)−1 P0 (ω) = 1/(λ − j −1ω). If λ ∈ −1R, x ∈/ L2 (R, g(ω)dω).
Thus there are no eigenvalues on the imaginary axis.
(ii) This follows from the fact that the essential spectrum is stable under the bounded√
perturbation and that there are no eigenvalues on σ(j −1M), see [8].
§ 3.3.

Eigenvalues of the operator Tj and the transition point formula

Our next task is to calculate roots of Eq.(3.9) to obtain eigenvalues of Tj =√√j −1M + j −1Kfj P. By putting λ = x + −1y, x, y ∈ R, Eq.(3.9) is rewrittenas
Im(fj )x
 x2 + (jω − y)2 g(ω)dω = − jK|fj |2 ,R
(3.12)
Re(fj )jω − yg(ω)dω = −.
2
2
jK|fj |2
R x + (jω − y)√In what follows, we suppose that Im(fj ) < 0. The case Im(fj ) ≥ 0 will be treated inSec.3.5. The next lemma is easily obtained.

Lemma 3.3.
(i) When Im(fj ) < 0, λ satisﬁes Re(λ) > 0 for any K > 0.
(ii) If K > 0 is suﬃciently large, there exists at least one eigenvalue λ near inﬁnity.
(iii) If K > 0 is suﬃciently small, there are no eigenvalues.

10

Hayato Chiba

Proof. Part (i) of the lemma immediately follows from the ﬁrst equation of Eq.(3.12).
To prove part (ii) of the lemma, note that if |λ| is large, Eq.(3.9) is rewritten as1
11.
+ O( 2 ) = √λλj −1Kfj√
Thus the Rouché’s theorem proves that Eq.(3.9) has a root λ ∼ j −1Kfj if K > 0
is suﬃciently large. To prove part (iii) of the lemma, we see that the left hand side of
the ﬁrst equation of Eq.(3.12) is bounded for any x, y ∈ R. To do so, let G(ω) be theprimitive function of g(ω) and ﬁx δ > 0 small. The left hand side of the ﬁrst equationof Eq.(3.12) is calculated asxg(ω)dωx2 + (jω − y)2

R
∞

=
=

y/j−δ

xg(ω)dω+2
2y/j+δ x + (jω − y)

−∞

∞

y/j−δ

xg(ω)dω+2
2y/j+δ x + (jω − y)−∞xg(ω)dω+2
x + (jω − y)2

y/j+δ
y/j−δ

xg(ω)dωx2 + (jω − y)2xg(ω)dωx2 + (jω − y)2x+ 2(G(y/j + δ) − G(y/j − δ)) +x + j 2 δ2

2jx(jω − y)

y/j+δ
y/j−δ

(x2 + (jω − y)2 )2G(ω)dω.

The ﬁrst three terms in the right hand side above are bounded for any x, y ∈ R. Since
G is continuous, there exists a number ξ such that the last term is rewritten as
y/j+δ
y/j−δ

2jx(jω − y)(x2 + (jω − y)2 )2G(ω)dω = 2jδ ·2xξ
(x2 + ξ 2 )2G(y/j + ξ/j).

This is bounded for any x, y ∈ R. Now we have proved that the left hand side of theﬁrst equation of Eq.(3.12) is bounded for any x > 0, although the right hand side tends
to inﬁnity as K → +0. Thus Eq.(3.9) has no roots if K is small.
Lemma 3.3 shows that if K > 0 is suﬃciently large, the trivial solution Zj = 0
of the system dZj /dt = Tj Zj is unstable because of the eigenvalues with positive real(j)parts. Our purpose in this subsection is to determine the bifurcation point Kc , which(j)is the minimum value of K such that if K < Kc , the operator Tj has no eigenvalues on
the right half plane. To calculate eigenvalues λ = λ(K) explicitly is diﬃcult in general.
However, note that since zeros of a holomorphic function do not vanish because of theargument principle, λ(K) disappears if and only if it is absorbed into the continuous√
spectrum σ(j −1M), on which D(λ) is not holomorphic. This fact suggests that to(j)determine Kc , it is suﬃcient to investigate Eq.(3.9) or Eq.(3.12) near the imaginary√√
axis. Since we are interested in λ(K) absorbed into σ(j −1M) ⊂ −1R, take theLinear stability of the incoherent solution for the Kuramoto-Daido model

11

y3
y2
y1

㱗1 (K)
supp (g)Figure 3. A schematic view of behavior of roots λ of Eq.(3.9) when K decreases.

limit x → +0 in Eq.(3.12):



lim

 x→+0(3.13)




 lim

x→+0

x
R

x2 + (jω − y)2Rx2 + (jω − y)2jω − yg(ω)dω = −Im(fj ),
jK|fj |2g(ω)dω = −Re(fj ).
jK|fj |2

These equations determine Kn and yn such that one of the eigenvalues λn (K) converges√
to −1yn as K → Kn + 0 (see Fig.3). To calculate them, we need the next lemma.
Lemma 3.4.
(i) Suppose that λn (K) →continuous at ω = yn .
(ii) If g(ω) is continuous at ω = y, then

(3.14)

x

lim

x→+0

R

x2 + (jω − y)2√−1yn as K → Kn . Then, g(ω) is

g(ω)dω = πg(y/j)/j.

Proof. To prove (i), suppose that g(ω) is discontinuous at ω = 0 without loss ofgenerality.
STEP 1: At ﬁrst, we suppose that g(ω) is piecewise continuous. Put g(+0) =h+ , g(−0) = h− and h+ = h− . In this case, for any ε > 0, there exists δ > 0 suchthat if −δ < ω < 0, then |g(ω) − h− | < ε and if 0 < ω < δ, then |g(ω) − h+ | < ε. For√
Eq.(3.9), we suppose |λ| = |x + −1y| < δ and y > 0. The case y < 0 is treated in a

12Hayato Chiba

similar manner. We calculate D(λ) as∞

D(λ) =
δ

g(ω)
√
dω +
λ − j −1ω
δ

+
0
∞

=
δ

+
δ

+
0

∞
δ

δ

(g(ω) − h+ ) dω +1√(h+ − h− ) dω + h−

λ − j −1ω0−δ λ − j1√−1ω(g(ω) − h− + h− ) dω

g(−ω)
√
dω
λ + j −1ω1√

λ − j −1ω

0

g(ω)
√
dω
−∞ λ − j −1ω(g(ω) − h+ + h+ ) dω +λ − j −1ω
g(ω)
√
dω +
λ − j −1ωδ

(3.15)

1
√

−δ

0

1
√

λ + j −1ωδ(g(−ω) − h− ) dω

dω
√
+ h−
0 λ − j −1ω

δ

dω
√
.
0 λ + j −1ω

Since |λ| < δ, there exists a positive number M , which is independent of λ, such that
∞
δ

g(±ω)
√
dω < M.
λ ∓ j −1ωThus |D(λ)| is estimated as
δ
δ
δ
dω
dω
dω
√
√
√
−ε
−ε
|D(λ)| ≥ |h+ − h− |0 |λ − j −1ω|0 |λ − j −1ω|0 |λ + j −1ω|

√

 h−
−1λ + jδ 

−
log √
− 2M.
j
−1λ − jδ √√
Since y > 0, |λ − j −1ω| < |λ + j −1ω|. This shows that

√

δ
 h−

dω
−1λ
+
jδ
 − 2M.
√log √(3.16) |D(λ)| ≥ (|h+ − h− | − 2ε)− j−1λ − jδ 
0 |λ − j −1ω|The right hand side tends to inﬁnity as λ → 0 if 2ε < |h+ − h− |. This proves thatEq.(3.9) has no roots at λ = 0 for positive K.
STEP 2: In general, since g(ω) is a non-negative measurable function, there exists
a monotonic increasing sequence {gn (ω)}∞
n=1 of non-negative simple functions such thatgn (ω) → g(ω) for each ω. In particular if g(ω) is discontinuous at ω = 0, we can choose{gn (ω)}∞n=1 so that gn (ω) is discontinuous at ω = 0 for any n ∈ N. Then, the proof isdone in the same way as STEP 1 by approximating g(ω) by gn (ω).
(ii) The formula Eq.(3.14) is proved in Ahlfors [1].
Let (y, K) be one of the solutions of Eq.(3.13). Since g(ω) is continuous at ω = y,substituting it into the ﬁrst equation of Eq.(3.13) yields(3.17)

πg(y/j) = −Im(fj ).
K|fj |2

Linear stability of the incoherent solution for the Kuramoto-Daido model

13

Substituting K = −Im(fj )/(π|fj |2 g(y/j)) obtained from the above into the second
equation of Eq.(3.13) results in(3.18)jω − y

lim

x→0

R

x2 + (jω − y)2g(ω)dω =πRe(fj )g(y/j).
j Im(fj )This equation for y determines imaginary parts to which λ(K) converges as Re(λ(K)) →+0. Let y1 , y2 , · · · be roots of Eq.(3.18). Then,(3.19)Kn =−Im(fj ), n = 1, 2, · · ·π|fj |2 g(yn /j)give the values such that Re(λ(K)) → 0 as K → Kn + 0. Now we obtain the next
theorem.

Theorem 3.5.

Suppose that Im(fj ) < 0. Let y1 , y2 , · · · be roots of Eq.(3.18).

Put
Kc(j) := inf Kn =(3.20)n−Im(fj ).
2
π|fj | supn g(yn /j)(j)

If 0 < K < Kc , the operator Tj has no eigenvalues, while if K is slightly larger than(j)
Kc , Tj has eigenvalues on the right half plane.
Note that inf n Kn is positive because of Lemma.3.3 (iii). As a corollary, we obtainthe transition point (bifurcation point to the partially locked state) conjectured by
Kuramoto [10]:Corollary 3.6 (Kuramoto’s transition point).
Suppose that the probability density function g(ω) is even and max g(ω) = g(0). If Re(f1 ) = 0 and Im(f1 ) = −1/2 (it(1)
corresponds to f (θ) = sin θ in Eq.(1.3)), then Kc deﬁned as above is given by
Kc(1) =

(3.21)

2
.
πg(0)

(1)

When K > Kc , the solution Z1 = 0 of Eq.(3.2) is unstable.
§ 3.4.

Semi-group generated by the operator T1 (Im(f1 ) < 0)

Since we are interested in the dynamics of the order parameter Z10 (t) = (Z1 , P0 ),in what follows, we consider only j = 1 while cases j = 2, 3, · · · are investigated in(1)
the same way. Theorem 3.5 shows that Kc is the least bifurcation point and the

14

Hayato Chiba(1)trivial solution Z1 (t, ω) = 0 of Eq.(3.2) is unstable if K is slightly larger than Kc . If√
(1)
0 < K < Kc , the spectrum of T1 is on the imaginary axis: σ(T1 ) = σ( −1M), andthus the dynamics of Z1 is nontrivial. In this subsection, we investigate the dynamics(1)
of Z1 and the order parameter for 0 < K < Kc . We will see that the order parameter
may decay exponentially even if the spectrum lies on the imaginary axis because ofexistence of resonance poles.
√√
Since −1M has the semi-group e −1Mt and since P is bounded, the operator√√
T1 = −1M+ −1Kf1 P also generates the semi-group (Kato [8]), say eT1 t . A solutionof Eq.(3.2) with an initial value q(ω) ∈ L2 (R, g(ω)dω) is given by eT1 t q(ω). The eT1 t iscalculated by using the Laplace inversion formula
T1 t(3.22)e= lim

y→∞ 2π1√√
x+ −1y−1√
x− −1yeλt (λ − T1 )−1 dλ,where x > 0 is chosen so that the contour is to the right of the spectrum of T1(Yosida [22]). At ﬁrst, let us calculate the resolvent (λ − T1 )−1 .

For any q(ω) ∈ L2 (R, g(ω)dω), the equality√((λ − −1M)−1 q, P0 )−1√
F0 (λ) := ((λ − T1 ) q, P0 ) =1 − −1Kf1 D(λ)

Lemma 3.7.
(3.23)
holds.
Proof. Put

R(λ)q := (λ − T1 )−1 q = (λ −

√

−1M −

√

−1Kf1 P)−1 q,

which yields(λ −√−1M)R(λ)q = q +
=q+

√
√

−1Kf1 PR(λ)q−1Kf1 (R(λ)q, P0 )P0 .

This is rearranged as
R(λ)q = (λ −√√√−1M)−1 q + −1Kf1 (R(λ)q, P0 )(λ − −1M)−1 P0 .

By taking the inner product with P0 , we obtain(R(λ)q, P0 ) = ((λ −This proves Eq.(3.23).

√−1M)−1 q, P0 ) +√−1Kf1 (R(λ)q, P0 )D(λ).

Linear stability of the incoherent solution for the Kuramoto-Daido model

15

Let Z10 (t) = (Z1 , P0 ) be the order parameter with the initial condition Z1 (0, ω) =q(ω). Eqs.(3.22) and (3.23) show that Z10 (t) is given by√√x+ −1y1((λ−−1M)−1 q, P0 )λt√√
(3.24)
Z10 (t) = (eT1 t q, P0 ) = limedλ.
y→∞ 2π −1 x−√−1y
1 − −1Kf1 D(λ)
One of the eﬀective way to calculate the integral above is to use the residue theorem.
Recall that the resolvent (λ − T1 )−1 is holomorphic on C\σ(T1 ). Since we assume that(1)
0 < K < Kc , T1 has no eigenvalues and the continuous spectrum lies on the imaginary√√
axis : σ(T1 ) = σ( −1M) = −1 · supp(g). Thus the integrand eλt F0 (λ) in Eq.(3.24)
is holomorphic on the right half plane and may not be holomorphic on σ(T1 ). However,under assumptions below, we can show that F0 (λ) has an analytic continuation F1 (λ)
through the line σ(T1 ) from right to left. Then, F1 (λ) may have poles on the left halfplane (the second Riemann sheet of the resolvent), which are called resonance poles [17].
The resonance pole µ aﬀects the integral in Eq.(3.24) through the residue theorem (seeFig.4). In this manner, the order parameter Z10 (t) can decay with the exponential rate
Re(µ). Such an exponential decay caused by resonance poles is well known in the theoryof Schrödinger operators [17], and for the Kuramoto model, it is investigated numericallyby Strogatz et al. [19] and Balmforth et al. [2].
At ﬁrst, we construct an analytic continuation of the function F0 (λ).

Lemma 3.8.
Suppose that the probability density function g(ω) and an initialcondition q(ω) are real analytic on R. If g(ω) and q(ω) have meromorphic continuationsg ∗ (λ) and q ∗ (λ) to the upper half plane, respectively, then the function F0 (λ) deﬁned
on the right half plane has the meromorphic continuation F1 (λ) to the left half plane,
which is given by√√√((λ − −1M)−1 q, P0 ) + 2πq ∗ (− −1λ)g ∗ (− −1λ)√√√.
(3.25)F1 (λ) =1 − −1Kf1 D(λ) − 2π −1Kf1 g ∗ (− −1λ)
Proof. By the formula (3.14), we obtain(3.26)√√lim ((λ− −1M)−1 q, P0 )− lim ((λ− −1M)−1 q, P0 ) = 2πq(Im(λ))·g(Im(λ)).

Re(λ)→+0

Re(λ)→−0√Thus the meromorphic continuation of ((λ − −1M)−1 q, P0 ) from right to left is givenby√((λ − −1M)−1 q, P0 )(Re(λ) > 0),√√√(3.27)−1∗∗((λ − −1M) q, P0 ) + 2πq (− −1λ)g (− −1λ) (Re(λ) < 0).
This proves Eq.(3.25).

16

Hayato Chiba

Poles of F1 (λ) (resonance poles) on the left half plane are given as roots of theequation√1, Re(λ) < 0(3.28)
D(λ) + 2πg ∗ (− −1λ) = √−1Kf1√
and poles of the function q ∗ (− −1λ). In the next theorem, we suppose for simplicity√
that q ∗ (− −1λ) has no poles. Now we calculate the order parameter Z10 (t).

Theorem 3.9.
For Eq.(3.2) with j = 1, suppose that
(1)(i) Im(f1 ) < 0 and 0 < K < Kc .
(ii) the probability density function g(ω) is real analytic on R and has a meromorphic
continuation g ∗ (λ) to the upper half plane.
(iii) an initial condition q(ω) is real analytic on R and has an analytic continuationq ∗ (λ) to the upper half plane.
(iv) there exists a positive number δ such that |F1 (λ)| → 0 as |λ| → ∞ in the angulardomains|arg(λ)| ≤ δ, |arg(λ) − π| ≤ δ.

(3.29)(v) there exist positive constants D and β such that
|F1 (λ)| ≤ Deβ|λ|(3.30)

in the angular domain π/2 + δ ≤ arg(λ) ≤ 3π/2 − δ.
Then, there exist resonance poles of T1 on the left half plane. Let α1 , α2 , · · · be resonancepoles such that |α1 | ≤ |α2 | ≤ · · · . Then, there exists a positive constant t0 such that theorder parameter is given by(3.31)

Z10 (t) = (eT1 t q, P0 ) =∞

pn (t)eαn t ,t > t0n=1

where pn (t) is a polynomial in t. In particular, Z10 (t) decays exponentially as t → ∞.
Proof. At ﬁrst, we prove the existence of resonance poles. Resonance poles areroots of Eq.(3.28), which is the analytic continuation of the equation (3.9) for j = 1.
Thus one of the resonance poles is obtained as a continuation of an eigenvalue λ(K).
(1)
Recall that λ(K) converges into the imaginary axis as K → Kc + 0. To prove that(1)
there exists a resonance pole on the left half plane when K < Kc , we have to show
(1)that λ(K) does not stay on the imaginary axis for K < Kc . Diﬀerentiating Eq.(3.9)
with respect to K, we obtain(3.32)λ (K)

1
1√√
g(ω)dω
=,−1ω)2−1K 2 f1
R (λ −

Linear stability of the incoherent solution for the Kuramoto-Daido model

C3

17

C2
r

C4

C1

r

d

C5

C6

Figure 4. The contour for the Laplace inversion formula.

which proves that λ (K) = 0. Further, roots y of Eq.(3.18), which determines eigenvalues on the imaginary axis, are isolated because both side of Eq.(3.18) are analytic withrespect to y. This means that λ(K) can not move along the imaginary axis. This provesthat an eigenvalue λ(K) gets across the imaginary axis from right to left as K decreases(1)from Kc , which gives a root of Eq.(3.28). Note that there may exist resonance poles
which are not continuations of eigenvalues (see Example 3.11).
Next, let us prove Eq.(3.31). Let d > 0 be a small number and r suﬃciently large
number. Take paths C1 to C6 as are shown in Fig.4:
C1 = {d +√−1y | − r ≤ y ≤ r},√
C2 = {x + −1r | 0 ≤ x ≤ d},√

C3 = {re −1θ | π/2 ≤ θ ≤ π/2 + δ},√

C4 = {re −1θ | π/2 + δ ≤ θ ≤ 3π/2 − δ},and C5 and C6 are deﬁned in a similar way to C3 and C2 , respectively. We put C =6
j=1 Cj .
Let α1 , α2 , · · · , αn(C) be resonance poles inside the closed curve C, where we assume
that there are no resonance poles on the curve C by deforming it slightly if necessary.
Let R1 (t), R2 (t), · · · , Rn(c) (t) be corresponding residues of eλt F1 (λ), respectively. Notethat if αj is a pole of F1 (λ) of order mj , Rj (t) is of the form Rj (t) = pj (t)eαj t with a18

Hayato Chiba

polynomial pj (t) of degree mj − 1. By the residue theorem, we have
√

2π −1



n(C)

Rj (t) =j=1

C6 +C1 +C2

eλt F0 (λ)dλ +C3 +C4 +C5eλt F1 (λ)dλ.

√
The integral C1 eλt F0 (λ)dλ/(2π −1) converges to Z10 (t) as r → ∞. It is easy to showthat the integrals along C2 , C3 , C5 , C6 tend to zero as r → ∞ because of the assumption(iv). We have to estimate the integral along C4 as
3π/2−δ
√
λt
rt cos θ
−1θ


e
F
(λ)dλ
≤
re
|F
(re
)|dθ
1
1


C4

π/2+δ

≤

√

max

π/2+δ≤θ≤3π/2−δ
π/2≤ Deβr|F1 (re −1θ )|3π/2−δ

rert cos θ dθ

π/2+δ2re−rt sin φ dφδπ/2

≤ Deβr2re−2rtφ/π dφ≤ Deπ  −2rtδ/π
−rt
.
e
·
−e
t
δ

βr

(3.33)

Thus if t > t0 := max {β, πβ/(2δ)}, this integral tends to zero as r → ∞.
Example 3.10.
If g(ω) is a rational function, the assumptions are satisﬁed whenq (λ) is bounded on the upper half plane. In this case, the number of resonance poles
is ﬁnite and thus Eq.(3.31) becomes ﬁnite sum. For example if g(ω) = 1/(π(1 + ω 2 )) is√
the Lorentzian distribution, a resonance pole is given by λ = −1Kf1 − 1 (a root of√
Eq.(3.28)). Therefore Z10 (t) decays with the exponential rates Re( −1Kf1 − 1).
∗

Example 3.11.
If g(ω) is the Gaussian distribution, the assumptions are sat∗
isﬁed when q (λ) is of exponential type; that is, there exist positive constants C andβ such that |q ∗ (λ)| ≤ Ceβ|λ| . Since the analytic continuation g ∗ (λ) has an essential
singularity at inﬁnity, there exist inﬁnitely many resonance poles and they accumulate
at inﬁnity.
§ 3.5.

Semi-group generated by the operator T1 (Im(f1 ) ≥ 0)

In Sec.3.1 and Sec.3.4, we investigate the semi-group generated by the operator√√
T1 = −1M + −1Kf1 P for the cases f1 = 0 and Im(f1 ) < 0, respectively. In thissubsection, we consider the case Im(f1 ) ≥ 0.
Theorem 3.12. Suppose that the assumptions (ii) to (v) of Thm.3.9 hold. If Im(f1 ) ≥0, for an arbitrarily ﬁxed K > 0, the order parameter Z10 (t) = (eT1 t q, P0 ) decaysLinear stability of the incoherent solution for the Kuramoto-Daido model

19

exponentially as t → ∞.

√We show an idea of the proof. If Im(f1 ) = 0, T1 / −1 = M + KRe(f1 )P is self-adjoint√
and a rank one perturbation of the multiplication M. By Theorem X-4.3 in [8], T1 / −1√
and M are unitarily equivalent. Since (e −1Mt q, P0 ) decays exponentially (see Sec.3.1),
we can prove that so is (eT1 t q, P0 ).
If Im(f1 ) > 0, change the parameter as K → −K. Then, the problem is reduced
to the case K < 0 and Im(f1 ) < 0, and Thm.3.12 is proved in a similar manner to theproof of Thm.3.9.

References[1] L. V. Ahlfors, Complex analysis. An introduction to the theory of analytic functions of
one complex variable, McGraw-Hill Book Co., New York, 1978.
[2] N. J. Balmforth, R. Sassi, A shocking display of synchrony, Phys. D, 143 (2000), 21–55.
[3] H. Chiba, D. Pazó, Stability of an [N/2]-dimensional invariant torus in the Kuramotomodel at small coupling, Phys. D, 238 (2009), 1068–1081.
[4] H. Daido, Order Function and Macroscopic Mutual Entrainment in Uniformly Coupled
Limit-Cycle Oscillators, Prog. Theor. Phys., 88 (1992), 1213–1218.
[5] H. Daido, Critical Conditions of Macroscopic Mutual Entrainment in Uniformly Coupled
Limit-Cycle Oscillators, Prog. Theor. Phys., 89 (1993), 929–934.
[6] H. Daido, Generic scaling at the onset of macroscopic mutual entrainment in limit-cycle
oscillators with uniform all-to-all coupling, Phys. Rev. Lett., 73 (1994), 760.
[7] H. Daido, Onset of cooperative entrainment in limit-cycle oscillators with uniform all-to-allinteractions: bifurcation of the order function, Phys. D, 91 (1996), 24–66.
[8] T. Kato, Perturbation theory for linear operators, Springer-Verlag, Berlin, 1995.
[9] Y. Kuramoto, Self-entrainment of a population of coupled non-linear oscillators, International Symposium on Mathematical Problems in Theoretical Physics, Lecture Notes inPhys., 39. Springer, Berlin, 1975.
[10] Y. Kuramoto, Chemical oscillations, waves, and turbulence, Springer Series in Synergetics,
19. Springer-Verlag, Berlin, 1984.
[11] Y. Maistrenko, O. Popovych, O. Burylko, P. A. Tass, Mechanism of desynchronization in
the ﬁnite-dimensional Kuramoto model, Phys. Rev. Lett., 93 (2004) 084102.
[12] Y. L. Maistrenko, O. V. Popovych, P. A. Tass, Chaotic attractor in the Kuramoto model,
Int. J. of Bif. and Chaos, 15 (2005) 3457–3466.
[13] E. A. Martens, E. Barreto, S. H. Strogatz, E. Ott ,P. So, T. M. Antonsen, Exact results
for the Kuramoto model with a bimodal frequency distribution, Phys. Rev. E, 79 (2009),026204.
[14] R. E. Mirollo, S. H. Strogatz, The spectrum of the locked state for the Kuramoto modelof coupled oscillators, Phys. D, 205 (2005), 249–266.
[15] R. Mirollo, S. H. Strogatz, The spectrum of the partially locked state for the Kuramotomodel, J. Nonlinear Sci., 17 (2007), 309–347.
[16] A. Pikovsky, M. Rosenblum, J. Kurths, Synchronization: A Universal Concept in Nonlinear Sciences, Cambridge University Press, Cambridge, 2001.

20

Hayato Chiba[17] M. Reed, B. Simon, Methods of modern mathematical physics IV. Analysis of operators,
Academic Press, New York-London, 1978.
[18] S. H. Strogatz, From Kuramoto to Crawford: exploring the onset of synchronization in
populations of coupled oscillators, Phys. D, 143 (2000), 1–20.
[19] S. H. Strogatz, R. E. Mirollo, P. C. Matthews, Coupled nonlinear oscillators below thesynchronization threshold: relaxation by generalized Landau damping, Phys. Rev. Lett.,
68 (1992), 2730–2733.
[20] S. H. Strogatz, R. E. Mirollo, Stability of incoherence in a population of coupled oscillators,
J. Statist. Phys., 63 (1991), 613–635.
[21] N. J. Vilenkin, Special functions and the theory of group representations, American Mathematical Society, 1968.
[22] K. Yosida, Functional analysis, Springer-Verlag, Berlin, 1995.