多重穩定性對生物擺動器的延遲、耦合強度和凹凸性之關係

國 立 交 通 大 學

應 用 數 學 系

碩

士

論

文

多重穩定性對生物擺動器的延遲、耦合強度

和凹凸性之關係

Delay, Coupling Strengths, Concavity and

Multistable Synchronization

of Biological Oscillators

研 究 生：陳怡樺

指導教授：莊重 教授

多重穩定性對生物擺動器的延遲、耦合強度

和凹凸性之關係

Delay, Coupling Strengths, Concavity and Multistable

Synchronization of Biological Oscillators

研 究 生：陳怡樺 Student：I-Hua Chen

指導教授：莊重 Advisor：Jonq Juang

國 立 交 通 大 學

應 用 數 學 系

碩 士 論 文

A Thesis

Submitted to Department of Applied Mathematics

College of Science

National Chiao Tung University

in partial Fulfillment of the Requirements

for the Degree of

Master

in

i

多重穩定性對生物擺動器的延遲、耦合強度

和凹凸性之關係

研究生:陳怡樺 指導老師:莊重 教授

國 立 交 通 大 學

應 用 數 學 系 碩 士 班

摘 要

我們研究了對於兩個 Mirollo-Strogatz-type 的擺動器系統中耦合強度係

數、延遲時間和其函數凹凸性之間的互相影響。對於延遲函數是凹向上的

擺動器，我們將可以得到一個關於激進耦合與壓抑耦合的完整的相圖。特

別地，我們證明出加了延遲時間與激進耦合在一個擺動器上，將會產生了

一個多重穩定的同步化。在另外一方面，對於加了延遲時間的壓抑耦合會

產生同相同步化的情形是不穩定的。

關鍵詞：生物擺動器、延遲時間、函數凹凸性、壓抑耦合、激進耦合。

中 華 民 國 九 十 七 年 六 月

Delay, Coupling Strengths, Concavity and Multistable

Synchronization of Biological Oscillators

Student：I-Hua Chen Advisor：Jonq Juang

Department of Applied Mathematics

National Chiao Tung University

ABSTRACT

We study the interaction between coupling strengths, delay and concavity for

a pair of two Mirollo-Strogatz-type oscillators. For a pair of two concave

oscillators with a nonzero delay, the complete phase diagrams with respect to

both inhibitory and excitatory coupling are given. In particular, we prove that the

delay and excitatory coupling induce multistable synchronization for such

system. On the other hand, the in-phase synchronization is shown to be unstable

for inhibitory coupling if the delay is not zero.

Key words: Biological Oscillators, Delay, Concavity, and Inhibitory and

Excitatory Coupling.

iii

誌 謝

光陰似箭，兩年轉眼間就過了。遙想兩年前的我，懷著既期待又害怕的心情入學。

非常興奮能考上交通大學，並且期待能跟大家一起學習、相處；但同時也害怕同學高手

如雲，自己沒能表現好。

這兩年上過很多老師的課。實變的吳培元老師、常微分方程的李明佳老師、動態系

統理論的林松山老師、高等機率的吳慶堂老師、應用數學方法的賴明治老師、動態系統

導論的莊重老師、生物資訊的許元春老師，每位都是認真教學、且富有知識與理念的老

師。十分感謝以上老師們給予我的指導，上老師的課對我而言都是寶貴的學習經驗。

當然，最感謝的就是我的指導老師莊重，除了對我們在課業上的教導外，在日常生

活中，當我們遇到了困難，他十分地關心我們，並且給予我們一些建議與幫助，是一位

不可多得的好老師。

另外，非常感謝我的博班學長梁育豪，以及學姊張郁泉在這兩年來對我們盡心的照

顧。當我們在寫論文遇到問題時，他們總是不吝於給予我們幫助與教導，能寫出這篇論

文，學長姐是在我背後的推手、大功臣。感謝之心，不足以用言語形容。也感謝我同老

闆的同學千茹，是跟我共同成長的好伙伴，有他的陪伴，讓我在研究的路上不孤單。還

有要感謝我的同班同學們，不只在課業上給我的幫助與扶持，更帶給我歡笑與樂趣，為

我的生活添加更多的色彩。

在這邊，要特別感謝我的父母，在我的背後鼓勵我、支持我。當我遇到困難時，聽

我訴苦、給我打氣；當我獲得小小的成就時，能分享我的喜悅，他們是我能繼續有動力

走下去的最大來源。再次感謝他們對我的支持與用心與信任，讓我能順利完成學業。

最後，我要感謝我的男朋友鈺傑在這兩年來盡心的照顧以及陪伴，感謝有你在我身

邊，也期待能繼續一起共同扶持走下去。感謝在交通大學認識的各位同學及學長姊、學

弟妹們，能與你們做朋友真的很高興，希望在畢業後能繼續和各位保持聯絡，也祝福自

己與各位有很美好的將來。

CONTENTS

Abstract (in Chinese)………i

Abstract (in English) ………ii

Acknowledgement……….……….iii

Contents……….……….iv

1. Introduction…………..………...1

2. Model………..………..3

3. Mathematical Analysis.………..5

4. Excitatory Couplings,

ε

> …...………7

0

4.1. Construction of the Firemap…...……….7

4.2. Dynamics in Cases……….………..10

4.3. Construction of Phase Diagrams….………..15

5. Inhibitory Couplings,

ε

< ……….………18

0

5.1. Construction of the Firemap……….………18

5.2. Dynamics in Cases……….………21

1. Introduction

Large assemblies of oscillator units can spontaneously evolve to a state of large scale organi-zation. Synchronization is the best known phenomenon of this kind, where after some transient regime a coherent oscillatory activity of the set of oscillators emerges. This interesting phenom-enon is quite common in many different disciplines such as engineering [28], physics [4, 13] and [24], chemistry [14], as well as biology [27]. For example, southeastern fireflies, where thousands of individuals gathered on trees flash in unison. Other examples of biological oscillators are the rhythmic activity of cells of the heart pacemaker [12, 18, 20] and [26], of cells of pancreas [22] and [23], and of neural networks [2, 8, 20, 21] and [25]. In the recent years, this topic has gained increasing attention as synchronous oscillations have been observed in the visual cortex [11, 6, 5], which were related to Gestalt properties of the stimulus. It has been pointed out that synchronous firing activity may be a part of higher brain functions and a method for integrating distributed information in an abstract representation [16, 17]. Abstracting from biophysical de-tails, neurons belong to an important class of oscillators characterized by a pulselike interaction, i.e., where the coupling consists in the transmission of a short pulse from an oscillator to its partners. For an understanding of the general principles underlying synchronization phenom-ena, it is useful to consider abstract oscillator models which contain various existing models under very general assumptions and can be treated conveniently. We begin with describing the Peskin’s model of n integrate-and-fire oscillators. Let the state of the i-th oscillator be denoted by xi, where xi is subject to the dynamics dx_dti = −rixi+ si, 0 ≤ xi ≤ 1, i = 1, 2, · · · , n with

input si> 0, a normalized threshold 1 and leakiness ri≥ 0. When xi = 1, the ith oscillator fires

and xi jumps back to zero. As a consequence of the firing of ith oscillator, the activation of any

other oscillator j is incremented by the coupling ǫj,i. Should no confusion arise, we write ǫj,i as

ǫji. This model was later generalized by Mirollo and Strogatz [19]. It was assumed that the state

variable xi evolves according to a map fi. When xi reaches the threshold, the oscillator fires and

xi jumps back instantly to zero, and the activation of any other oscillator j is incremented by

the positive coupling ǫij. Specifically, xi evolve according to xi = fi(φi), where fi: [0, 1] → [0, 1]

is smooth, and strictly increasing, i.e., f_i′ > 0 on (0, 1). Here φi is a phase variable so that (i)

dφi

dt = 1 Ti

, where Ti is the cycle period for oscillator xi when evolving freely, (ii) φi = 0 when

the oscillator is at its lowest state xi= 0, and (iii) φi≡ 1 at the end of cycle when the oscillator

reaches the threshold xi = 1. Therefore, fi satisfy fi(0) = 0, fi(1) = 1. These maps fi are to

be called evolution maps. The inverses of fi are to be denoted by gi. If fi ≡ f , gi ≡ g, Ti ≡ T

and ǫij ≡ ǫ for all i, j, then the corresponding system is called identical. Otherwise, it is called

nonidentical. Moreover, if f_i′′> 0(resp., f_i′′ < 0) for all i, then the system is said to consist of the concave(resp., convex) oscillators. Assuming the identical system of convex oscillators, Mirollo and Strogatz [19] proved rigorously that the globally pulse-coupled oscillators always synchronize with zero phase difference. Their results solved the first conjecture of Peskin. Recently, Chang and Juang [3] prove that for ”nearly” identical system of convex oscillators, it will fire in unison. That solves the Peskin’s second conjecture. Those results pertain only to excitatory couplings, in realistic applications, however, one is also confronted with inhibitory couplings, which are

abundant, e.g., in the central nervous systems and whose importance for synchronization was pointed out recently [15]. Furthermore, in most applications the transmission of a pulse requires a finite propagation time. It is an important question, how synchronization over long distances can be achieved when such temporal delays prevail (like in the visual cortex [7]). That is to say every biological system has to deal with substantial delays that seem, heuristically speaking, to constrain the process of synchronization.

In [9, 10], the study of general mechanisms of synchronization of the system of convex oscilla-tors in cases where delay and also inhibitory couplings are present was given. In particular, they presented a complete mathematical analysis for pairs of two Mirollo-Strogatz-type oscillators for a wide range of delays τ and coupling strengths ǫ. Specifically, they showed that for inhibitory couplings, the presence of delays can lead to stable in-phase synchronization. For excitatory couplings, they showed that no stable in-phase synchronization exists. However, it was shown numerically in [1] that globally coupled oscillators with pulse interaction can synchronize under broader conditions. In particular, they demonstrated even the nonidentical system of concave oscillators can synchronize provided that the concavity of the system is not too large. Such numerical observation is also recently proved in [3].

The purpose of this thesis is to study how the roles of the concavity of the oscillators, the presence of delays and excitatory / inhibitory couplings play out in reaching or not reaching synchronization, and then compare the results with those from the system consisting of convex oscillators [9, 10].

We study the interaction between coupling strengths, delay and concavity for a pair of two Mirollo-Strogatz-type oscillators. For a pair of two concave oscillators with a nonzero delay, the complete phase diagrams with respect to both inhibitory and excitatory coupling are given. In particular, we prove that the delay and excitatory coupling induce multistable synchronization for such system. On the other hand, the in-phase synchronization is shown to be unstable for inhibitory coupling if the delay is not zero.

Summing up the earlier results and ours here, we conclude the following. For system of two identical concave oscillators, it will not reach synchrony without delay. With delay, the system will acquire in-phase synchronization with excitatory coupling on.

2. Model

The network consists of N relaxation oscillators, which are caricatures of real pulse-coupled neurons in biological systems [19]. Each oscillator i may be described by a smooth function f (φi), which is concave up and monotonically increasing [f′ > 0, f′′ > 0, f (0) = 0, f(1)=1]. f

plays the role of an amplitude (e.g. the membrane potential) and φi ∈ [0, 1] is a phase, which

in the case of vanishing input from other oscillators corresponds to the normalized time elapsed since the last firing of i. We assume from here on that the speed of each oscillator is one. When f reaches the threshold fs := 1, the oscillator fires and φi and f are reset to zero. After a

time delay t = τ, τ ∈ (0, 0.5) , the spike reaches all the other oscillators (no self-interaction) and raises (excitatory couplings) or lowers (inhibitory couplings) their amplitudes by an amount ǫ = ǫ(N − 1)−1, where ǫ denotes the normalized coupling strength (ǫ ∈ (0, 1]). The coupling to the oscillators j may be represented equivalently by an increase or decrease in phase ∆φj

φj+ ∆φj = g(min[f (φj) + ǫ, 1]) := F+(φj, ǫ), where g = f−1 (2.1)

φj+ ∆φj = g(max[f (φj) + ǫ, 0]) := F−(φj, ǫ), where g = f−1 (2.2)

where Eqs. (2.1) and (2.2) refer to excitatory and inhibitory coupling, respectively. We point out that the concavity of f is responsible for the dependence of ∆φj on φj, the larger the phase

φj, the smaller the phase shift ∆φj.

1 f D 2 f D e e e 3 f D φ f (φ)

Figure 2.1. Function f (φ) and the dependence of the phase shift on φ. With excitatory couplings, an increase of ǫ in the amplitude f corresponds to a shift of phase that is larger when starting with a smaller phase (∆φ2 < ∆φ1). A negative

phase shift ∆φ3 occurs with inhibitory couplings.

In this work, we consider the identical system of two concave oscillators with the delay. Before we treat a pair of these oscillators in a mathematical analysis, we note some simple properties of the functions F− and F+ introduced in Eqs. (2.1) and (2.2) that we will need in the next

paragraph:

A1: F+(φ, ǫ) > φ, for φ < 1.

A2: F−(φ, ǫ) < φ, for φ > 0.

A3: F+(c + φ, ǫ) − F+(c − φ, ǫ) < (c + φ) − (c − φ) = 2φ, for c − φ > 0 and 3

F+(c + φ, ǫ) < 1.

A4: F−(c + φ, ǫ) − F−(c − φ, ǫ) > (c + φ) − (c − φ) = 2φ, for c + φ < 1 and

F−(c − φ, ǫ) > 0. A5: f (φ2) − f (φ1) < f (φ2+ a) − f (φ1+ a), if φ1 < φ2, a > 0, f′> 0, f′′> 0. A6: f (φ′1+ a′) − f (φ′1− a′) < f (φ′2+ a′) − f (φ′2− a′), for φ′1 = φ1+ a 2, φ ′ 2 = φ2+ a 2, a′₌ a 2 and φ1 < φ2. A7: F+(φ2, ǫ) − φ2< F+(φ1, ǫ) − φ1, if 0 < φ1< φ2, f′ > 0, f′′> 0, F (φ2, ǫ) < 1.

3. Mathematical Analysis

In this section, we will derive phase diagrams that allow one to determine if and how two oscillators synchronize their activities. We hereby consider a system S of two oscillators A and B, both either inhibitorily or excitatorily coupled together with time delay τ .

To study the dynamics of S, we begin with assuming that the oscillator A just reaches the threshold and is reset to φA = 0 and φ = φB > 0. We further assume that if φB < τ , then

φB must have fired φB time earlier. This assumption is not necessary in a mathematical sense,

but makes the analysis easier by reducing the number of case distinctions. It should also be noted that since the speed of the oscillators is assumed to be one, the phase φ and the time φ is interchangeable. As the system S evolves, the phase positions of φA and φB at time t are to

be denoted by φA(t) and φB(t), respectively. We next define a firemap and a return map for

the system S of two oscillators A and B. Let t = tp,i denote the time when oscillator i has just

fired its pth time and its phase is reset to zero. In this situation the system S is said to reach a ground state. The firemap and the return map are defined as follows, similarly as in [10]:

(i) Firemap h(φj(tp,i)) = φl(tq,k), where j 6= i, l 6= k and

tq,k =      min r∈N,k∈{A,B}{tr,k|tr,k > tp,i+ τ } if φ ∈ I1, min

r∈N,k∈{A,B}{tr,k|tr,k > tp,i} otherwise .

(ii)Return map R(φj(tp,i)) = φj(tp+1,i) for i 6= j.

Suppose the system S just reaches a ground state in oscillator i, then the firemap takes the phase position of the non-ground state oscillator into the phase position of the non-ground state oscillator when the system reaches the immediate next ground state. The return map takes the phase position of the non-ground state oscillator into the phase position of the non-ground state oscillator when the system reaches the next ground state in the same oscillators.

Before we begin the analysis of the different cases or configurations of the dynamics, we want to illustrate our motivation for our choice of intervals in the subspace of initial phase differences φ. Let us consider that S is a GS with oscillator A just being reset to φA= 0 such that φ = φB.

In a first interval I1, both oscillators have fired, but their spikes did not reach their destination

yet. Therefore, the consequences of the two pulses being received have to be evaluated. In a second interval I2, only the spike of oscillator A did not reach B and has to be taken into

account. In a third interval I3, oscillator B will reach the threshold before the spike of A can

be received. These considerations lead to the following definitions for I1, I2, and I3:

I1 : φB ∈ [0, τ ),

I2 : φB ∈ [τ, 1 − τ ],

I3 : (Ex0, In0) : φB ∈ (1 − τ, 1].

On the one hand, the dynamics is very simple if we are looking at domain I3. After a time

t = 1 − φ (from here on, we identify t with the time elapsed since S has been in a GS), S will be in a GS with φB = 0, which leads to a firemap h(φ) = 1 − φ. Additional case distinctions

are not required here, and I3 will be referenced as region Ex0 (for excitatory couplings) or In0

(for inhibitory couplings) for reasons of conformity.

On the other hand, the detailed analysis of I1 and I2following in the next section requires one

to distinguish between excitatory and inhibitory coupling. In each case, we first heuristically describe the temporal development of S to motivate the mathematical notation of the dynamics that will follow afterwards. For brevity, we denote a spike originating from oscillator A as ”spike A” and the oscillator A simply as ”A”. Each domain will be partitioned into smaller regions, which we denote ExN or InN with successive numbering for excitatory and inhibitory coupling, respectively.

4. Excitatory Couplings, ǫ > 0 4.1. Construction of the FireMap.

Configuration I1

Ex1 : ǫ < 1 − f (τ + φ), ǫ < f (1 − φ) − f (τ − φ) and φ ∈ I1 (For fixed ǫ, φ ∈ [0, min{g(1 −

ǫ) − τ, u1b, τ }), where u1b satisfies ǫ = f (1 − φ) − f (τ − φ).) time t φA φB 0 0 φ τ − φ τ − φ → F+(τ − φ, ǫ) τ τ F+(τ − φ, ǫ) + φ < 1 τ + φ → F+(τ + φ, ǫ) < 1 h(φ) = 1 − |F+(τ + φ, ǫ) − F+(τ − φ, ǫ) − φ| (4.1)

Using relation A3, we know that 0 ≤ F+(τ + φ, ǫ) − F+(τ − φ, ǫ) ≤ 2φ for all φ ∈ Ex1, and

both the left and right equalities hold just as φ = 0. Then the firemap satisfies h(φ) = 1 − |F+(τ + φ, ǫ) − F+(τ − φ, ǫ) − φ|

≥ 1 − |2φ − φ|

= 1 − φ > 1 − τ. (4.2)

Moreover, h(φ) = 1 − φ just as φ = 0.

Ex2a : 1 − f (φ + τ ) ≤ ǫ < f (1 − φ) − f (τ − φ), and φ ∈ I1 (For fixed ǫ, φ ∈

[g(1 − ǫ) − τ, min{u1b, τ }).) time t φA φB 0 0 φ τ − φ τ − φ → F+(τ − φ, ǫ) τ τ F+(τ − φ, ǫ) + φ < 1 τ + φ → F+(τ + φ, ǫ) = 1 h(φ) = F+(τ − φ, ǫ) + φ (4.3)

Since F+(τ − φ, ǫ) > F+(τ + φ, ǫ) − 2φ ∀φ ∈ Ex2a, and F+(τ + φ, ǫ) = 1,

h(φ) = F+(τ − φ, ǫ) + φ

> F+(τ + φ, ǫ) − 2φ + φ

≥ 1 − 2φ + φ = 1 − φ

≥ 1 − τ. (4.4)

Ex2b : f (1 − φ) − f (τ − φ) ≤ ǫ < 1 − f (τ + φ), and φ ∈ I1 (For fixed ǫ, φ ∈ [u1b, min{g(1 − ǫ) − τ, τ }).) time t φA φB 0 0 φ τ − φ τ − φ → F+(τ − φ, ǫ) < 1 τ τ − φ + 1 − F+(τ − φ, ǫ) 1 → 0 τ + 1 − F+(τ − φ, ǫ) τ F+(τ − φ, ǫ) + φ − 1 ≥ 0 τ + φ → F+(τ + φ, ǫ) < 1 h(φ) = F+(τ − φ, ǫ) − F+(τ + φ, ǫ) + φ (4.5)

Since F+(τ − φ, ǫ) − F+(τ + φ, ǫ) < 0 ∀φ ∈ Ex2b, the firemap satisfies

h(φ) = F+(τ − φ, ǫ) − F+(τ + φ, ǫ) + φ < 0 + φ = φ. (4.6) Ex3 : f (1 − φ) − f (τ − φ) ≤ ǫ < 1 − f (τ − φ), 1 − f (τ + φ) ≤ ǫ < 1 − f (τ − φ) and φ ∈ I1. time t φA φB 0 0 φ τ − φ τ − φ → F+(τ − φ, ǫ) < 1 τ τ − φ + 1 − F+(τ − φ, ǫ) 1 → 0 τ + 1 − F+(τ − φ, ǫ) τ φ − 1 + F+(τ − φ, ǫ) ≥ 0 τ + φ → F+(τ + φ, ǫ) = 1 h(φ) = φ − 1 + F+(τ − φ, ǫ) (4.7)

Since F+(τ − φ, ǫ) < 1, the firemap satisfies

h(φ) = φ − 1 + F+(τ − φ, ǫ)

< φ − 1 + 1 = φ. (4.8)

Ex4 : ǫ ≥ 1 − f (τ − φ) and φ ∈ I1 (For fixed ǫ, φ ∈ [0, τ − g(1 − ǫ)).)

time t φA φB

0 0 φ

τ − φ τ − φ → F+(τ − φ, ǫ) = 1 τ

τ φ τ + φ → F+(τ + φ, ǫ) = 1

h(φ) = φ (4.9)

Configuration I2

Ex5 : ǫ < 1 − f (τ + φ) and φ ∈ I2 (For fixed ǫ, φ ∈ [τ, g(1 − ǫ) − τ ).)

time t φA φB

0 0 φ

τ τ τ + φ → F+(τ + φ, ǫ) < 1

h(φ) = 1 − F+(τ + φ, ǫ) + τ (4.10)

Since τ + φ < F+(τ + φ, ǫ) < 1 and τ ≤ φ, the firemap satisfies

h(φ) = 1 − F+(τ + φ, ǫ) + τ

< 1 − τ − φ + φ = 1 − τ, (4.11)

and

h(φ) = 1 − F+(τ + φ, ǫ) + τ

> 1 − 1 + τ = τ. (4.12)

Ex6 : ǫ ≥ 1 − f (φ + τ ) and φ ∈ I2 (For fixed ǫ, φ ∈ [max{τ, g(1 − ǫ) − τ }, 1 − τ ].)

time t φA φB 0 0 φ τ τ τ + φ → F+(τ + φ, ǫ) = 1 h(φ) = τ (4.13) Configuration I3 Ex0 : φ ∈ I3 time t φA φB 0 0 φ 1 − φ 1 − φ 1 → 0 h(φ) = 1 − φ (4.14) 9

4.2. Dynamics in Cases.

Configuration I1

Dynamics in Ex1 : ǫ < 1 − f (τ + φ), ǫ < f (1 − φ) − f (τ − φ) and φ ∈ I1 (For fixed ǫ,

φ ∈ [0, min{g(1 − ǫ) − τ, u1b, τ }), where u1b satisfies ǫ = f (1 − φ) − f (τ − φ).) h(φ) = 1 − |F+(τ + φ, ǫ) − F+(τ − φ, ǫ) − φ|

By Eq. (4.2), the firemap h : Ex1 7→ Ex0. Specifically, when φ satisfies F+(τ − φ, ǫ) + φ =

F+(τ + φ, ǫ), h(φ) = 1, which means φA and φB synchronize. Moreover, the return map R

satisfies

R(φ) = |F+(τ + φ, ǫ) − F+(τ − φ, ǫ) − φ| ≤ φ. (4.15)

The equality holds just as φ = 0. Since φ ∈ Ex1 and R(φ) < φ, R(φ) ∈ Ex1.

Proposition 4.1. For each φ ∈ Ex1, Rk_{(φ) → 0, as k → ∞.}

Proof. Suppose not, i.e., ∃φ ∈ Ex1, Rk_{(φ) 9 0, as k → ∞. Then R}k_{(φ) 6= 0 ∀k ∈ N. Otherwise,}

Rm(φ) = 0 for large enough m ∈ N. By Eq. (4.15), Rk(φ) ∈ Ex1 ∀k ∈ N and R(φ) < φ. Inductively, φ > R(φ) > · · · > Rk−1(φ) > Rk_{(φ) > 0. Since {R}k_{(φ)} is a monotone and}

bounded sequence in Ex1, and R(φ) is continuous in Ex1. {Rk_{(φ)} is a convergent sequence,}

i.e., Rk_{(φ) → φ}∗_{, for some φ}∗ _{∈ (0, τ ), and such φ}∗ _{is a fixed point. It is a contradiction since}

no fixed point in Ex1\{0}, by Eq. (4.15).



Corollary 4.1. For each φ ∈ Ex1, φA and φB synchronize.

Dynamics in Ex2a: 1 − f (τ + φ) ≤ ǫ < f (1 − φ) − f (τ − φ) and φ ∈ I1 (For fixed ǫ,

φ ∈ [g(1 − ǫ) − τ, min{u1b, τ }).)

h(φ) = F+(τ − φ, ǫ) + φ

Remark 4.1. φ = 0 6∈ Ex2a.

By Eq. (4.4), the firemap h : Ex2a7→ Ex0. Then the return map R satisfies

R(φ) = h(h(φ)) = 1 − h(φ) = 1 − F+(τ − φ, ǫ) − φ

Rk−1(φ) > Rk(φ) > 0. Thus, Rk(φ) → φ∗ for some φ∗ ∈ Ex2a, and such φ∗ is a fixed point. It

makes a contradiction with Eq. (4.16).



Corollary 4.2. For each φ ∈ Ex2a, there is a k ∈ N, depending on φ, such that Rk(φ) ∈ Ex1.

Dynamics in Ex2b : f (1 − φ) − f (τ − φ) ≤ ǫ < 1 − f (τ + φ) and φ ∈ I1 (For fixed ǫ,

φ ∈ [u1b, min{g(1 − ǫ) − τ, τ }).)

h(φ) = F+(τ − φ, ǫ) − F+(τ + φ, ǫ) + φ

Remark 4.2. φ = 0 6∈ Ex2b.

By Eq. (4.6), the firemap h : Ex2b 7→ Ex1 ∪ Ex2b.

Proposition 4.3. There is no φ ∈ Ex2b such that hk(φ) ∈ Ex2b for all k = 0, 1, 2, · · · .

Proof. Suppose not, i.e., ∃φ ∈ Ex2b such that hk(φ) ∈ Ex2b for all k = 0, 1, 2 · · · . Note that

h(φ) is continuous in Ex2b. Moreover, h(φ) < φ, ∀φ ∈ Ex2b. Inductively φ > R(φ) > · · · >

Rk−1(φ) > Rk(φ) > 0. Thus, hk(φ) → φ∗ for some φ∗ ∈ Ex2b, and such φ∗ is a fixed point. It

makes a contradiction with Eq. (4.6).



Corollary 4.3. For each φ ∈ Ex2b, there is a k ∈ N, depending on φ, such that hk(φ) ∈ Ex1.

Dynamics in Ex3 : 1 − f (τ + φ) ≤ ǫ < 1 − f (τ − φ), f (1 − φ) − f (τ − φ) ≤ ǫ < 1 − f (τ − φ) and φ ∈ I1.

h(φ) = φ − 1 + F+(τ − φ, ǫ)

Remark 4.3. φ = 0 6∈ Ex3.

By Eq.(4.8), the firemap h : Ex3 7→ I1. More precisely,

(i) If ǫ < 1 − f (τ ), h : Ex3 7→ Ex1 ∪ Ex2a∪ Ex2b∪ Ex3, since the slope of the lower bound of

Ex3 is negative (−f′(τ + φ) < 0 and −f′(1 − φ) + f′(τ − φ) < 0.

(ii) If ǫ ≥ 1 − f (τ ), h : Ex3 7→ Ex3 ∪ Ex4 , since the slope of the upper bound of Ex3 is positive (f′(τ − φ) > 0).

Proposition 4.4. (i) If ǫ < 1 − f (τ ), there is no φ ∈ Ex3 such that hk_{(φ) ∈ Ex3 for all}

k = 0, 1, 2, · · · . (ii) If ǫ ≥ 1 − f (τ ), for each φ ∈ Ex3, hk_{(φ) ∈ Ex3 for all k ∈ N. Moreover,}

hk_{(φ) → τ − g(1 − ǫ) for all φ ∈ Ex3.}

Proof.

(i) Suppose not, i.e., ∃φ ∈ Ex3 such that hk_{(φ) ∈ Ex3, for all k = 0, 1, 2, · · · . Since as}

ǫ < 1 − f (τ ), Ex3 = [g(1 − ǫ), τ ), or Ex3 = [u1b, τ ) and the sequence {h

k_{(φ)} must satisfy}

hk_{(φ) → φ}∗_{, where φ}∗ _{is a fixed point for h in Ex3. It makes a contradiction with Eq. (4.8).} 11

(ii) If ǫ ≥ 1 − f (τ ), Ex3 = (τ − g(1 − ǫ), τ ), where τ − g(1 − ǫ) ≥ 0, and the equality holds just as ǫ = 1 − f (τ ). Because of h(τ−) = τ − 1 + F+(τ − τ, ǫ) = τ − 1 + g(ǫ) < τ, (4.17) and h((τ − g(1 − ǫ))+) = τ − g(1 − ǫ) − 1 + F+(τ − τ + g(1 − ǫ), ǫ) = τ − g(1 − ǫ) − 1 + 1 = τ − g(1 − ǫ). (4.18) Also, for any φ ∈ Ex3,

h′(φ) = 1 + F+′(τ − φ, ǫ)

= 1 − g′(f (τ − φ) + ǫ)f′(τ − φ) > 1 − g′(f (τ − φ))f′(τ − φ)

= 1 − 1 = 0. (4.19)

We obtain τ − g(1 − ǫ) = h((τ − g(1 − ǫ))+

) < h(φ) < h(τ−) < τ . Therefore, for each φ ∈ Ex3, hk_{(φ) ∈ Ex3 for all k ∈ N. Moreover, by Eqs. (4.8), (4.17), (4.18), (4.19), the sequence {h}k_(φ)}

is decreasing. Hence, hk(φ) → φ∗. Since there is no fixed point in Ex3, φ∗∈ bd(Ex3). It follows that φ∗ = τ − g(1 − ǫ) which is also the boundary of Ex4.

 Specifically, if ǫ = 1 − f (τ ), we obtain φ∗ _{= 0, i.e., φ}

A and φB synchronize.

Corollary 4.4. (i) If ǫ < 1 − f (τ ), for each φ ∈ Ex3, hk_{(φ) ∈ Ex1 ∪ Ex2}

a∪ Ex2b for some

k ∈ N. (ii) If ǫ ≥ 1 − f (τ ), for each φ ∈ Ex3, hk(φ) → τ − g(1 − ǫ) as k → ∞. Specifically, if ǫ = 1 − f (τ ), for each φ ∈ Ex3, φA and φB synchronize.

Dynamics in Ex4 : ǫ ≥ 1 − f (τ − φ) and φ ∈ I1 (For fixed ǫ, φ ∈ [0, τ − g(1 − ǫ)).)

h(φ) = φ

Thus, h : Ex4 7→ Ex4. Specifically, if φ = 0 ∈ Ex4, hk_{(φ) = 0 ∈ Ex4, ∀k ∈ N.}

Configuration I2

Dynamics in Ex5 : ǫ < 1 − f (τ + φ) and φ ∈ I2 (For fixed ǫ, φ ∈ [τ, g(1 − ǫ) − τ ).)

h(φ) = 1 − F+(τ + φ, ǫ) + τ

By Eqs. (4.11), (4.12), we have τ < h(φ) < 1 − τ . Hence h : Ex5 7→ I2.

Lemma 4.1. For each φ ∈ Ex5, hk_{(φ) ∈ Ex5 for k = 0, 1, 2, · · · .}

Proof. For each φ ∈ Ex5, since ǫ < 1 − f (τ + φ) and φ > τ , ǫ < 1 − f (2τ ). Let φ0 = g(1 − ǫ),

∆φ0 = g(f (φ0)+ǫ)−ǫ and ∆2τ = g(f (2τ )+ǫ)−ǫ. Since g(f (φ0)+ǫ) = 1 and f′′ > 0, ∆φ0 < ∆2τ ,

i.e., F+(φ0, ǫ) − φ0 < F+(2τ, ǫ) − 2τ . Clearly, F+(τ + φ, ǫ) ≥ F+(2τ, ǫ) for φ ∈ [τ, 1 − τ ]. Then

we have

h(φ) ≤ h(τ ) = 1 − F+(2τ, ǫ) + τ

< 1 + φ0− F+(φ0, ǫ) − 2τ + τ

= 1 + g(1 − ǫ) − g(f (g(1 − ǫ)) + ǫ) − τ = g(1 − ǫ) − τ.

Thus, h(φ) ∈ (τ, g(1 − ǫ) − τ ), i.e., h(φ) ∈ Ex5. Inductively, hk_{(φ) ∈ Ex5 for k = 0, 1, 2, · · · .}

Therefore, the lemma holds.

 Lemma 4.2. −1 < h′(φ) < 0 and 0 < R′(φ) < 1 ∀φ ∈ Ex5.

Proof. From direct computation, 0 > h′(φ) = −g′(f (τ + φ) + ǫ) · f′(τ + φ) > −g′(f (τ + φ)) · f′_{(τ + φ) = −1, since g}′′ _{< 0 and g}′_{(f (τ + φ)) · f}′_{(τ + φ) = 1. Since R}′_{(φ) = h}′_{(h(φ)) · h}′_{(φ), we}

obtain 0 < R′(φ) < 1.

 Proposition 4.5. There exists a unique fixed point for h in Ex5, and it is an attractor. Proof. Define A(φ) = h(φ) − φ. Clearly, A′_{(φ) = h}′_{(φ) − 1. From Lemma 4.2, we have h}′_{(φ) < 0,}

and then A′_{(φ) = h}′_{(φ) − 1 < 0. By the definition of Ex5, Ex5 = [τ, g(1 − ǫ) − τ ). Denote}

δ ≡ g(1 − ǫ) − τ . Now, let’s check the values of A(δ−) and A(τ ).

(i) If φ = δ, we figure out that F+(δ + τ, ǫ) = 1. It follows that h(δ−) = τ , and

A(δ−) = h(δ−) − δ = τ − δ < 0. (ii) If φ = τ , we have

A(τ ) = h(τ ) − τ

= 1 − F+(2τ, ǫ) + τ − τ

= 1 − F+(2τ, ǫ) > 0.

The last inequality follows directly from the definition of Ex5.

Hence h has a unique fixed point φ∗ ∈ Ex5, i.e., h(φ∗_{) = φ}∗_{. Since R(φ}∗_{) = h}2

(φ∗) and 0 < R′(φ) < 1, such φ∗ is also the unique fixed point for R in Ex5. Moreover, since

φ∗ < R(φ) < φ if φ > φ∗,

and

φ∗ > R(φ) > φ if φ < φ∗, the fixed point φ∗ for R and hence for h is an attractor.

 Corollary 4.6. For each φ ∈ Ex5, hk_{(φ) converges to the fixed point φ}∗ _{∈ Ex5, which is an}

attractor.

Dynamics in Ex6 : ǫ ≥ 1 − f (τ + φ) and φ ∈ I2 (For fixed ǫ, φ ∈ [max{τ, g(1 − ǫ) −

τ }, 1 − τ ].)

h(φ) = τ

Thus, h : Ex6 7→ I2.

Corollary 4.7. For each φ ∈ Ex6, the firemap maps φ to the line φ = τ .

Configuration I3

Dynamics in Ex0 : φ ∈ I3

h(φ) = 1 − φ

Thus, h : Ex0 7→ I1.

4.3. Construction of Phase Diagrams.

Peskin modeled the pacemaker as a network of N ”integrate-and-fire” oscillators, each char-acterized by a voltagelike state variable xi, subject to the dynamics

dxi

dt = −rxi+ s, 0 ≤ xi≤ 1, i = 1, 2, · · · N. (4.20)

Mirollo and Strogatz [19] proposed a new model and they assume that x evolves according to x = f (φ), where f : [0, 1] → [0, 1] is smooth, monotonic increasing (f′ > 0), and concave down (f′′_{< 0). Here φ ∈ [0, 1] is a phase variable such that (i) dφ/dt = 1/T , where T is a cycle}

period. (ii) f (0) = 0. (iii) f (1) = 1. Then they generate the function f and its inverse function g which are given in the following, respectively.

f (φ) = s r − s r( s − r s ) φ (4.21) g(x) = ln( s−rx s ) ln(s−r_s )

In contrast, we assume there exists a function f which is smooth, monotonic increasing (f′ _{> 0), and concave up (f}′′_{> 0). So, we switch the function f and g of Eq. (4.21) to be the}

model of our phase diagrams.

We fixed s = 1, and change the variable r, then we have three different kinds of figures in Figure 4.1:

(a) When r = 0.4 and τ = 0.2, f (1 − φ) − f (τ − φ) > 1 − f (τ + φ) ∀φ ∈ I1.

(b) When r = 0.7 and τ = 0.2, f (1 − φ) − f (τ − φ) < 1 − f (τ + φ) ∀φ ∈ I1.

(c) When r = 0.72 and τ = 0.35, f (1 − φ) − f (τ − φ) and 1 − f (τ + φ) have an intersection in I1.

Moreover, we add Figure 4.2 which is a phase diagram we make from Ernst [10] as a contrast. And it’s energy function f and inverse function g is exactly the Mirollo-Strogatz-type (Eq. 4.21) oscillators we just mentioned.

And there are some notations for the following figures:

l1 : ǫ = 1 − f (τ + φ)

l2 : ǫ = f (1 − φ) − f (τ − φ)

l3 : ǫ = 1 − f (τ − φ)

In the following tables, we fixed an arbitrary ǫ > 0, and shift φ in [0,1]. Therefore, we clearly see which area (ExN ) does the oscillator begin and figure out where it ends.

(I) f (φ) = ln( s−rφ s ) ln(s−r_s ) , f ′ _{> 0, f}′′_{> 0}

Ex1

Ex2

a

Ex3

Ex4

Ex5

Ex0

Ex6

l1

l

2

l3

1-f( )τ f(1- )τ 1-f(2 )τ

Ex0

₄

Ex0

₃ (a) Figure 4.1. (a) s = 1, r = 0.4, τ = 0.2. (b) s = 1, r = 0.7, τ = 0.2. (c) s = 1, r = 0.72, τ = 0.35. Ex1 1-f( )τ f(1- )τ 1-f(2 )τ l1 l2 l3 Ex3 Ex6 Ex5 Ex0 Ex2 b Ex4 (b) Ex4 Ex2_b Ex2_a Ex5 Ex0 Ex6 Ex1 Ex3 1-f( )τ f(1- )τ 1-f(2 )τ f e l1 l3 l2 (c) ǫ Dynamics Multistability 0 < ǫ < 1 − f (2τ ) φ → 0, if φ ∈ I1, I3 Complete Syn. φ → φ∗_Ex5, if φ ∈ I2 Lag Syn. with Lag φ∗_Ex5

(II) f (φ) = s r − s r( s − r s ) φ_{, f}′ _{> 0, f}′′ _{< 0}

Ex4

l3 l2 l1 1-f( )τ f(1- )τ 1-f(2 )τ

Ex1

Ex3

Ex2

Ex1

Ex1

Ex5

Ex6

Ex0

Figure 4.2. _{s := 1, r := 0.4, τ = 0.2, where s > r > 0} ǫ Dynamics Multistability 0 < ǫ < 1 − f (2τ ) φ → φ ∗

Ex5∗, if φ ∈ I1, Ex5∗, I3 Lag Syn. with Lag φ∗Ex5∗

φ →Period(τ, h(τ )), if φ ∈ I2\Ex5∗ Lag Syn. with Lag τ

1 − f (2τ ) ≤ ǫ < 1 − f (τ ) φ → φ

∗

Ex4 or τ , if φ ∈ I1, I3 Depends on the behavior of Ex2

φ → τ , if φ ∈ I2 Lag Syn. with Lag τ

1 − f (τ ) ≤ ǫ ≤ 1 φ → φ

∗

Ex4, if φ ∈ I1, I3 Lag Syn. with Lag φ∗Ex4

φ → τ , if φ ∈ I2 Lag Syn. with Lag τ

5. Inhibitory Couplings, ǫ < 0 5.1. Construction of the FireMap.

Configuration I1

In1 : |ǫ| < f (τ − φ) and φ ∈ I1

⇒ ǫ > −f (τ − φ) and φ ∈ I1 (For fixed ǫ, φ ∈ [0, τ − g(−ǫ)).)

time t φA φB

0 0 φ

τ − φ τ − φ → F−(τ − φ, ǫ) > 0 τ

τ F−(τ − φ, ǫ) + φ τ + φ → F−(τ + φ, ǫ) > 0

h(φ) = 1 − F−(τ + φ, ǫ) + F−(τ − φ, ǫ) + φ (5.1)

Using relation A4, F−(τ + φ, ǫ) ≥ F−(τ − φ, ǫ) + 2φ, we have φB(τ ) ≥ φA(τ ). Since F−(τ −

φ, ǫ) > 0, the firemap satisfies

h(φ) = 1 − F−(τ + φ, ǫ) + F−(τ − φ, ǫ) + φ

> 1 − (τ + φ) + F−(τ − φ, ǫ) + φ

= 1 − τ + F−(τ − φ, ǫ)

> 1 − τ. (5.2)

In2 : f (τ − φ) ≤ |ǫ| < f (τ + φ) and φ ∈ I1

⇒ −f (τ − φ) ≥ ǫ > −f (τ + φ) and φ ∈ I1 (For fixed ǫ, φ ∈ [0, τ ).)

time t φA φB

0 0 φ

τ − φ τ − φ → F−(τ − φ, ǫ) = 0 τ

τ φ τ + φ → F−(τ + φ, ǫ) > 0

h(φ) = 1 − |F−(τ + φ, ǫ) − φ| (5.3)

Since F−(τ + φ, ǫ) < τ + φ, the firemap satisfies

h(φ) = 1 − |F−(τ + φ, ǫ) − φ|

> 1 − |τ + φ − φ|

= 1 − τ. (5.4)

(ii) If F−(τ + φ, ǫ) < φ,

Define In2c : f (φ) − f (τ + φ) > ǫ > −f (τ + φ) (For fixed ǫ, φ ∈ [g(−ǫ) − τ, min{l2b, τ }).) R(φ) = h(φ) = 1 + F−(τ + φ, ǫ) − φ

(iii) If F−(τ + φ, ǫ) = φ,

Define In2d: ǫ = f (φ) − f (τ + φ) (For fixed ǫ, φ = l2b.) h(φ) = 1

In3 : f (τ + φ) ≤ |ǫ| and φ ∈ I1

⇒ ǫ ≥ −f (τ + φ) and φ ∈ I1 (For fixed ǫ, φ ∈ [0, min{g(−ǫ) − τ, τ }).)

time t φA φB 0 0 φ τ − φ τ − φ → F−(τ − φ, ǫ) = 0 τ τ φ τ + φ → F−(τ + φ, ǫ) = 0 h(φ) = R(φ) = 1 − φ (5.5) Configuration I2 In4 : |ǫ| < f (τ + φ) − f (2τ ) and φ ∈ I2

⇒ ǫ > f (2τ ) − f (τ + φ) and φ ∈ I2 (For fixed ǫ, φ ∈ (l4, 1 − τ ], where l4 = g(f (2τ ) − ǫ) − τ .)

time t φA φB

0 0 φ

τ τ τ + φ → F−(τ + φ, ǫ) > 0

h(φ) = 1 − F−(τ + φ, ǫ) + τ (5.6)

We assume that |ǫ| is small enough such that F−(τ +φ, ǫ) > 2τ . Since 2τ < F−(τ +φ, ǫ) < τ +φ,

the firemap satisfies

h(φ) = 1 − F−(τ + φ, ǫ) + τ < 1 − 2τ + τ = 1 − τ, (5.7) and h(φ) = 1 − F−(τ + φ, ǫ) + τ > 1 − τ − φ + τ = 1 − φ ≥ τ. (5.8) 19

In5 : f (τ + φ) − f (2τ ) ≤|ǫ| < f (τ + φ) and φ ∈ I2 ⇒ f (2τ ) − f (τ + φ) ≥ ǫ > −f (τ + φ) and φ ∈ I2 time t φA φB 0 0 φ τ τ τ + φ → F−(τ + φ, ǫ) > 0 h(φ) = 1 −|F−(τ + φ, ǫ) − τ | (5.9)

We assume that |ǫ| is bigger than that in In4 such that 0 < F−(τ + φ, ǫ) ≤ 2τ , then the

firemap satisfies

h(φ) = 1 − |F−(τ + φ, ǫ) − τ |

≥ 1 − |2τ − τ |

= 1 − τ. (5.10)

In6 : f (τ + φ) ≤ |ǫ| and φ ∈ I2

⇒ −f (τ + φ) ≥ ǫ and φ ∈ I2 (For fixed ǫ, φ ∈ [τ, g(−ǫ) − τ ].)

time t φA φB 0 0 φ τ τ τ + φ → F−(τ + φ, ǫ) = 0 h(φ) = R(φ) = 1 − τ (5.11) Configuration I3 In0 : φ ∈ I3 time t φA φB 0 0 φ 1 − φ 1 − φ 1 → 0 h(φ) = 1 − φ (5.12)

5.2. Dynamics in Cases.

Configuration I1

Dynamics in In1 : ǫ > −f (τ − φ) and φ ∈ I1 (For fixed ǫ, φ ∈ [0, τ − g(−ǫ)).)

h(φ) = 1 − F−(τ + φ, ǫ) + F−(τ − φ, ǫ) + φ

By Eq. (5.2), thus h : In1 7→ In0. Then the return map satisfies R(φ) = h(h(φ)) = 1 − h(φ)

= F−(τ + φ, ǫ) − F−(τ − φ, ǫ) − φ (5.13)

≥ 2φ − φ = φ.

The equality holds just as φ = 0. Since φ ∈ In1 and R(φ) > φ, R(φ) ∈ In1 ∪ In2. Specifically, if φ = 0 ∈ In1, Rk_{(φ) = 0 ∈ In1.}

Proposition 5.1. There is no φ ∈ In1\{0} such that Rk_{(φ) ∈ In1\{0} for all k = 0, 1, 2, · · · .}

Proof. Suppose not, i.e., ∃φ ∈ In1\{0} s.t. Rk_{(φ) ∈ In1\{0} for all k = 0, 1, 2, · · · . Note that}

R(φ) is continuous in (0, τ − g(−ǫ)], since

R(φ) = (

F−(τ + φ, ǫ) − F−(τ − φ, ǫ) − φ, φ ∈ In1;

F−(τ + φ, ǫ) − φ φ ∈ In2a∨ In2b. (5.14, 5.18)

The return maps of In2a and In2b will be proved in Eqs. (5.14), (5.18).

⇒ (

R((τ − g(−ǫ))−_{) = F}

−(2τ − g(−ǫ)) − (τ − g(−ǫ));

R(τ − g(−ǫ)) = F−(2τ − g(−ǫ)) − (τ − g(−ǫ)).

Moreover, R(φ) > φ, ∀φ ∈ In1\{0}. Thus, Rk_{(φ) → φ}∗ _{for some φ}∗ _{∈ (0, τ − g(−ǫ)], and}

such φ∗ is a fixed point. That is R(φ∗) = φ∗. It follows that φ∗ = τ − g(−ǫ), or ǫ = −f (τ − φ∗), then F−(τ − φ∗, ǫ) = 0. Hence, R(φ∗) = φ∗ ⇒ F−(τ + φ∗, ǫ) − F−(τ − φ∗, ǫ) − φ∗ = φ∗ ⇒ F−(τ + φ∗, ǫ) − φ∗ = φ∗ ⇒ F−(τ + φ∗, ǫ) = 2φ∗ ⇒ f (τ + φ∗) + ǫ = f (2φ∗) ⇒ f (τ + φ∗) − f (τ − φ∗) = f (2φ∗) − f (0) = f (φ∗+ φ∗) − f (φ∗− φ∗)

By relation A6, we know f (τ + φ∗) − f (τ − φ∗) > f (φ∗+ φ∗) − f (φ∗− φ∗), since τ > φ∗. →←  Corollary 5.1. For each φ ∈ In1\{0}, there is a k ∈ N, depending on φ, such that Rk_{(φ) ∈ In2.}

Dynamics in In2a : −f (τ − φ) ≥ ǫ ≥ −f (τ ) and φ ∈ I1 (For fixed ǫ, φ ∈ [l2a, τ ], where l2a = τ − g(−ǫ).)

h(φ) = 1 − F−(τ + φ, ǫ) + φ

By Eq. (5.4), thus h : In2a7→ In0. Then the return map satisfies

R(φ) = h(h(φ)) = 1 − h(φ)

= F−(τ + φ, ǫ) − φ (5.14)

< τ + φ − φ = τ.

Thus R(φ) ∈ In1 ∪ In2.

Lemma 5.1. For each fixed ǫ, there exists a fixed point for R in In2a.

Proof. Note that R(φ) is continuous in [τ − g(−ǫ), τ ]. If φ = l2a, we know that ǫ = −f (τ − l2a). The return map

R(l2a) = F−(τ + l2a, ǫ) − l2a

= g(f (τ + l2a) − f (τ − l2a)) − l2a

> g(f (l2a+ l2a) − f (l2a − l2a)) − l2a (A6)

= 2l2a − l2a = l2a. (5.15)

If φ = τ , the return map

R(τ ) = F−(2τ, ǫ) − τ

< 2τ − τ = τ. (5.16)

By Eqs. (5.15), (5.16), there exists a fixed point for R in In2a.

 Since R(φ) = φ, R(φ) = φ ⇒ F−(τ + φ, ǫ) − φ = φ ⇒ F−(τ + φ, ǫ) = 2φ ⇒ g(f (τ + φ) + ǫ) = 2φ ⇒ ǫ = f (2φ) − f (τ + φ). (5.17)

then the fixed point for R in In2a is unique.

For any φ ∈ In2a,

R′(φ) = F₋′(τ + φ, ǫ) − 1

= g′(f (τ + φ) − |ǫ|)f′(τ + φ) − 1 > g′(f (τ + φ))f′(τ + φ) − 1 > 1 − 1 = 0.

Thus, l2a < R(l2a) < R(φ) < R(τ ) < τ . Hence, In2a is a trapping region. Let φ∗ be the fixed point of R(φ), i.e., R(φ∗) = φ∗.

(i)As φ > φ∗_,

R(φ) < φ (Otherwise, there exists another fixed point in [φ, τ ].) ⇒ Rk(φ) < · · · < R2(φ) < R(φ) < φ, ∀k = 1, 2, · · ·

and Rk(φ) > Rk(φ∗) = φ∗, ∀k = 1, 2, · · · (∵ R is increasing.) ⇒ φ∗ < Rk(φ) < Rk−1(φ) < · · · < φ, ∀k = 1, 2, · · ·

Thus, the sequence {Rk_{(φ)} is decreasing, and R}k_{(φ) → φ}∗_.

(ii)As φ < φ∗_,

Similarly, the sequence {Rk(φ)} is increasing, and Rk(φ) → φ∗. The unique fixed point φ∗ is an attractor.

 Corollary 5.2. For each φ ∈ In2a, Rk(φ) converges to the fixed point φ∗. Specifically, if φ = 0,

then the fixed point φ∗ _{= 0, i.e., φ}

A and φB synchronize.

Dynamics in In2b : −f (τ ) > ǫ > f (φ) − f (τ + φ) and φ ∈ I1 (For fixed ǫ, φ ∈ [l2b, τ ), where l2b satisfies ǫ = f (φ) − f (τ + φ).)

h(φ) = 1 − F−(τ + φ, ǫ) + φ

Remark 5.1. φ = 0 6∈ In2b.

By Eq. (5.4), thus h : In2b7→ In0. Then the return map satisfies

R(φ) = h(h(φ)) = 1 − h(φ)

= F−(τ + φ, ǫ) − φ (5.18)

< τ + φ − φ = τ.

Since −f (τ ) > ǫ, R(φ) 6∈ In1 ∪ In2a. It follows that R(φ) ∈ In2\In2a∪ In3.

Proposition 5.3. There is no φ ∈ In2b such that Rk(φ) ∈ In2b for all k = 0, 1, 2, · · · .

Proof. Note that R(φ) is continuous in [l2b, τ ]. Assume that there exists a fixed point for R in In2b. Since the return map of In2b is the same as that in In2a. It follows that

R(φ) = φ

⇔ ǫ = f (2φ) − f (τ + φ) > f (2φ − 2φ) − f (τ + φ − 2φ) (5.17, A5) ⇔ ǫ = f (2φ) − f (τ + φ) > −f (τ − φ) > −f (τ ) (For φ 6= 0 ∈ In2b.)

⇔ ǫ > −f (τ ).

This contradicts the region of the ǫ in In2b. Then there is no fixed point for R in In2b.

Since R(l2b) = 0 and R(τ ) = F−(2τ, ǫ) − τ < τ , it follows that R(φ) < φ. Otherwise, there exists fixed points in [l2b, τ ]. Thus, φ is decreasing in In2b.

Suppose not, i.e., ∃φ ∈ In2b, s.t. Rk(φ) ∈ In2b, for all k = 0, 1, 2, · · · . Since φ is decreasing,

Rk(φ) → φ∗, for some φ∗ ∈ [l2b, τ ]. It follows that φ

∗ _{= l}

2b is a fixed point in In2a. It makes contradiction.

 Corollary 5.3. For each φ ∈ In2b, there is a k ∈ N, depending on φ, such that Rk(φ) ∈

In2c∪ In2d∪ In3.

Dynamics in In2c : f (φ) − f (τ + φ) > ǫ > −f (τ + φ) and φ ∈ I1 (For fixed ǫ,

φ ∈ [g(−ǫ) − τ, min{l2b, τ }).)

R(φ) = h(φ) = 1 + F−(τ + φ, ǫ) − φ

Remark 5.2. φ = 0 6∈ In2c.

By Eq. (5.4), thus R, h : In2c 7→ In0. We iterate the return map and derive

h(R(φ)) = 1 − R(φ) = 1 − h(φ)

= φ − F−(τ + φ, ǫ) (5.19)

< φ.

Since R(φ) ∈ In0, h(R(φ)) < φ, h(R(φ)) ∈ In2c∪ In3.

Note that h(R(φ)) is continuous in In2c, and h(R(g(−ǫ) − τ )) = g(−ǫ) − τ . Clearly, we know

that h′(R(φ)) = 1 − F₋′(τ + φ, ǫ) < 0. For each φ ∈ In2c, since g(−ǫ) − τ < φ, we have

h(R(g(−ǫ) − τ )) > h(R(φ)) ⇒ g(−ǫ) − τ > h(R(φ))

Dynamics in In2d: ǫ = f (φ) − f (τ + φ) and φ ∈ I1 (For fixed ǫ, φ = l2b.) h(φ) = 1

Corollary 5.5. For each φ ∈ In2d, φA and φB synchronize.

Dynamics in In3 : ǫ ≤−f (τ + φ) and φ ∈ I1 (For fixed ǫ, φ ∈ [0, min{g(−ǫ) − τ, τ }).)

h(φ) = R(φ) = 1 − φ

Thus, h, R : In3 7→ In0. Specifically, if φ = 0, φA and φB synchronize.

Corollary 5.6. For each φ ∈ In3, φ and h(φ) is a Period.

Configuration I2

Dynamics in In4 : ǫ > f (2τ ) − f (τ + φ) and φ ∈ I2 (For fixed ǫ, φ ∈ (l4, 1 − τ ], where

l4 = g(f (2τ ) − ǫ) − τ .)

h(φ) = 1 − F−(τ + φ, ǫ) + τ

By Eqs. (5.7), (5.8), thus τ < h(φ) ≤ 1 − τ such that h : In4 7→ In4 ∪ In5 ∪ In6. Note that h is continuous in [l4, 1 − τ ]. Since,

h(φ) = ( 1 − F−(τ + φ, ǫ) + τ, φ ∈ In4; 1 − |F−(τ + φ, ǫ) − τ | φ ∈ In5. (5.9) ⇒ h(l + 4) = 1 − τ ; h(l4) = 1 − τ (5.20)

Lemma 5.2. There exists a unique fixed point for h in In4.

Proof. From Eq. (5.20), we know that h(l4) = 1 − τ > l4. It follows that

h(φ) = 1 − F−(τ + φ, ǫ) + τ ⇒ h′(φ) = −F₋′ (τ + φ, ǫ) < −1 (5.21) ⇒ Z 1−τ l4 h′(φ)dφ < Z 1−τ l4 −1dφ ⇒ h(1 − τ ) − h(l4) < −[(1 − τ ) − l4] ⇒ h(1 − τ ) < h(l4) − [(1 − τ ) − l4] = l4 < 1 − τ. (5.22)

Since h(l4) > l4 and h(1 − τ ) < 1 − τ , h(φ) has fixed points.

Now, we claim that the fixed point for h in In4 is unique.

Suppose not, i.e., ∃φ1< φ2 ∈ In4 s.t. h(φ1) = φ1 and h(φ2) = φ2. By Mean Value Theorem,

h′(φ′) = h(φ2) − h(φ1) φ2− φ1

for some φ′ _{∈ (φ} 1, φ2)

⇒ h′(φ′) = 1 for some φ′ ∈ (φ1, φ2)

This is a contradiction with Eq. (5.21). Thus, there exists a unique fixed point φ∗ _{for h in In4.}



Since there exists a fixed point for h,

h(φ) = φ

⇒ 1 − F−(τ + φ, ǫ) + τ = φ

⇒ 1 + τ − φ = g(f (τ + φ) + ǫ)

⇒ ǫ = f (1 + τ − φ) − f (τ + φ). (5.23)

Thus, the fixed point in In4 must satisfy ǫ(φ) = f (1 + τ − φ) − f (τ + φ).

Proposition 5.4. Except for φ = φ∗, there is no φ ∈ In4 such that hk_{(φ) ∈ In4 for all}

k = 0, 1, 2, · · · . i.e., The unique fixed point φ∗ for h in In4 is a repellor.

Proof. Suppose not, i.e., ∃φ ∈ In4 s.t. hk_{(φ) ∈ In4, for all k = 0, 1, 2, · · · . Let}

S1:= {φ : h1(φ) ∈ In4} ⇒ S1 is an interval, S1 6= ∅. (∵ φ∗ ∈ S1 and h′(φ) < −1)

S2:= {φ : h2(φ) ∈ In4} ⇒ S2 is an interval, S2 6= ∅. (∵ φ∗ ∈ S2 ⊆ S1 and

d dφh 2 (φ) > 1) .. . ... ⇒ S := {φ : hk_{(φ) ∈ In4 ∀k} =} ∞T k=1 Sk is an interval (noempty).

By Eq. (5.21), R′(φ) = h′(h(φ))h′(φ) > 1. Let H(φ) = R(φ) − φ, then we have H′(φ) = R′(φ) − 1 > 0. 1◦ (i) If φ < φ∗, H(φ) < H(φ∗) = 0 ⇒ R(φ) < φ. (ii) If φ > φ∗, H(φ) > H(φ∗_{) = 0 ⇒ φ < R(φ).} Then (i) Rk(φ) < · · · < R2 (φ) < R(φ) < φ < φ∗, ∀k ∈ N. and (ii) Rk_{(φ) > · · · > R}2 (φ) > R(φ) > φ > φ∗, ∀k ∈ N. 2◦ (i) If φ < φ∗_, H(R(φ)) < H(φ) ⇒ R(R(φ)) − R(φ) < R(φ) − φ ⇒ |R(φ) − φ| < |R2 (φ) − R(φ)| (ii) If φ > φ∗_, H(φ) < H(R(φ)) ⇒ |R(φ) − φ| < |R2 (φ) − R(φ)|

Corollary 5.7. For each φ ∈ In4, there exists a unique fixed point φ∗ which is a repellor. Also there is a k ∈ N, depending on φ, such that hk_{(φ) ∈ In5 ∪ In6.}

Dynamics in In5 : f (2τ ) − f (τ + φ) ≥ ǫ > −f (τ + φ) and φ ∈ I2

h(φ) = 1 − |F−(τ + φ, ǫ) − τ |

By Eq. (5.10), thus h : In5 7→ In0 ∪ l1−τ, where l1−τ is the line φ = 1 − τ . The firemap

h(φ) = 1 − τ ∈ l1−τ just as φ = l4. Also, if φA(τ ) > φB(τ ), R(φ) = h(φ). Furthermore, h(φ) = 1 ⇔ |F−(τ + φ, ǫ) − τ | = 0 ⇔ g(f (τ + φ) + ǫ) = τ ⇔ ǫ = f (τ ) − f (τ + φ). Therefore, φA and φB synchronize just as ǫ(φ) = f (τ ) − f (τ + φ).

Corollary 5.8. For each φ ∈ In5, h : In5 7→ In0 ∪ l1−τ. Specifically, when ǫ(φ) = f (τ ) − f (τ +

φ), φA and φB synchronize.

Dynamics in In6 : −f (τ + φ) ≥ ǫ and φ ∈ I2 (For fixed ǫ, φ ∈ [τ, g(−ǫ) − τ ].)

h(φ) = R(φ) = 1 − τ

Thus h : In6 7→ In4 ∪ In5.

Corollary 5.9. For each φ ∈ In6, h : In6 7→ In4 ∪ In5.

Configuration I3

Dynamics in In0 : φ ∈ I3

h(φ) = 1 − φ

Thus, h : In0 7→ I1.

Corollary 5.10. In0 is the same as Ex0 and we have h : In0 7→ I1.

5.3. Construction of Phase Diagrams.

Since for all φ ∈ In6, h(φ) = 1 − τ , and the line 1 − τ belongs to In4, In5, In6. Let’s find out how the oscillators run when φ = 1 − τ . We find that different τ leads to three different cases. Theorem 5.1. Case 1: _{τ ≤ g(}1 2) − 1 2 (Figure 5.1: τ = 0.25) ⇒ ∃G1, G2 s.t. f (2τ ) − 1 < G1 ≤ G2 < −f (2τ ), and G1 = G2 when τ = g( 1 2) − 1 2. (i) If ǫ ∈ [G1, G2], then h(1 − τ ) ∈ In6.

(ii)If ǫ ∈ [f (2τ ) − 1, G1) ∪ (G2, 0), then h(1 − τ ) ∈ In5.

(iii)If ǫ ∈ (−1, f (2τ ) − 1), then h(1 − τ ) ∈ In0. Case 2: _g(1 2) − 1 2 < τ ≤ 1 2g( 1 2) (Figure 5.2: τ = 0.3) ⇒ (i) If ǫ ∈ [f (2τ ) − 1, 0), then h(1 − τ ) ∈ In5.

(ii)If ǫ ∈ (−1, f (2τ ) − 1), then h(1 − τ ) ∈ In0. Case 3: _{τ >} 1

2g( 1

2) (Figure 5.3: τ = 0.4) ⇒ (i) If ǫ ∈ [f (2τ ) − 1, 0), then h(1 − τ ) ∈ In5.

(ii)If ǫ ∈ (−1, f (2τ ) − 1), then h(1 − τ ) ∈ In0.

Proof. We clearly know that the line 1 − τ ∈ In4, In5, In6. If φ ∈ In4, then ǫ > f (2τ ) − f (τ + φ) and φ ∈ I2, i.e., if φ ∈ In4, then ǫ > f (2τ ) − 1. On the other hand, If φ ∈ In6, then ǫ ≤ −f (τ + φ)

and φ ∈ I2, i.e., if φ ∈ In6, then ǫ < −f (2τ ).

Case 1: ∀φ = 1 − τ ∈ In4, by Eq. (5.20), h(1 − τ ) < l4 ⇒ h(1 − τ ) ∈ In5 ∨ In6. (5.24) (i) h(1 − τ ) ∈ In6 ⇔ ( f (2τ ) − 1 < ǫ < −f (2τ ) |ǫ| ≥ f (τ + h(1 − τ )).

On the one hand, since f (2τ ) − 1 < ǫ < −f (2τ ), f (2τ ) − 1 < −f (2τ ), or equivalently, τ < 1 2g(

1 2). On the other hand, since h(1 − τ ) = 1 − g(1 − ǫ) + τ ,

|ǫ| ≥ f (τ + h(1 − τ )) ⇔ −ǫ ≥ f (2τ + 1 − g(1 + ǫ)) ⇔ g(−ǫ) ≥ 2τ + 1 − g(1 + ǫ)

Here K′_{(ǫ) = 0 as ǫ = −}1

2, and the maximum value of K is K(− 1 2) = 2g( 1 2). Since K(− 1 2 + x) = K(−1 2 − x) i.e., K is symmetric at ǫ = − 1 2 and [f (2τ ) − 1] + [−f (2τ )] = 2 · (− 1 2), for a fixed τ , ∃ǫ ∈ (f (2τ ) − 1, −f (2τ )) s.t. g(−ǫ) + g(1 + ǫ) ≥ 2τ + 1, if and only if 2g(1 2) ≥ 2τ + 1, i.e., τ ≤ g(1 2) − 1 2 ≤ 1 2g( 1

2). Moreover, the collection of such ǫ is an interval.

Thus, there is a ǫ such that h(1 − τ ) ∈ In6 with 1 − τ ∈ In4 if and only if τ ≤ g(1 2) − 1 2, and as τ ≤ g(1 2) − 1

2, set S = {1 − τ ∈ In4, h(1 − τ ) ∈ In6} = [G1, G2] 6= ∅ for some G1, G2 ∈ (f (2τ ) − 1, −(f τ )).

(ii)If ǫ ∈ [f (2τ )− 1, G1)∪ (G2, 0), then 1− τ ∈ In4. By Eq. (5.24) and (i), we obtain h(1− τ ) ∈ In5.

(iii)If ǫ ∈ (−1, f (2τ ) − 1), then 1 − τ ∈ In5 by the definition of In5. Thus, h(1 − τ ) ∈ In0. Case 2:

If τ ≤ 1 2g(

1

2), then f (2τ ) − 1 ≤ −f (2τ ). Hence, ∄[G1, G2] ∈ (f (2τ ) − 1, −f (2τ )). (i) The proof is the same as Case 1-(ii).

(ii)The proof is the same as Case 1-(iii). Case 3:

If τ > 1 2g(

1

2), thenf (2τ ) − 1 > −f (2τ ). (i) The proof is the same as Case 1-(ii).

(ii)The proof is the same as Case 1-(iii). 

In the following Figure 5.1 ∼ 5.3, we use the particular energy function f and inverse function g which are used in Section 4.3 (The inverse of Eq. (4.21)). Also, we fix s := 1, r := 0.9, and three different τ in Figure 5.1 ∼ 5.3.

Moreover, we add Figure 5.4 which is the phase diagram we make from Ernst [10] as a contrast. And it’s energy function f and inverse function g is exactly the Mirollo-Strogatz-type (Eq. 4.21) oscillators we mentioned before. Also, we fixed s := 1, r := 0.4, and τ = 0.2 in Figure 5.4.

And there are some notations for the following figures: l1 : ǫ = −f (τ − φ)

l2 : ǫ = −f (τ + φ)

l3 : ǫ = f (2τ ) − f (τ + φ)

D1 : ∀φ ∈ D1, φ is an attractor.

D2 : ∀φ ∈ D2, the complete synchronization occurs.

D3 : ∀φ ∈ D3, φ is a repellor.

D4 : ∀φ ∈ D4, the complete synchronization occurs.

In the following tables, we fixed an arbitrary ǫ < 0, and shift φ in [0,1]. Therefore, we clearly see which area (InN ) does the oscillator begin and figure out where it ends.

(I) f (φ) = ln( s−rφ s ) ln(s−r_s ) , f ′_{> 0, f}′′_{> 0} Case 1: τ ≤ g(1 2) − 1 2 In2_a In2c

In6

In4

=In2_d l2 l3 D3 D2 In2b

In 1

l1 D1 D4

In5

In0

In3

-f( )

τ

f(2 )-1

τ

-f(2 )

τ

-f(2 ) f( )

τ

+

τ

G

₂

G

₁ Figure 5.1. _{s := 1, r := 0.9, τ = 0.25, where s > r > 0} ǫ Dynamics Multistability −f (τ ) ≤ ǫ < 0 φ→ φ ∗

I n2a, if φ ∈ I1, I2\D3, I3 Lag Syn. with Lag φ ∗ I n2a φ→ φ, if φ ∈ D3 Lag Syn. with Lag φ

G2< ǫ <−f (τ )

φ→Period(φ∗

I n3,1 − φ∗I n3), if φ ∈ I1\D2, I2\{D3∪ D4}, I3 Lag Syn. with Lag φ∗I n3

φ→ 0, if φ ∈ D2, D4 Complete Syn.

φ→ φ, if φ ∈ D3 Lag Syn. with Lag φ

φ→Period(φ∗

I n3,1 − φ∗I n3), if φ ∈ I1, I3 Lag Syn. with Lag φ∗I n3

φ→Period(φ∗

I n3,1 − φ∗I n3), if φ ∈ In4\{l1−τ ∪ D3}, In5\D4 Lag Syn. with Lag φ∗I n3

G1 ≤ ǫ ≤ G2 φ→Period(1 − τ, h(1 − τ )), if φ ∈ l_1−τ, In6 Lag Syn. with Lag 1 − τ

φ→ φ, if φ ∈ D3 Lag Syn. with Lag φ

φ→ 0, if φ ∈ D4 Complete Syn.

φ→Period(φ∗

Case 2: g(1 2) − 1 2 < τ ≤ 1 2g( 1 2)

In1

In0

In3

In2

_a

In6

In4

In5

In2

c

In2

_b

= In2

_d

-f( )

τ

-f(2 )+

τ

f( )

τ

-f(2 )

τ

f(2 )-1

τ

l1

l2

l3

D3

D2

D4

D

1 Figure 5.2. _{s := 1, r := 0.9, τ = 0.3, where s > r > 0} ǫ Dynamics Multistability −f (τ ) ≤ ǫ < 0 φ→ φ ∗

I n2a, if φ ∈ I1, I2\D3, I3 Lag Syn. with Lag φ ∗ I n2a φ→ φ, if φ ∈ D3 Lag Syn. with Lag φ

f(2τ ) − 1 < ǫ < −f (τ ), −1 ≤ ǫ < f (2τ ) − 1

φ→Period(φ∗

I n3,1 − φ∗I n3), if φ ∈ I1\D2, I2\{D3∪ D4}, I3 Lag Syn. with Lag φ∗I n3

φ→ 0, if φ ∈ D2, D4 Complete Syn.

φ→ φ, if φ ∈ D3 Lag Syn. with Lag φ

ǫ= f (2τ ) − 1

φ→ 1 − τ , if φ ∈ l1−τ, In6 Lag Syn. with Lag 1 − τ

φ→Period(φ∗

I n3,1 − φ∗I n3), if φ ∈ I1, I2\{D4∪ l_1−τ}, I3 Lag Syn. with Lag φ∗I n3

φ→ 0, if φ ∈ D4 Complete Syn.

Case 3: τ > 1 2g( 1 2)

In1

In0

In3

In2

_a

In6

In4

In5

In2

_c

In2

_b

-f( )

τ

-f(2 )+

τ

f( )

τ

-f(2 )

τ

f(2 )-1

τ

l1

l2

l3

D

3

= In2

_d

D2

D4

D1

Figure 5.3. _{s := 1, r := 0.9, τ = 0.4, where s > r > 0} ǫ Dynamics Multistability −f (τ ) ≤ ǫ < 0 φ→ φ ∗

I n2a, if φ ∈ I1, I2\D3, I3 Lag Syn. with Lag φ ∗ I n2a φ→ φ, if φ ∈ D3 Lag Syn. with Lag φ

f(2τ ) − 1 < ǫ < −f (τ ), −1 ≤ ǫ < f (2τ ) − 1

φ→Period(φ∗

I n3,1 − φ∗I n3), if φ ∈ I1\D2, I2\{D3∪ D4}, I3 Lag Syn. with Lag φ∗I n3

φ→ 0, if φ ∈ D2, D4 Complete Syn.

(II) f (φ) = s r − s r( s − r s ) φ_{, f}′ _{> 0, f}′′_{< 0} l3

In1

In3

In5

-f( )τ f(2 )-1τ -f(2 )τ

In4

In0

In6

l2

In2

l1 Figure 5.4. _{s := 1, r := 0.4, τ = 0.2, where s > r > 0} ǫ Dynamics Multistability −f (τ ) ≤ ǫ < 0 φ → 0, if φ ∈ I1, In5, I3 Complete Syn.

φ → φ∗

I n4, if φ ∈ In4 Lag Syn. with Lag φ∗I n4

f (2τ ) − 1 < ǫ < −f (τ ) φ → φ

∗

I n3, if φ ∈ I1, In5, I3 Lag Syn. with Lag φ∗I n3

φ → φ∗

I n4, if φ ∈ In4 Lag Syn. with Lag φ∗I n4

−1 ≤ ǫ < f (2τ ) − 1 φ → φ∗

I n3, if φ ∈ I1, I2, I3 Lag Syn. with Lag φ∗I n3

6. Conclusion

It was numerically demonstrated in [10] that with N ≫ 2 convex oscillators, the system reveal multistable phase clustering for inhibitory couplings. For N ≫ 2 concave oscillators the corresponding system also has multistable phase clustering for excitatory couplings. Since for N = 2, the corresponding system has stable in-phase synchronization. It will be interesting to treat N as a parameter so as to see how the system evolves from synchronization to clustering synchronization as N increases.

References

[1] S. Bottani, Pulse-Coupled Relaxation Oscillators: From Biological Synchronization to Self-Organized Critical-ity, Phys. Rev. Lett., 74, 4189 (1995).

[2] J. Buck, Synchronous Rhythmic Flashing of Fireflies. II., Q. Rev. Biol., 63, 265 (1988).

[3] Y. C. Chang, and J. Juang, Stable Synchrony In Globally-Coupled Integrate-And-Fire Oscillators, preprint. [4] H. Daido, Lower Critical Dimension for Populations of Oscillators with Randomly Distributed Frequencies: A

Renormalization-Group Analysis, Phys. Rev. Lett., 61, 231 (1988).

[5] J. Deppisch, H. U. Bauer, T. Schillen, P. K¨onig, K. Pawelzik, and T. Geisel, Network 4, 243 (1993).

[6] R. Eckhorn, R. Bauer, W. Jordan, M. Brosch, W. Kruse, M. Munk, and H. J. Reitboeck, Biol. Cybern. 60, 121 (1988).

[7] A. K. Engel, P. K¨onig, A. K. Kreiter, and W. Singer, Science 252, 1177 (1991).

[8] G. Ermentrout, An adaptive model for synchrony in the firefly Pteroptyx malaccae, J. Math. Biol., 29, 571 (1991).

[9] U. Ernst, K. Pawelzik, and T. Geisel, Phys. Rev. Lett. 74, 1570 (1995).

[10] U. Ernst, K. Pawelzik, and T. Geisel, Delay-induced multistable synchronization of biological oscillators, Phys. Rev. E, 57, 2150 (1998).

[11] C. M. Gray, P. Kr¨onig, A. K. Engel, and W. Singer, Nature (London) 338, 334 (1989).

[12] J. Jalife, Mutual entrainment and electrical coupling as mechanisms for synchronous firing of rabbit sinoatrial pacemaker cells, J. Physiol., 356, 221 (1984).

[13] Y. Kuramoto and I. Nishikawa, Statistical macrodynamics of large dynamical systems. Case of a phase tran-sition in oscillator communities, J. Stat. Phys., 49, 596 (1987).

[14] Y. Kuramoto, Collective synchronization of pulse-coupled oscillators and excitable units, Physica D, 50, 15 (1991).

[15] W. W. Lytton and T. J. Sejnowski, J. Neurophys. 66,1059 (1991); C. van Vreeswijk, L. F. Abbott, and G. B. Ermentrout, J. Comp. Neurosci. (to be published); A. Nischwitz, H. Gl¨under, A. von Oertzen, and P. Klausner, ANN Proceedings (Springer-Verlag, Berlin, 1993), p. 622.

[16] C. v. d. Malsburg, Technical Report No. 81-2, MPI for Biophysical Chemistry, 1981. [17] C. v. d. Malsburg and W. Schneider, Biol. Cybern. 54, 29 (1986).

[18] D. C. Michaels, E. P. Matyas, J. Jalife, Mechanisms of sinoatrial pacemaker synchronization: a new hypoth-esis, Circulation Res., 61, 704 (1987).

[19] R. E. Mirollo and S. H. Strogatz, Synchronization of pulse-coupled biological oscillators, SIAM J. Appl. Math., 50_{, 1645 (1990).}

[20] C. S. Peskin, Mathematical Aspects of Heart Physiology, Courant Inst. Math. Sci., New York Univ., 268 (1975).

[21] J. Rubin and D. Terman, Geometric analysis of population rhythms in synaptically coupled neuronal networks, Neural Computation, 12, 597 (2000).

[22] A. Sherman, J. Rinzel, J. Keizer, Emergence of organized bursting in clusters of pancreatic beta-cells by channel sharing, Biophys. J., 54, 411 (1988).

[23] A. Sherman and J. Rinzel., Collective properties of insulin secreting cells, in “Cell to Cell Signalling: From Experiments to Theoretical Models”, A. Goldbeter, ed., Academic Press, London, 61 (1989).

[24] S. H. Strogatz, C. M. Marcus, R. M. Westervelt, and R. E. Mirollo, Simple Model of Collective Transport with Phase Slippage, Phys. Rev. Lett., 61, 2380 (1988).

[25] S. Strogatz, Lecture Notes in Biomathematics, Springer, Berlin, 100, (1993).

[26] Vincent Torre, A theory of synchronization of heart pace-maker cells, J. Theoret. Biol., 61, 55 (1976). [27] A. T. Winfree, The Geometry of Biological Time (Springer-Verlag, New York, 1980).

[28] C. W. Wu, Synchronization in Coupled Chaotic Circuits and Systems, World Scientific series on nonlinear science, Series A, 41, World Scientific, Singapore, (2002).