Emergence of synchronization in Kuramoto model with frustration under general network topology

1.   Introduction

Synchronized behavior in complex systems is ubiquitous and has been extensively investigated in various academic communities such as physics, biology, engineering [2,7,29,36,37,39,40,42], etc. Recently, sychronization mechanism has been applied in control of robot systems and power systems [12,13,34]. The rigorous mathematical treatment of synchronization phenomena was started by two pioneers Winfree [43] and Kuramoto [27,28] several decades ago, who introduced different types of first-order systems of ordinary differential equations to describe the synchronous behaviors. These models contain rich emergent behaviors such as synchronization, partially phase-lcoking and nonlinear stability, etc., and have been extensively studied in both theoretical and numerical level [3,5,11,14,15,16,17,19,26,32,39].

In this paper, we address the synchronous problem of Kuramoto model on a general graph under the effect of frustration. To fix the idea, we consider a digraph $\mathcal{G} = (\mathcal{V}, \mathcal{E})$ consisting of a finite set $\mathcal{V} = \{1,\ldots, N\}$ of vertices and a set $\mathcal{E} \subset \mathcal{V} \times \mathcal{V}$ of directed arcs. We assume that Kuramoto oscillators are located at vertices and interact with each other via the underlying network topology. For each vertex $i$, we denote the set of its neighbors by $\mathcal{N}_i$, which is the set of vertices that directly influence vertex $i$. Now, let $\theta_i = \theta_i(t)$ be the phase of the Kuramoto oscillator at vertex $i$, and define the $(0,1)$-adjacency matrix $(\chi_{ij})$ as follows:

Then, the set of neighbors of $i$-th oscillator is actually $\mathcal{N}_i : = \{j: \chi_{ij} >0\}$. In this setting, the Kuramoto model with frustration on a general network is governed by the following ODE system:

where $\Omega_i, K, N$ and $\alpha \in (0, \frac{\pi}{2})$ are the natural frequency of the $i$th oscillator, coupling strength, the number of oscillators and the uniform frustration between oscillators, respectively. For the case of nonpositive frustration, we can reformulate such a system into (1) form by taking $\hat{\theta}_i = - \theta_i$ for $i = 1,2,\ldots, N$. Note that the well-posedness of system (1) is guaranteed by the standard Cauchy-Lipschitz theory since the vector field on the R.H.S of (1) is analytic.

Comparing to the original Kuramoto model, there are two additional structures, i.e., frustration and general digraph. The frustration was introduced by Sakaguchi and Kuramoto [38], due to the observation that a pair of strongly coupled oscillators eventually oscillate with a common frequency that deviates from the average of their natural frequencies. On the other hand, the original all-to-all symmetric network is an ideal setting, thus it is natural to further consider general digraph case. Therefore, the frustration model with general digraph is more realistic in some sense. Moreover, these two structures also lead to richer phenomenon. For instance, the author in [6] observed that the frustration is common in disordered interactions, and the author in [44] found that frustration can induce the desynchronization through varying the value of $\alpha$ in numerical simulations. For more information, please refer to [4,10,21,23,24,25,30,33,35,41].

However, mathematically, for the Kuramoto model, the frustration and general digraph structures generate a lot of difficulties in rigorous analysis. For instance, the conservation law and gradient flow structure are lost, and thus the asymptotic states and dissipation mechanism become non-trivial. For all-to-all and symmetric case with frustration, in [20], the authors provided sufficient frameworks leading to complete synchronization under the effect of uniform frustration. In their work, they required initial configuration to be confined in half circle. Furthermore, the authors in [31] dealt with the stability and uniqueness of emergent phase-locked states. In particular, the authors in [22] exploited order parameter approach to study the identical Kuramoto oscillators with frustration. They showed that an initial configuration whose order parameter is bounded below will evolve to the complete phase synchronization or the bipolar state exponentially fast. On the other hand, for non-all-to-all case without frustration, the authors in [9] lifted the Kuramoto model to second-order system such that the second-order formulation enjoys several similar mathematical structures to that of Cucker-Smale flocking model [8]. But this method only works when the size of initial phases is less than a quarter circle, as we know the cosine function becomes negative if $\frac{\pi}{2}< \theta<\pi$. In [45], for the Kuramoto model on general network whose initial configuration is distributed in half circle, the authors exploited the idea of node decomposition in [25] and showed that frequency synchronization emerges under sufficient frameworks. To the best knowledge of the authors, there is few work on the Kuramoto model over general digraph with frustration. The authors in [18] studied the Kuramoto model with frustrations on a complete graph which is a small perturbation of all-to-all network, and provided synchronization estimates in half circle.

Our interest in this paper is studying the system (1) with uniform frustration on a general digraph. As far as the authors know, when the ensemble is distributed in half circle, the dissipation structure of the Kuramoto type model with general digragh is still unclear. The main difficulty comes from the loss of uniform coercive inequality, which is due to the non-all-to-all and non-symmetric interactions. Thus we cannot expect to capture the dissipation from Gronwall-type inequality of phase diameter. For example, the time derivative of the phase diameter may be zero at some time for general digraph case. To this end, we switch to apply the idea of node decomposition in [25] to gain the dissipation through hypo-coercivity. Due to the lack of monotonicity of sine function in half circle, we follow the delicate constructions and estimates of the convex combinations in [45]. Eventually, we have the following main theorem.

Theorem 1.1. Suppose the network topology $(\chi_{ij})$ in system (1) contains a spanning tree, $D^\infty$ is a given positive constant such that $D^\infty <\frac{\pi}{2}$, and all oscillators are initially confined in half circle, i.e.,

Then for sufficient large coupling strength $K$ and small frustration $\alpha$, there exists a finite time $t_*>0$ such that

and

where $C_1$ and $C_2$ are positive constants.

Remark 1. The first part of Theorem 1.1 claims that all oscillators confined initially in half circle will enter into a region less than quarter cirlce after some finite time. It is natural to ask how large $K$ and how small $\alpha$ we need to guarantee the first part of Theorem 1.1. In fact, according to the proof in later sections, we have the following explicit constraints on $K$ and $\alpha$,

where $d$ is the number of general nodes which is smaller than $N$ (see Lemma 2.7), $D(\Omega)$ is the diameter of natural frequency, and the other parameters $\zeta$, $\gamma$, $\eta$, $\beta$ and $c$ are positive constants which satisfy the following properties,

where $A(2N,j)$ denotes the number of permutations. It's obvious that we can find admissible parameters satisfying (3) since $D(\theta(0)) < \pi$. Once the parameters are fixed, we immediately conclude (2) holds for small $\alpha$ and large $K$.

Remark 2. For $t \in [0,t_*)$, the phase diameter $D(\theta(t))$ satisfies $D(\theta(t)) < \zeta$ which follows from $t_* < \bar{t}$ in Lemma 4.5 and Lemma 2.8. After $t_*$, all oscillators are confined in a region less than quarter circle, and Kuramoto model (1) will be equivalent to Cucker-Smale type model with frustration (see (53)). Therefore, we can directly apply the methods and results in [9] to conclude the emergence of frequency synchronization for large coupling and small frustration (see Lemma 4.6). Therefore, to guarantee the emergence of synchronization, it suffices to show the detailed proof of the first part of Theorem 1.1.

The rest of the paper is organized as follows. In Section 2, we recall some concepts on the network topology and provide a priori local-in-time estimate on the phase diameter of whole ensemble with frustration. In Section 3, we consider a strongly connected ensemble with frustration for which the initial phases are distributed in a half circle. We show that for large coupling strength and small frustration, the phase diameter is uniformly bounded by a value less than $\frac{\pi}{2}$ after some finite time. In Section 4, we study the general network with a spanning tree structure under the effect of uniform frustration. We use the inductive argument and show that Kuramoto oscillators will concentrate into a region less than quarter circle in finite time, which eventually leads to the emergence of synchronization exponentially fast. In Section 5, we present some simulations to illustrate the main results in Theorem 1.1. Section 6 is devoted to a brief summary.

2.   Preliminaries

In this section, we first introduce some fundamental concepts such as synchronization, spanning tree and node decomposition of a general network (1). Then, we will provide some necessary notations and a priori estimate that will be frequently used in later sections.

For the Kuramoto-type model, we recall the definition of synchronization as follows.

Definition 2.1. Let $\theta(t) : = (\theta_1(t), \theta_2(t), \ldots, \theta_N(t))$ be a time-dependent Kuramoto phase vecor. The configuration $\theta(t)$ exhibits (asymptotic) complete frequency synchronization if and only if the relative frequencies tend to zero asymptotically:

2.1.   Directed graph

Let the network topology be registered by the neighbor set $\mathcal{N}_i$ which consists of all neighbors of the $i$th oscillator. Then, for a given set of $\{\mathcal{N}_i\}_{i = 1}^N$ in system (1), we have the following definition.

Definition 2.2. (1)The Kuramoto digraph $\mathcal{G} = (\mathcal{V}, \mathcal{E})$ associated to system (1) consists of a finite set $\mathcal{V} = \{1,2,\ldots,N\}$ of vertices, and a set $\mathcal{E} \subset \mathcal{V} \times \mathcal{V}$ of arcs with ordered pair $(j,i) \in \mathcal{E}$ if $j \in \mathcal{N}_i$.

(2)A path in $\mathcal{G}$ from $i_1$ to $i_k$ is a sequence $i_1,i_2,\ldots,i_k$ such that

If there exists a path from $j$ to $i$, then vertex $i$ is said to be reachable from vertex $j$.

(3)The Kuramoto digraph contains a spanning tree if we can find a vertex such that any other vertex of $\mathcal{G}$ is reachable from it.

In order to guarantee the emergence of synchronization, we will always assume the existence of a spanning tree throughout the paper. Now we recall the concepts of root and general root introduced in [25]. Let $l,k \in \mathbb{N}$ with $1 \le l \le k \le N$, and let $C_{l,k} = (c_l,c_{l+1},\ldots,c_k)$ be a vector in $\mathbb{R}^{k-l+1}$ such that

For an ensembel of $N$-oscillators with phases $\{\theta_i\}_{i = 1}^N$, we set $\mathcal{L}_l^k(C_{l,k})$ to be a convex combination of $\{\theta_i\}_{i = l}^k$ with the coefficient $C_{l,k}$:

Note that each $\theta_i$ is a convex combination of itself, and particularly $\theta_N = \mathcal{L}_N^N(1)$ and $\theta_1 = \mathcal{L}_1^1(1)$.

Definition 2.3. (Root and general root)

1. We say $\theta_k$ is a root if it is not affected by the rest oscillators, i.e., $j \notin \mathcal{N}_k$ for any $j \in \{1,2,\ldots,N\} \setminus \{k\}$.

2. We say $\mathcal{L}_l^k(C_{l,k})$ is a general root if $\mathcal{L}_l^k(C_{l,k})$ is not affected by the rest oscillators, i.e., for any $i\in \{l, l+1,\ldots,k\}$ and $j \in \{1,2,\ldots,N\} \setminus \{l,l+1,\ldots, k\}$, we have $j \notin \mathcal{N}_i$.

Lemma 2.4.[25] The following assertions hold.

1. If the network contains a spanning tree, then there is at most one root.

2. Assume the network contains a spanning tree. If $\mathcal{L}_k^N(C_{k,N})$ is a general root, then $\mathcal{L}_1^l(C_{1,l})$ is not a general root for each $l \in \{1,2,\ldots,k-1\}$.

2.2.   Node decomposition

In this part, we will recall the concept of maximum node introduced in [25]. Then, we can follow node decomposition introduced in [25] to represent the whole graph $\mathcal{G}$ (or say vertex set $\mathcal{V}$) as a disjoint union of a sequence of nodes. The key point is that the node decomposition shows a hierarchical structure, then we can exploit this advantage to apply the induction principle. Let $\mathcal{G} = (\mathcal{V}, \mathcal{E}), \mathcal{V}_1 \subset \mathcal{V}$, and a subgraph $\mathcal{G}_1 = (\mathcal{V}_1, \mathcal{E}_1)$ is the digraph with vertex set $\mathcal{V}_1$ and arc set $\mathcal{E}_1$ which consists of the arcs in $\mathcal{G}$ connecting agents in $\mathcal{V}_1$. For a given digraph $\mathcal{G} = (\mathcal{V}, \mathcal{E})$, we will identify a subgraph $\mathcal{G}_1 = (\mathcal{V}_1, \mathcal{E}_1)$ with its vertex set $\mathcal{V}_1$ for convenience. Now we first present the definition of nodes below.

Definition 2.5. [25] (Node) Let $\mathcal{G}$ be a digraph. A subset $\mathcal{G}_1$ of vertices is called a node if it is strongly connected, i.e., for any subset $\mathcal{G}_2$ of $\mathcal{G}_1$, $\mathcal{G}_2$ is affected by $\mathcal{G}_1\setminus\mathcal{G}_2$. Moreover, if $\mathcal{G}_1$ is not affected by $\mathcal{G}\setminus\mathcal{G}_1$, we say $\mathcal{G}_1$ is a maximum node.

Intuitively, a node can be understood through a way that a set of oscillators can be viewed as a "large" oscillator. The concept of node can be exploited to simplify the structure of the digraph, which indeed helps us to catch the attraction effect more clearly in the network topology.

Lemma 2.6.[25] Any digraph $\mathcal{G}$ contains at least one maximum node. A digraph $\mathcal{G}$ contains a unique maximum node if and only if $\mathcal{G}$ has a spanning tree.

Lemma 2.7.[25] (Node decomposition)Let $\mathcal{G}$ be any digraph. Then we can decompose $\mathcal{G}$ to be a union as $\mathcal{G} = \bigcup_{i = 0}^d (\bigcup_{j = 1}^{k_i} \mathcal{G}_i^j)$ such that

1. $\mathcal{G}_0^j$ are the maximum nodes of $\mathcal{G}$, where $1\le j\le k_0$.

2. For any $p,q$ where $1 \le p\le d$ and $1\le q \le k_p$, $\mathcal{G}_p^q$ are the maximum nodes of $\mathcal{G} \setminus (\bigcup_{i = 0}^{p-1} (\bigcup_{j = 1}^{k_i} \mathcal{G}_i^j))$.

Proof. As $\mathcal{G}$ is assumed to be any digraph, from Lemma 2.2, we see that $\mathcal{G}$ contains at least one maximum node. Therefore, we can find all maximum nodes of $\mathcal{G}$ and denote them by $\mathcal{G}^j_0$ with $1\le j \le k_0$, where $k_0$ is the number of maximum nodes. Next, we get rid of $\bigcup_{j = 1}^{k_0} \mathcal{G}_0^j$ and collect all maximum nodes in the remaining digraph $\mathcal{G}\setminus (\bigcup_{j = 1}^{k_0} \mathcal{G}_0^j)$. Denote all maximum nodes of the remainder $\mathcal{G}\setminus (\bigcup_{j = 1}^{k_0} \mathcal{G}_0^j)$ by $\mathcal{G}_1^j$ with $1 \le j \le k_1$, provided that there are $k_1$ maximum nodes for $\mathcal{G}\setminus (\bigcup_{j = 1}^{k_0} \mathcal{G}_0^j)$. Then we can repeat the same process and find the maximum nodes $\mathcal{G}^q_p$ of $\mathcal{G}\setminus (\bigcup_{i = 0}^{p-1}(\bigcup_{j = 1}^{k_i}\mathcal{G}_i^j))$ with $1 \le q \le k_p$. After $d$ steps, we can construct $\mathcal{G} = (\bigcup_{i = 0}^{d}(\bigcup_{j = 1}^{k_i}\mathcal{G}_i^j))$ as $\mathcal{G}$ consists of finite $N$ oscillators.

Remark 3. Lemma 2.7 shows a clear hierarchical structure on a general digraph. For the convenience of later analysis, we make some comments on important notations and properties that are used throughout the paper.

1. From the definition of maximum node, for $1\le q \ne q' \le k_p$, we see that $\mathcal{G}_p^q$ and $\mathcal{G}_p^{q'}$ do not affect each other. Actually, $\mathcal{G}_p^q$ will only be affected by $\mathcal{G}_0$ and $\mathcal{G}_i^j$, where $1 \le i \le p-1, \ 1 \le j \le k_i$. Thus in the proof of our main theorem (see Theorem 1.1), without loss of generality, we may assume $k_i = 1$ for all $1 \le i \le d$. Hence, the decomposition can be further simplified and expressed by

where $\mathcal{G}_p$ is a maximum node of $\mathcal{G}\setminus (\bigcup_{i = 0}^{p-1} \mathcal{G}_i)$.

2.For a given oscillator $\theta_i^{k+1} \in \mathcal{G}_{k+1}$, we denote by $\bigcup_{j = 0}^{k+1} \mathcal{N}_i^{k+1}(j)$ the set of neighbors of $\theta_i^{k+1}$, where $\mathcal{N}_i^{k+1}(j)$ represents the neighbors of $\theta_i^{k+1}$ in $\mathcal{G}_j$. Note that the node decomposition and spanning tree structure in $\mathcal{G}$ guarantee that $\bigcup_{j = 0}^k \mathcal{N}_i^{k+1}(j) \ne \emptyset$.

2.3.   Notations and local estimates

In this part, for notational simplicity, we introduce some notations, such as the extreme phase, phase diameter of $\mathcal{G}$ and the first $k+1$ nodes, natural frequency diameter, and cardinality of subdigraph:

Finally, we provide a priori local-in-time estimate on the phase diameter to finish this section, which states that all oscillators can be confined in half circle in short time.

Lemma 2.8. Let $\theta_i$ be a solution to system (1) and suppose the initial phase diameter satisfies

where $\zeta$ and $\gamma$ are positive constants less than $\pi$.Then there exists a finite time $\bar{t} >0$ such that the phase diameter of whole ensemble remains less than $\zeta$ before $\bar{t}$, i.e.,

Proof. From system (1), we see that the dynamics of extreme phases is given by the following equations

where $\theta_M$ and $\theta_m$ denote the maximum and minimum phases. We combine the above equations to estimate the dynamics of phase diameter

where $\Omega_M$ and $\Omega_m$ are the natural frequencies of the extreme phases $\theta_M$ and $\theta_m$ respectively, and we use the formula

When the phase diameter satisfies $D(\theta(t)) \le \zeta < \pi$, it is obvious that

Then we see from (4) that

That is to say, when $D(\theta(t)) \le \zeta < \pi$, the growth of phase diameter is no greater than the linear growth with positive slope $D(\Omega) + 2NK \sin \alpha$. Now we integrate on both sides of (5) from $0$ to $t$ to have

Therefore, it yields that there exists a finite time $\bar{t}>0$ such that

where $\bar{t}$ is given as below

3.   Strongly connected case

We will first study the special case, i.e., the network is strongly connected. Without loss of generality, we denote by $\mathcal{G}_0$ the strongly connected digraph. According to Definition 2.5, Lemma 2.6 and Lemma 2.7, the strongly connected network consists of only one maximum node. Then in the present special case, according to Remark 2, we only need to show that the phase diameter is uniformly bounded by a value less than $\frac{\pi}{2}$ after some finite time, which ultimately yields the emergence of frequency synchronization. We now introduce an algorithm to construct a proper convex combination of oscillators, which can involve the dissipation from interaction of general network. More precisely, the algorithm for $\mathcal{G}_0$ consists of the following three steps:

Step 1. For any given time $t$, we reorder the oscillator indexes to make the oscillator phases from minimum to maximum. More specifically, by relabeling the agents at time $t$, we set

In order to introduce the following steps, we first provide the process of iterations for $\bar{\mathcal{L}}_k^{N_0}(\bar{C}_{k,N_0})$ and $\underline{\mathcal{L}}_1^l(\underline{C}_{1,l})$ as follows:

$\bullet$($\mathcal{A}_1$): If $\bar{\mathcal{L}}_k^{N_0}(\bar{C}_{k,N_0})$ is not a general root, then we construct

$\bullet$($\mathcal{A}_2$): If $\underline{\mathcal{L}}_1^l(\underline{C}_{1,l})$ is not a general root, then we construct

Step 2. According to the strong connectivity of $\mathcal{G}_0$, we immediately know that $\bar{\mathcal{L}}_1^{N_0}(\bar{C}_{1,N_0})$ is a general root, and $\bar{\mathcal{L}}_k^{N_0}(\bar{C}_{k,N_0})$ is not a general root for $k >1$. Therefore, we may start from $\theta^0_{N_0}$ and follow the process $\mathcal{A}_1$ to construct $\bar{\mathcal{L}}_k^{N_0}(\bar{C}_{k,N_0})$ until $k = 1$.

Step 3. Similarly, we know that $\underline{\mathcal{L}}_1^{N_0}(\underline{C}_{1,N_0})$ is a general root and $\underline{\mathcal{L}}_1^l(\underline{C}_{1,l})$ is not a general root for $l < N_0$. Therefore, we may start from $\theta^0_1$ and follow the process $\mathcal{A}_2$ until $l = N_0$.

It is worth emphasizing that the order of the oscillators may change along time $t$, but at each time $t$, the above algorithm does work. For convenience, the algorithm from Step 1 to Step 3 will be referred to as Algorithm $\mathcal{A}$. Then, based on Algorithm $\mathcal{A}$, we will provide a priori estimates on a monotone property about the sine function, which will be crucially used later.

Lemma 3.1. Let $\{\theta^0_i\}$ be a solution to system (1) with strongly connected network $\mathcal{G}_0$. Moreover at time $t$, we also assume that the oscillators are well-ordered as (6), and the phase diameter and quantity $\eta$ satisfiy the following conditions:

where $\zeta, \gamma$ are given in the condition (3).Then at time $t$, the following assertions hold

Proof. For the proof of this lemma, please see [45] for details.

Recall the strongly connected ensemble $\mathcal{G}_0$, and denote by $\theta_i^0\ (i = 1,2,\ldots,N_0)$ the members in $\mathcal{G}_0$. For the oscillators in $\mathcal{G}_0$, based on a priori estimates in Lemma 3.1, we will design proper coefficients of convex combination which helps us to capture the dissipation structure. Now we assume that at time $t$, the oscillators in $\mathcal{G}_0$ are well-ordered as follows,

Then we apply the process $\mathcal{A}_1$ from $\theta^0_{N_0}$ to $\theta^0_1$ and the process $\mathcal{A}_2$ from $\theta^0_1$ to $\theta^0_{N_0}$ to respectively construct

where $N_0$ is the cardinality of $\mathcal{G}_0$ and $\eta$ is given in the condition (3). By induction criteria, we can deduce explict expressions about the constructed coefficients:

Note that $\bar{a}^0_{N_0+1-i} = \underline{a}^0_i, \ i = 1,2\ldots,N_0$. And we set

We define a non-negative quantity $Q^0 = \bar{\theta}_0 - \underline{\theta}_0$ where $\bar{\theta}_0 = \bar{\theta}^0_1$ and $\underline{\theta}_0 = \underline{\theta}^0_{N_0}$. Note that $Q^0(t)$ is Lipschitz continuous with respect to $t$. We then establish the comparison relation between $Q^0$ and the phase diameter $D_0(\theta)$ of $\mathcal{G}_0$ in the following lemma.

Lemma 3.2. Let $\{\theta^0_i\}$ be a solution to system (1) with strongly connected digraph $\mathcal{G}_0$. Assume that for the group $\mathcal{G}_0$, the coefficients $\bar{a}_k^0$'s and $\underline{a}_k^0$'s satisfy the scheme (7). Then at each time $t$, we have the following relation

where $\eta$ satisfies the condition (3).

Proof. As we choose the same design for coefficients of convex combination as that in [45], the proof of this lemma is same as that in [45], please see [45] for details.

From Lemma 3.2, we see that the quantity $Q^0$ can control the phase diameter $D_0(\theta)$, which play a key role in analysing the bound of phase diameter. Based on Algorithm $\mathcal{A}$ and Lemma 3.1, we first study the dynamics of the constructed quantity $Q^0$.

Lemma 3.3. Let $\{\theta^0_i\}$ be the solution to system (1) with strongly connected digraph $\mathcal{G}_0$. Moreover, for a given positive constant $D^\infty < \min\left\{\frac{\pi}{2},\zeta\right\}$, assume the following conditions hold,

where $\zeta, \gamma$ are constants, $\beta$ is given in Lemma 3.2 and$c = \frac{\left(\sum_{j = 1}^{N_0-1} \eta^j A(2N_0,j) + 1\right)\gamma}{\sin \gamma}.$Then the phase diameter of the graph $\mathcal{G}_0$ is uniformly bounded by $\gamma$:

and the dynamics of $Q^0(t)$ is controlled by the following differential inequality

Proof. The proof is similar to [45] under the assumption that the frustration $\alpha$ is sufficiently small. However, due to the presence of frustration, there are some slight differences in the process of analysis, thus we put the detailed proof in the Appendix A.

Lemma 3.3 states that the phase diameter of the digraph $\mathcal{G}_0$ remains less than $\pi$ and provides the dynamics of $Q^0$. We next exploit the dynamics of $Q^0$ and find some finite time such that all oscillators in $\mathcal{G}_0$ will be trapped into a region of quarter circle after the time.

Lemma 3.4. Let $\{\theta^0_i\}$ be a solution to system (1) with strongly connected digraph $\mathcal{G}_0$, and suppose the assumptions in Lemma 3.3 hold. Then there exists time $t_0\ge0$ such that

where $t_0$ can be estimated as below and bounded by $\bar{t}$ given in Lemma 2.8

Remark 4. According to Lemma 2.8, we see that $D_0(\theta(t)) < \zeta$ for $t \in [0,t_0)$ since $t_0 < \bar{t}$.

Proof. From Lemma 3.3, we see that the dynamics of quantity $Q^0(t)$ is governed by the following inequality

We next show that there exists some time $t_0$ such that the quantity $Q^0$ in (11) is uniformly bounded after $t_0$. There are two cases we consider separately.

$\diamond$ Case 1. We first consider the case that $Q^0(0) > \beta D^\infty$. When $Q^0(t) \in [\beta D^\infty,Q^0(0)]$, from (66), we have

That is to say, when $Q^0(t)$ is located in the interval $[\beta D^\infty,Q^0(0)]$, $Q^0(t)$ will keep decreasing with a rate bounded by a uniform slope. Then we can define a stopping time $t_0$ as follows,

And based on the definition of $t_0$, we see that $Q^0$ will decrease before $t_0$ and has the following property at $t_0$,

Moerover, from (12), it is easy to see that the stopping time $t_0$ satisfies the following upper bound estimate,

Now we study the upper bound of $Q^0$ on $[t_0, +\infty)$. In fact, we can apply (12), (13) and the same arguments in (64) to derive

$\diamond$ Case 2. For another case that $Q^0(0) \le \beta D^\infty$. We can apply the similar analysis in (64) to obtain

Then in this case, we directly set $t_0 = 0$.

Therefore, from (15), (16), and Lemma 3.2, we derive the upper bound of $D_0(\theta)$ on $[t_0, +\infty)$ as below

On the other hand, in order to verify (10), we further study $t_0$ in (17). Combining (14) in Case 1 and $t_0 = 0$ in Case 2, we see that

Here, we use the truth that $Q^0(0) \le D_0(\theta(0)) < \zeta$. Then from the assumption about $K$ in (9), i.e.,

it yields the following estimate about $t_0$,

where in this special strongly connected case, it's clear that $N_0 = N$ and $D_0(\theta) = D(\theta)$ in Lemma 2.8.

Thus, combining (17), (18) and (19), we derive the desired results.

4.   General network

In this section, we investigate the general network with a spanning tree structure, and prove our main result Theorem 1.1, which state that synchronization will emerge for Kuramoto model with frustration. According to Definition 2.5 and Lemma 2.6, we see that the digraph $\mathcal{G}$ associated to system (1) has a unique maximum node if it contains a spanning tree structure. And From Remark 3, without loss of generality, we assume $\mathcal{G}$ is decomposed into a union as $\mathcal{G} = \bigcup_{i = 0}^d \mathcal{G}_i$, where $\mathcal{G}_p$ is a maximum node of $\mathcal{G}\setminus (\bigcup_{i = 0}^{p-1} \mathcal{G}_i)$.

We have studied the situation $d = 0$ in Section 3, and we showed that the phase diameter of the digraph $\mathcal{G}_0$ is uniformly bounded by a value less than $\frac{\pi}{2}$ after some finite time, i.e., the oscillators of $\mathcal{G}_0$ will concentrate into a region of quarter circle. However, for the case that $d > 0$, $\mathcal{G}_k$'s are not maximum nodes in $\mathcal{G}$ for $k \ge 1$. Hence, the methods in Lemma 3.3 and Lemma 3.4 can not be directly exploited for the situation $d>0$. More precisely, the oscillators in $\mathcal{G}_i$ with $i < k$ perform as an attraction source and affect the agents in $\mathcal{G}_k$. Thus when we study the behavior of agents in $\mathcal{G}_k$, the information from $\mathcal{G}_i$ with $i < k$ can not be ignored.

From Remark 3 and node decomposition, the graph $\mathcal{G}$ can be represented as

and we denote the oscillators in $\mathcal{G}_k$ by $\theta^k_i$ with $1 \le i \le N_k$. Then we assume that at time $t$, the oscillators in each $\mathcal{G}_k$ are well-ordered as below:

For each subdigraph $\mathcal{G}_k$ with $k \ge 0$ which is strongly connected, we follow the process in Algorithm $\mathcal{A}_1$ and $\mathcal{A}_2$ to construct $\bar{\mathcal{L}}_{l-1}^{N_k}(\bar{C}_{l-1,N_k})$ and $\underline{\mathcal{L}}_1^{l+1}(\underline{C}_{1,l+1})$ by redesigning the coefficients $\bar{a}^k_l$ and $\underline{a}_l^k$ of convex combination as below:

By induction principle, we deduce that

Note that $\bar{a}^k_{N_k+1-i} = \underline{a}^k_i, \ i = 1,2\ldots,N_k$. By simple calculation, we have

And we further introduce the following notations,

Due to the analyticity of the solution, $Q^k(t)$ is Lipschitz continuous. Similar to Section 3, we will first establish the comparison between the quantity $Q^k(t)$ and phase diameter $D_k(\theta(t))$ of $\bigcup_{i = 0}^k \mathcal{G}_i$, which plays a crucial role in later analysis.

Lemma 4.1. Let $\theta = \{\theta_i\}$ be a solution to system (1), and assume that the network contains a spanning tree and for each subdigraph $\mathcal{G}_k$, the coefficients $\bar{a}^k_l$ and $\underline{a}^k_l$ of convex combination in Algorithm $\mathcal{A}$ satisfy the scheme (20). Then at each time $t$, we have the following relation

where $D_k(\theta) = \max\limits_{0 \le i \le k} \max\limits_{1 \le j \le N_i} \{\theta_j^i\} - \min\limits_{0 \le i \le k} \min\limits_{1 \le j \le N_i} \{\theta_j^i\}$ and $\eta$ satisfies the condition (3).

Proof. As we adopt the same construction of coefficients of convex combination in [45] which deals with the Kuramoto model without frustration on a general network, thus for the detailed proof of this lemma, please see [45].

Now we are ready to prove our main Theorem 1.1. To this end, we will follow similar arguments in Section 3 to complete the proof. Actually, we will investigate the dynamics of the constructed quantity $Q^k(t)$ that involves the influences from $\mathcal{G}_i$ with $i <k$, which yields the hypo-coercivity of the phase diameter. Applying similar arguments in Lemma 3.3 and Lemma 3.4, we have the following estimates for $\mathcal{G}_0$.

Lemma 4.2. Suppose that the network topology contains a spanning tree, and let $\theta = \{\theta_i\}$ be a solution to (1). Moreover, assume that the initial data and the quantity $\eta$ satisfy

where $\zeta, \gamma$ are positive constants. And for a given sufficiently small $D^\infty < \min\{\frac{\pi}{2},\zeta\}$, assume the frustration $\alpha$ and coupling strength $\kappa$ satisfy

where $d$ is the number of general nodes and$c = \frac{(\sum_{j = 1}^{N-1} \eta^j A(2N,j) + 1)\gamma}{\sin \gamma}.$Then the following two assertions hold for the maximum node $\mathcal{G}_0$:

1. The dynamics of $Q^0(t)$ is governed by the following equation

2. there exists time $t_0 \ge 0$ such that

where $t_0$ can be estimated as below and bounded by $\bar{t}$ given in Lemma 2.8

Next, inspiring from Lemma 4.2, we make the following reasonable ansatz for $Q^k(t)$ with $0 \le k \le d$.

Ansatz:

1. The dynamics of quantity $Q^k(t)$ in time interval $[0,\infty)$ is governed by the following differential inequality,

where we assume $D_{-1}(\theta(t)) = 0$.

2. There exists a finite time $t_k \ge 0$ such that, the phase diameter $D_k(\theta(t))$ of $\bigcup_{i = 0}^k \mathcal{G}_i$ is uniformly bounded after $t_k$, i.e.,

where $t_k$ subjects to the following estimate,

In the subsequence, we will split the proof of the ansatz into two lemmas by induction criteria. More precisely, based on the results in Lemma 4.2 as the initial step, we suppose the ansatz holds for $Q^k$ and $D_k(\theta)$ with $0 \le k \le d-1$, and then prove that the ansatz also holds for $Q^{k+1}$ and $D_{k+1}(\theta)$.

Lemma 4.3. Suppose the assumptions in Lemma 4.2 are fulfilled, and the ansatz in (27), (28) and (29) holds for some $k$ with $0 \le k \le d-1$. Then the ansatz (27) holds for $k+1$.

Proof. We will use proof by contradiction criteria to verify the ansatz for $Q^{k+1}$. To this end, define a set

From Lemma 2.8, we see that

It is clear that $\bar{t} \in \mathcal{B}_{k+1}$. Thus the set $\mathcal{B}_{k+1}$ is not empty. Define $T^* = \sup \mathcal{B}_{k+1}$. We will prove by contradiction that $T^* = +\infty$. Suppose not, i.e., $T^* < +\infty$. It is obvious that

As the solution to system (1) is analytic, in the finite time interval $[0,T^*)$, $\bar{\theta}_i$ and $\bar{\theta}_j$ either collide finite times or always stay together. Similar to the analysis in Lemma 3.3, without loss of generality, we only consider the situation that there is no pair of $\bar{\theta}_i$ and $\bar{\theta}_j$ staying together through all period $[0,T^*)$. That means the order of $\{\bar{\theta}_i\}_{i = 0}^{k+1}$ will only exchange finite times in $[0,T^*)$, so does $\{\underline{\theta}_i\}_{i = 0}^{k+1}$. Thus, we divide the time interval $[0,T^*)$ into a finite union as below

such that in each interval $J_l$, the orders of both $\{\bar{\theta}_i\}_{i = 0}^{k+1}$ and $\{\underline{\theta}_i\}_{i = 0}^{k+1}$ are preseved, and the order of oscillators in each subdigraph $\mathcal{G}_i$ with $0 \le i \le k+1$ does not change. In the following, we will show the contradiction via two steps.

$\star$ Step 1. In this step, we first verify the Ansatz (27) holds for $Q^{k+1}$ on $[0, T^*)$, i.e.,

As the proof is slightly different from that in [45] and rather lengthy, we put the detailed proof in Appendix B.

$\star$ Step 2. In this step, we will study the upper bound of $Q^{k+1}$ in (31) in time interval $[t_k,T^*)$, where $t_k$ is given in Ansatz (28) for $D_k(\theta)$. For the purposes of discussion, we rewrite the equation (31) in the interval $[0,T^*)$ as below

where $c$ is expressed by the following equation

For the term $D_k(\theta)$ in (32), under the assumption of induction criteria, the Ansatz (28) holds for $D_k(\theta)$, i.e., there exists time $t_k$ such that

And from the condition (26), it is obvious that

This directly yields that

where $0\le k\le d-1, \beta< 1, c>1.$ Then for the purposes of analysing the last two terms in the bracket of (32), we add the esimates in (34) and (35) to get

From Lemma 2.8, we have $t_k < \bar{t} \le T^*$, thus it makes sense when we consider the time interval $[t_k, T^*)$. Now based on the above estiamte (36), we apply the differential equation (32) and study the upper bound of $Q^{k+1}$ on $[t_k, T^*)$. We claim that

Suppose not, then there exists some $\tilde{t} \in (t_k, T^*)$ such that $Q^{k+1}(\tilde{t}) > M_{k+1}$. We construct a set

Since $Q^{k+1}(t_k) \le M_{k+1}$, the set $\mathcal{C}_{k+1}$ is not empty. Define $t^* = \sup \mathcal{C}_{k+1}$. Then it is easy to see that

From the construction of $M_{k+1}$, (36) and (38), it is clear that for $t \in (t^*, \tilde{t}]$

Wen apply the above inequality and integrate on both sides of (32) from $t^*$ to $\tilde{t}$ to get

which contradicts to the truth $Q^{k+1}(\tilde{t}) - M_{k+1} >0$. Thus we complete the proof of (37).

$\star$ Step 3. In this step, we will construct a contradiction to (30). From (37), Lemma 2.8 and the fact that

we directly obtain

From Lemma 4.1 and the condition (25), it yields that

Since $D_{k+1}(\theta(t))$ is continuous, we have

which obviously contradicts to the assumption $D_{k+1}(\theta(T^*)) = \gamma$ in (30).

Thus, we combine all above analysis to conclude that $T^* = +\infty$, that is to say,

Then for any finite time $T>0$, we apply (39) and repeat the analysis in Step 1 to obtain that the differential inequality (27) holds for $Q^{k+1}$ on $[0,T)$. Thus we obtain the dynamics of $Q^{k+1}$ in whole time interval $[0,+\infty)$ as below:

Therefore, we complete the proof of the Ansatz (27) for $Q^{k+1}$.

Lemma 4.4. Suppose the conditions in Lemma 4.2 are fulfilled, and the ansatz in (27), (28) and (29) holds for some $k$ with $0 \le k \le d-1$. Then the ansatz (28) and (29) holds for $k+1$.

Proof. From Lemma 4.3, we know the dynamic of $Q^{k+1}$ is governed by (40). For the purposes of discussion, we rewrite the differential equation (40) as below and discuss it on $[t_k, + \infty)$, i.e.,

where $c$ is given in (33). In the subsequence, we will apply (41) to find a finite time $t_{k+1}$ such that the quantity $Q^{k+1}$ in (41) is uniformly bounded by some value less than $\frac{\pi}{2}$ after $t_{k+1}$. We split into two cases to discuss.

$\bullet$ Case 1. We first consider the case that $Q^{k+1}(t_k) > \frac{\beta^{d-k}D^\infty}{ [4(2N+1)c]^{d-k-1}}$. In this case, When $Q^{k+1}(t) \in [\frac{\beta^{d-k}D^\infty}{ [4(2N+1)c]^{d-k-1}},Q^{k+1}(t_k)]$, we combine (36) and (41) to have

That is to say, when $Q^{k+1}(t)$ is located in the interval $[\frac{\beta^{d-k}D^\infty}{ [4(2N+1)c]^{d-k-1}},Q^{k+1}(t_k)]$, $Q^{k+1}(t)$ will keep decreasing with a rate bounded by a uniform slope. Therefore, we can define a stopping time $t_{k+1}$ as follows,

Then, based on (42) and the definition of $t_{k+1}$, we see that $Q^{k+1}$ will decrease before $t_{k+1}$ and has the following property at $t_{k+1}$,

Moreover, from (42), it yields that the stopping time $t_{k+1}$ satisfies the following upper bound estimate,

Now we study the upper bound of $Q^{k+1}$ on $[t_{k+1}, +\infty)$. In fact, we can apply (42), (43) and the same arguments in (37) to derive

On the other hand, in order to verify (29), we further study $t_{k+1}$ in (44). For the first part on the right-hand side of (44), from Lemma 2.8 and the fact that

we have the following estimates

where the denominator on the right-hand side of above inequality is positive from the conditions about $K$ and $\alpha$ in (26). For the term $t_k$ in (44), based on the assumption (29) for $t_k$, we have

Thus it yields from (44), (46) and (47) that the time $t_{k+1}$ satisfies

Moreover, from (26), it is easy to see that the coupling strength $K$ satisfies the following inequality

Thus we combine (48) and (49) to verify the Ansatz (29) for $k+1$ in the first case, i.e., the time $t_{k+1}$ subjects to the following estimate,

$\bullet$ Case 2. For another case that $Q^{k+1}(t_k) \le \frac{\beta^{d-k}D^\infty}{ [4(2N+1)c]^{d-k-1}}$. Similar to the analysis in (37), we apply (42) to conclude that

In this case, we directly set $t_{k+1} = t_k$. Then, from (47), it yields that the inequalities (48) and (50) also hold, which finish the verification of the Ansatz (29) in the second case.

Finally, we are ready to verify the ansatz (28) for $k+1$. Actually, we can apply (45), (51) and Lemma 4.1 to have the upper bound of $D_{k+1}(\theta)$ on $[t_{k+1}, +\infty)$ as below

Then we combine (48), (50) and (52) in Case 1 and similar analysis in Case 2 to conclude that the Ansatz (28) and (29) is true for $D^{k+1}(\theta)$.

Now, we are ready to prove our main result.

Lemma 4.5. Let $\theta = (\theta_1, \theta_2, \ldots, \theta_N)$ be a solution to system (1), and suppose that the network contains a spanning tree and the assumptions in Lemma 4.2 are fulfilled. Then there exists a finite time $t_*\ge0$ such that

where $t_* < \bar{t}$ and $\bar{t}$ is given in Lemma 2.8.

Proof. Combining Lemma 4.2, Lemma 4.3 and Lemma 4.4, we apply inductive criteria to conclude that the Ansatz (27) –(29) hold for all $0\leq k\leq d$. Then, it yields from (28) and (29) that there exists a finite time $0 \le t_d < \bar{t}$ such that

Thus we set $t_* = t_d$ and derive the desired result.

Remark 5. For the Kuramoto model with frustration, in Lemma 4.5, we show the phase diameter of whole ensemble will be uniformly bounded by a value $D^\infty$ less than $\frac{\pi}{2}$ after some finite time. Under the assumption that $\alpha$ is sufficiently small such that $D^\infty + \alpha < \frac{\pi}{2}$, the interaction function $\cos x$ in the dynamics of frequency is positive after the finite time. Thus, we can lift system (1) to the second-order formulation, which enjoys the similar form to Cucker-Smale model with the interaction function $\cos x$.

More precisely, we can introduce phase velocity or frequency $\omega_i(t) : = \dot{\theta}_i(t)$ for each oscillator, and directly differentiate (1) with respect to time $t$ to derive the equivalent second-order Cucker-Smale type model as below

Now for the second-order system (53), we apply the results in [9] for Kuramoto model without frustration on a general digraph and present the frequency synchronization for Kuramoto model with frustrations.

Lemma 4.6. Let $(\theta(t),\omega(t))$ be a solution to system (53), and suppose the network contains a spanning tree and the assumptions in Lemma 4.2 are fulfilled.Then there exist positive constants $C_1$ and $C_2$ such that

where $D(\omega(t)) = \max_{1\le i \le N} \{\omega_i(t)\} - \min_{1\le i \le N} \{\omega_i(t)\}$ is the diameter of frequency.

Proof. According to Lemma 4.5 and the condition $D^\infty + \alpha < \frac{\pi}{2}$ in (26), we see that there exists time $t_* \ge 0$ such that

Therefore, we can apply the methods and results in the work of Dong et al. [9] for Kuramoto model without frustration to yield the emergence of exponentially fast synchronization. As the proof is almost the same as that in [9], we omit its details.

Proof of Theorem 1.1: Combining Lemma 4.5 and Lemma 4.6, we ultimately derive the desired result in Theorem 1.1.

5.   Numerical simulations

In this section, we present several numerical simulations to illustrate the main results in Theorem 1.1, which state that in sufficiently large coupling strength and small frustration regimes, synchronous behavior will emerge for the Kuramoto model with frustration in half circle case.

For the simulation, we use the fourth-order Runge-Kutta method and employ the parameter $N = 10$. The natural frequencies $\Omega_i$ are randomly chosen from the interval $(-1,1)$, and the initial configuration $\theta_{i0}$ are randomly chosen from the interval $(0,\frac{5\pi}{6})$ which are confined in half circle. The interaction network we exploit are presented in Figure 1, which contains a spanning tree structure. And the corresponding $(0,1)$-adjacency matrix $(\chi_{ij})$ is given as below:

For the fixed digraph in Figure 1, we choose large $K = 20$ and small $\alpha = 0.01$ in Figure 2. Figure 2a shows that the phase diameter is uniformly bounded by a value less than $\frac{\pi}{2}$ after some finite time, and we observe that the frequency diameter $D(\omega(t))$ decays to zero exponentially fast in Figure 2b. This is consistent with the analytical results in Theorem 1.1.

In Figure 3, we fix the coupling strength $K = 1$, and respectively employ $\alpha = 0, 0.01, 0.1, 0.3$ to observe the effects of frustration. Figure 3a–3d shows the dynamics of frequencies and we observe that in small frustration regime (i.e., $\alpha = 0, 0.01$), the asymptotic behavior evolves to frequency synchronization, while the state of desynchronization occurs for large frustration value (i.e., $\alpha = 0.1, 0.3$). Figure 3e displays the dynamics of frequency diameter for two cases $\alpha = 0.01, 0.3$.

In Figure 4, we fix the frustration $\alpha = 0.1$ and respectively choose $K = 0.8,1,1.2,1.4$ for simulation. Figure 4a–4d presents that the state of desynchronization occurs when $K$ is small (i.e., $K = 0.8,1$), while frequency synchronization asymptotically emerges for large $K$ (i.e., $K = 1.2,1.4$). For comparison, the dynamics of frequency diameter when $K = 0.8, 1.4$ is displayed in Figure 4e.

Moreover, for $K = 1$, we observe in Figure 3a that frequency synchronization emerges asymptotically with zero frustration (i.e., $\alpha = 0$), whereas desynchronization occurs when $\alpha = 0.1$ in Figure 3c. This implies that the frustration hinders synchronization. And in Figure 4c–4d, synchronization emerges when $K = 1.2, 1.4$ for $\alpha = 0.1$. This indicates that a larger coupling strength $K$ is needed to guarantee the emergence of synchronization than that in zero frustration case.

6.   Summary

In this paper, under the effect of frustration, we provide sufficient frameworks leading to the complete synchronization for the Kuramoto model with general network containing a spanning tree. To this end, we follow a node decomposition introduced in [25] and construct hypo-coercive inequalities through which we can study the upper bounds of phase diameters. When the initial configuration is confined in a half circle, for sufficiently small frustration and sufficiently large coupling strength, we show that the relative differences of Kuramoto oscillators adding a phase shift will be confined into a region of quarter circle in finite time, thus we can directly apply the methods and results in [9] to prove that the complete frequency synchronization emerges exponentially fast. And we provide some numerical simulations to illustrate the main results.

Acknowledgments

We really appreciate the editors and reviewers for their thorough reviews and insightful suggestions.

Appendix A.   Proof of Lemma 3.3

We will split the proof into six steps. In the first step, we suppose by contrary that the phase diameter of $\mathcal{G}_0$ is bounded by $\gamma$ in a finite time interval. In the second, third and forth steps, we use induction criteria to construct the differential inequality of $Q^0(t)$ in the finite time interval. In the last two steps, we exploit the derived differential inequality of $Q^0(t)$ to conclude that phase diameter of $\mathcal{G}_0$ is bounded by $\gamma$ on $[0, +\infty)$, and thus the differential inequality of $Q^0(t)$ obtained in the forth step also holds on $[0, +\infty)$.

$\bullet$ Step 1. Define a set

From Lemma 2.8 where $N = N_0$ in the present section, the set $\mathcal{B}_0$ is non-empty since

which directly yields that $\bar{t} \in \mathcal{B}_0$. Define $T^* = \sup \mathcal{B}_0$. And we claim that $T^* = +\infty$. Suppose not, i.e., $T^* < +\infty$, then we apply the continuity of $D_0(\theta(t))$ to have

In particular, we have $\bar{t} \le T^*$. The analyticity of the solution to system (1) is guaranteed by the standard Cauchy-Lipschitz theory. Therefore, in the finite time interval $[0,T^*)$, any two oscillators either collide finite times or always stay together. If there are some $\theta_i$ and $\theta_j$ always staying together in $[0,T^*]$, we can view them as one oscillator and thus the total number of oscillators that we need to study can be reduced. This is a more simpler situation, and we can similarly deal with it. Therefore, we only consider the case that there is no pair of oscillators staying together in $[0,T^*)$. For this case, there are only finite many collisions occurring through $[0,T^*)$. Thus, we divide the time interval $[0,T^*)$ into a finite union as below

where the end point $t_l$ denotes the collision instant. It is easy to see that there is no collision in the interior of $J_l$. Now we pick out any time interval $J_l$ and assume that

$\bullet$ Step 2. According to the notations in (8), we follow the process $\mathcal{A}_1$ and $\mathcal{A}_2$ to construct $\bar{\theta}^0_n$ and $\underline{\theta}^0_n, \ 1\le n\le N_0$, respecively. We first study the dynamics of $\bar{\theta}_{N_0}^0 = \theta_{N_0}^0$, i.e.,

For the dynamics of $\bar{\theta}^0_{N_0-1}$, according to the process $\mathcal{A}_1$ and $\bar{a}^0_{N_0-1} = \eta(N_0 + 2)$ in (7), we apply (55) and estimate the derivative of $\bar{\theta}^0_{N_0-1}$ as follows,

where we use

Next we show the term $\mathcal{I}_1$ is non-positive. We only consider the situation $\gamma > \frac{\pi}{2}$, and the case $\gamma \le \frac{\pi}{2}$ can be similarly dealt with. It is clear that

Note that $\bar{k}_{N_0} \le N_0-1$ since $\bar{\mathcal{L}}^{N_0}_{N_0}(\bar{C}_{N_0,N_0})$ is not a general root. Therefore, if $0 \le \theta^0_{N_0}(t) - \theta^0_{\bar{k}_{N_0}}(t) \le \frac{\pi}{2}$, we immediately obtain that

which implies that

On the other hand, if $\frac{\pi}{2} < \theta^0_{N_0}(t) - \theta^0_{\bar{k}_{N_0}}(t) < \gamma$, we use the fact

to conclude that $\eta \sin(\theta^0_{\bar{k}_{N_0}} - \theta^0_{N_0}) \le -1$. Hence, in this case, we still obtain that

Thus, for $t\in J_l$, we combine above analysis to conclude that

Then combining (56) and (57), we derive that

$\ $

$\bullet$ Step 3. Now we apply the induction principle to cope with $\bar{\theta}^0_n$ in (8), which is construced in the iteration process $\mathcal{A}_1$. We will prove for $1 \le n \le N_0$ that,

In fact, it is known that (59) already holds for $n = N_0, N_0-1$ from (55) and (58). Then by induction criteria, suppose (59) holds for $n$. Next we verify that (59) still holds for $n-1$. According to the Algorithm $\mathcal{A}_1$ and similar calculations in (56), the dynamics of the quantity $\bar{\theta}^0_{n-1}(t)$ subjects to the following estimates

Moreover, we can prove the term $\mathcal{I}_2$ is non-positive. As the proof is very similar as that in the previous step, we omit the details and directly claim that $\mathcal{I}_2 \le 0$, which together with (60) verifies (59).

$\bullet$ Step 4. Now, we set $n = 1$ in (59) and apply Lemma 3.1 to have

where $\bar{k}_{1} = \min_{j \in \bigcup_{i = 1}^{N_0} \mathcal{N}^0_{i}(0)} j = 1$ due to the strong connectivity of $\mathcal{G}_0$. Similarly, we can follow the process $\mathcal{A}_2$ to construct $\underline{\theta}^0_k$ in (8) until $k = N_0$. Then, we can apply the similar argument in (59) to obtain that,

Then we recall the notations $\bar{\theta}_0 = \bar{\theta}^0_1$ and $\underline{\theta}_0 = \underline{\theta}^0_{N_0}$, and combine (61) and (62) to obtain that

where we use the property

As the function $\frac{\sin x}{x}$ is monotonically decreasing in $(0, \pi]$, we apply (54) to obtain that

Moreover, due to the fact $Q^0(t) \le \theta^0_{N_0}(t) - \theta^0_1(t)$, we have

Note that the constructed quantity $Q^0(t) = \bar{\theta}_0(t) - \underline{\theta}_0(t)$ is Lipschitz continuous on $[0,T^*)$. Moreover, the above analysis does not depend on the time interval $J_l, \ l = 1,2,\ldots,r$, thus the differential inequality (63) holds almost everywhere on $[0,T^*)$.

$\bullet$ Step 5. Next we study the upper bound of $Q^0(t)$ in the period $[0,T^*)$. Define

We claim that

Suppose not, then there exists some $\tilde{t} \in [0,T^*)$ such that $Q^0(\tilde{t}) > M_0$. We construct a set

Since $0 \in \mathcal{C}_0$, the set $\mathcal{C}_0$ is not empty. Then we denote $t^* = \sup \mathcal{C}_0$. It is easy to see that

For the given constant $D^\infty < \min\{\frac{\pi}{2},\zeta\}$, based on the assumptions about the frustration and the coupling strength in (9), it is clear that

where $c = \frac{\left(\sum_{j = 1}^{N_0-1} \eta^j A(2N_0,j) + 1\right)\gamma}{\sin \gamma}.$ Thus combing the construction of $M_0$, (65) and (66), we obtain that for $t \in (t^*, \tilde{t}]$, the following estimate holds,

Then, we apply the above inequality and integrate on the both sides of (63) from $t^*$ to $\tilde{t}$ to get

which obviously contradicts to the fact $Q^0(\tilde{t}) - M_0 >0$, and verifies (64).

$\bullet$ Step 6. Now we are ready to show the contradiction to (54), which implies that $T^* = +\infty$. In fact, from the fact that $\beta < 1, D^\infty < \zeta$ and $Q^0(0) \le D_0(\theta(0)) < \zeta$, we see

Then we apply the relation $\beta D_0(\theta(t)) \le Q^0(t)$ given in Lemma 3.2 and the assumption $\eta > \frac{2}{1 - \frac{\zeta}{\gamma}}$ in (9) to obtain that

As $D_0(\theta(t))$ is continuous, we have

which contradicts to the situation that $D_0(\theta(T^*)) = \gamma$ in (54). Therefore, we conclude that $T^* = +\infty$, which implies that

Then for any finite time $T>0$, we apply (67) and repeat the same argument in the second, third, forth steps to obtain the dynamics of $Q^0(t)$ in (63) holds on $[0,T)$. Thus we obtain the following differential inequality of $Q^0$ on the whole time interval:

Appendix B.   Proof of Step 1 in Lemma 4.3

We will show the detailed proof of Step 1 in Lemma 4.3. Now we pick out any interval $J_l$ with $1\le l \le r$, where the orders of both $\{\bar{\theta}_i\}_{i = 0}^{k+1}$ and $\{\underline{\theta}_i\}_{i = 0}^{k+1}$ are preseved and the order of oscillators in each subdigraph $\mathcal{G}_i$ with $0\le i \le k+1$ will not change in each time interval. Then, we consider four cases depending on the possibility of relative position between $\bigcup_{i = 0}^k \mathcal{G}$ and $\mathcal{G}_{k+1}$.

Figure 5 shows the four possible relations between $\bigcup_{i = 0}^k \mathcal{G}$ and $\mathcal{G}_{k+1}$ at any time $t$. Case $1$ and Case $4$ are similar and relative simple, while the analysis on Case $2$ and Case $3$ are similar but much more complicated. Therefore, we will only show the detailed proof of Case $2$ for simplicity. In this case, we have from Figure 5 that

Without loss of generality, we assume that

$\bullet$ Step 1. Similar to (59), we claim that for $1 \le n \le N_{k+1}$, the following inequalities hold

where $S_k = \sum_{i = 0}^k N_i$. In the following, we will prove the claim (68) via induction principle.

$\bullet$ Step 1.1. As an initial step, we first verify that (68) holds for $n = N_{k+1}$. In fact, we have

where we use

$\diamond$ Estimates on $\mathbf{\mathcal{I}_{11}}$ in (69). We know that $\theta^{k+1} _{N_{k+1}}$ is the largest phase among $\mathcal{G}_{k+1}$, and all the oscillators in $\bigcup_{i = 0}^{k+1} \mathcal{G}_{i}$ are confined in half circle before $T^*$. Therefore, it is clear that

Then we immediately have

$\diamond$ Estimates on $\mathbf{\mathcal{I}_{12}}$ in (69). For $\theta^l_j$ which is the neighbor of $\theta^{k+1}_{N_{k+1}}$ in $\mathcal{G}_l$ with $0 \le l \le k$, i.e., $j\in \mathcal{N}_{N_{k+1}}^{k+1}(l)$, we consider two possible orderings between $\theta^l_j$ and $\theta^{k+1}_{N_{k+1}}$:

If $\theta^l_j \le \theta^{k+1}_{N_{k+1}}$, we immediately have

If $\theta^l_j > \theta^{k+1}_{N_{k+1}}$, from the fact that

we immediately obtain

Thus we use the property of $\sin x \le x, \ x \ge 0$ and (71) to get

Therefore, combining the above discussion, we have

From (69), (70) and (72), it yields that (68) holds for $n = N_{k+1}$.

$\bullet$ Step 1.2. Next, we assume that (68) holds for $n$ with $2 \le n \le N_{k+1}$, and we will show that (68) holds for $n-1$. Following the process $\mathcal{A}_1$ and similar analysis in (69), we have

Next we do some estimates about the terms $\mathcal{I}_{21}$ and $\mathcal{I}_{22}$ in (73) seperately.

$\diamond$ Estimates on $\mathbf{\mathcal{I}_{21}}$ in (73). Without loss of generality, we only deal with $\mathcal{I}_{21}$ under the situation $\gamma > \frac{\pi}{2}$. We first apply the strong connectivity of $\mathcal{G}_{k+1}$ and Lemma 3.1 to obtain that

where $\bar{k}_{n} = \min_{j \in \bigcup_{i = n}^{N_{k+1}} \mathcal{N}^{k+1}_{i}(k+1)} j \le n-1$. According to (74), we consider two cases depending on comparison between $\theta^{k+1}_{N_{k+1}} - \theta^{k+1}_{\bar{k}_{n}}$ and $\frac{\pi}{2}$.

(i) For the first case that $0 \le \theta^{k+1}_{N_{k+1}} - \theta^{k+1}_{\bar{k}_{n}} \le \frac{\pi}{2}$, we immediately obtain that for $j \in \mathcal{N}^{k+1}_{n-1}(k+1), \ j>n-1$,

Then it yieldst from (74), (75) and $\eta > 2$ that

(ii) For the second case that $\frac{\pi}{2} < \theta^{k+1}_{N_{k+1}} - \theta^{k+1}_{\bar{k}_{n}} < \gamma$, it is known that

which yields $\eta \sin(\theta^{k+1}_{\bar{k}_{n}} - \theta^{k+1}_{N_{k+1}} ) \le -1$. Thus we immediately derive that

Then, we combine the above arguments in (i) and (ii) to obtain

$\diamond$ Estimates on $\mathbf{\mathcal{I}_{22}}$ in (73). For the term $\mathcal{I}_{22}$, there are three possible relations between $\theta^{k+1}_{n-1}$ and $\theta^l_j$ with $0 \le l \le k$:

(i) If $\theta^l_j \le \theta^{k+1}_{n-1}$, we immediately have $\sin (\theta^{l}_j - \theta^{k+1}_{n-1}) \le 0$.

(ii) If $\theta^{k+1}_{n-1} < \theta^l_j \le \theta^{k+1}_{N_{k+1}}$, we consider two cases separately:

(a) For the case that $0 \le \theta^{k+1}_{N_{k+1}} - \theta^{k+1}_{\bar{k}_{n}} \le \frac{\pi}{2}$, it is clear that

Thus from the above inequality and (74), we have

(b) For another case that $\frac{\pi}{2} < \theta^{k+1}_{N_{k+1}} - \theta^{k+1}_{\bar{k}_{n}} < \gamma$, it is yields from (76) that

Hence, combining the above arguments in (a) and (b), we obtain that

$\ $

(iii) If $\theta^l_j > \theta^{k+1}_{N_{k+1}}$, we exploit the concave property of sine function in $[0, \pi]$ to get

For the second part on the right-hand side of above inequality (78), we apply the same analysis in (ii) to obtain

For the first part on the right-hand side of (78), the calculation is the same as (72), thus we have

Therefore, we combine the above estimates to obtain

Then from (73), (77), and (79), it yields that

This means that the claim (68) does hold for $n-1$. Therefore, we apply the inductive criteria to complete the proof of the claim (68).

$\bullet$ Step 2. Now we are ready to prove (27) on $J_l$ for Case 2. In fact, we apply Lemma 3.1 and the strong connectivity of $\mathcal{G}_{k+1}$ to have

From the notations in (22) and (23), it is known that

Thus, we exploit the above inequality and set $n = 1$ in (68) to obtain

We further apply the similar arguments in (80) to derive the differential inequality of $\underline{\theta}_{k+1}$ as below

Due to the monotone decreasing property of $\frac{\sin x}{x}$ in $(0, \pi]$ and from (30), it is obvious that

Then we combine the above inequality, (80), (81) and (21) to get

where we use the fact that $Q^{k+1}(t) \le \theta^{k+1}_{N_{k+1}}(t) - \theta^{k+1}_1(t)$ and (21). Eventually, for Case 2, we obtain the dynamics for $Q^{k+1}(t)$ in (27) on $J_l$, i.e.,

As our analysis does not depend on the choice of $\mathcal{J}_l$, the differential inequality (82) holds on the interval $[0,T^*)$.