Dynamics of a stochastic SIRS epidemic model with standard incidence and vaccination

1.   Introduction

The theoretical study of weakly coupled limit cycle oscillators is being actively developed in several research fields. For example, in statistical physics, various models are being developed, whereas in network science, synchronization on complex networks is attracting attention. As we state later, mathematical analysis of this field is being promoted, especially by those who are concerned with functional equations. Furthermore, this phenomenon is applied to various areas of engineering, including neural networks, bio-sciences, and network engineering[9].

Limit cycle oscillators, which are also called nonlinear oscillators, are different from harmonic oscillators, whose limit cycle is vulnerable to the perturbation forces from outside. The synchronization phenomenon is another feature of coupled limit cycle oscillators under specific conditions.

As is well known, such a phenomenon was first discovered by Huygens in the 17th century, who devised the pendulum clock for navigation officers. The physical formulation of this phenomenon was rigorously discussed later in the 1960's. In 1967, Winfree[37] proposed an attractive formulation, in which he successfully and rigorously defined the phase of the oscillator. He also revealed the region of synchronization in the sense of the phase space.

Based on Winfree's contribution, Kuramoto[19] proposed a simplified but more sophisticated model, which is called the Kuramoto model. The features of his approach are phase reduction and mean field approximation. His discussion begins with the usual dynamical system:

where X is the $n$-dimensional vector, and ${\bf{F}}$ is the $n$-dimensional vector valued function. Based on Winfree's theory, Kuramoto defined the phase of this system as that satisfying

where $\omega = 2\pi/T_0$ is the frequency of the oscillator with period $T_0$. From this definition, he derived the temporal behavior of the phase as

which is the start point of phase reduction. In an $N$-oscillator system, he described the temporal evolution of the phase (or disturbance of phase) of each oscillator:

where $\phi_j \; (j=1,2,\ldots,N)$ is the phase or phase disturbance of the $j$-th oscillator, $\omega_j$, the natural frequency (for the definition of the natural frequency, see p.67 in[19], for instance), and $K_{jk}$, the coupling strength between the $j$th and $k$th oscillators.

Kuramoto[17][19] further sophisticated (1.1) by applying the mean field approximation and the assumption $K_{jk}=K/N$, which means that all oscillators couple with uniform strength (Kuramoto called this model as global coupling[20]):

where $r=r(t)$ is the order parameter, and $\eta$ is the average phase. Note that $i=\sqrt{-1}$ hereafter. As is well known, $r \in [0,1]$ measures the coherence strength of the oscillators.

Since this model is sufficiently simple for rigorously analyzing and simulating on computers, numerous investigations have been conducted on it. For instance, Daido[8] derived a model that replaces the term $\sin \theta$ by a more general function, and Bonilla[4] revealed stability with the aid of asymptotic expansion. Recently, Li[26] considered the Kuramoto-type model with intrinsic frustrations.

As for the stability analysis, the work by Strogatz and Mirollo[34] is the first one, which concerns the spectrum of the incoherent state as the number of oscillators tends to infinity. Later, Crawford[6][7] applied the center manifold reduction to verify the stability of incoherence in detail.

Due to the limitation of space, we refer the reader to the survey by Acebrón[3] of this research area.

From the perspective of network science, it is interesting to generalize the network topology and distribution of coupling strength. Ichinomiya[14] proposed a model on random networks, and Nakao[29] numerically analyzed a model on a complex network.

Although (1.2) is a system of ordinary differential equations, adding white noise to it makes it possible to apply partial differential equation-based analysis. We consider

where $\{\xi_j(t)\}_{j=1}^N$ are the independent Wiener processes satisfying

$< \cdot>$ stands for the mean value, $D>0$, the diffusion coefficient, $\delta(\cdot)$, the Dirac's delta function, and $\delta_{jk}$, the Kronecker's delta. In virtue of the theory of stochastic processes (see Remark 1 below), the Fokker-Planck equation, which describes the temporal behavior of the probability density of particles subject to (1.3)[10][30][33], is written as

where $\varrho(\theta,t;\omega)$ is the probability distribution function of the oscillators' phase $\theta$ at time $t$ with natural frequency $\omega$, and $g(\omega)$ is the probability distribution function of $\omega$. The independent variables $\theta$ and $t$ represent the phase of oscillators and time, respectively, and we regard $\omega$ as a parameter.

Remark 1. Let $X(t)$ be a Markovian process defined on $t>0$, and $\{K_j(x)\}_{j=1}^\infty$ be a series of its intensity defined by

Then, the probability distribution function $w(x,t)$ of $X(t)$ at time $t$ satisfies[33]

A Markovian process satisfying

is called to be continuous. One example of this type of stochastic process is the one driven by the Brownian motion, like (1.3). In this case, (1.5) becomes

which is the Fokker-Planck equation. Applying this argument to (1.3), where $K_1 = V_1$ and $K_2 = 2D$, the corresponding Fokker-Planck equation is written as

which corresponds to (1.4)$_1$ (hereafter, we represent the $i$th equation of (a.b) as ${\rm (a.b)}_i$).

By substituting (1.4)$_2$ into (1.4)$_1$, we have a single integro-differential equation with suitable boundary and initial conditions:

which is the so-called Kuramoto-Sakaguchi equation[24]. This approach enables macroscopic analysis when a system consists of numerous oscillators. As we state later, (1.6) with $D=0$ has also been discussed in past arguments.

Kuramoto also presented a direction to take the spatial validity of coupling strength into account in his original model, which he called the non-local coupling model. In it, the strength of the connection between oscillators depends on the distance between them[21]. Due to the varying strength of connections, it was shown that characteristic patterns, such as a chimera pattern, emerge[21][22]. Numerical studies of the chimera state have also been conducted, mainly by Abrams and Strogatz[1][2].

In spite of numerous contributions concerning numerical simulations, there have been few studies regarding the mathematical analysis of this model. In this paper, we discuss the existence and uniqueness of the global-in-time solution to this model. We also discuss the nonlinear stability of the trivial stationary solution and existence of the vanishing diffusion limit. We note that the vanishing diffusion limit is discussed in the function spaces of higher order derivative than that by Ha and Xiao's[11] by applying a different approach.

The remainder of this paper is organized as follows. In the next section, we formulate the problem of the phase reduction of non-local coupling oscillators, including both with and without diffusion. In Section 3, we introduce related work. In Sections 4 and 5, we are concerned with the local and global-in-time solvability of the problem and the stability of the stationary incoherent solution, respectively. We discuss the existence of the solution in the vanishing diffusion limit in Section 6, and provide a conclusion and discuss remaining issues in Section 7.

2.   Formulation of problem

In this section, we formulate the problem to be considered.We begin the discussion with the temporal evolution of the phase with additive noise under non-local coupling [15][20][31]:

where $x$ and $y$ stand for the location of each oscillator; $\Gamma(\cdot)$, the phase coupling function which is periodic with respect to its argument and depend only on the difference of each oscillator's phase, and $G(\cdot)$, its strength; $\omega$, the natural frequency of each oscillator, and $\xi(x,t)$, the Wiener process satisfying:

with $D$ being the diffusion coefficient. Kuramoto mentioned that the non-local coupling model encompasses the original Kuramoto model as a specific case of (2.1) when $G(x)$ is a constant (p.132 in[20], [21]).

It is also to be noted that $\sin(\phi_j-\phi_k)$ in (1.1) is the first Fourier mode approximation of the periodic function $\Gamma(\phi_j-\phi_k)$[19]. It is possible to adopt $\sin(\cdot)$ for the non-local coupling case as in[16], but in many cases, they apply $\Gamma(\cdot)$ for the non-local coupling model. By regarding $\omega$ and $\phi$ as random variables, the average of the second term in the right-hand side is written as[15][16][18]

where $g(\cdot)$ is the probability distribution function of the natural frequency, and $\varrho=\varrho(\cdot,t;x,\omega)$ is the probability distribution function of the phase at time $t$, with location $x$ and natural frequency $\omega$ of each oscillator. We replace the second term in the right-hand side of (2.1) with the representation above to obtain:

If we consider the case that $G(x)$ is a constant $K$ in (2.1), the average of the second term in the right-hand side is

and (2.2) becomes

which corresponds to (1.3) as a infinite population limit of it.

Note that (2.2) is the mean field approximation to (2.1), in the sense that we replace the second term of the right-hand side of (2.1) with its average. However, while (1.4)$_2$ is derived by taking the infinite limit of oscillators' population, this representation of the average is theoretically exact.

By tracing the same argument as we derived (1.4), the Fokker-Planck equation corresponding to (2.2) is written as (hereafter we denote the phase by $\theta$)

Together with suitable initial and boundary conditions, the problem corresponding to (2.2) is written as[15][31].

where

It is to be noted that, in[15], [21] and[31], they deal with problem (2.3) with $F[\varrho_1,\varrho_2]$ replaced with

which corresponds to the case $g(\omega^\prime) = \delta(\omega^\prime-\omega).$ This means that all the oscillators have the same natural frequency, and $\varrho$ does not depend on $\omega$. In other words, $\omega$ is not regarded as a random variable. However, it is generally reasonable to apply (2.3)-(2.4) as Kawamura[16] did, as in the original Kuramoto model ((1.4) or (1.6)). Therefore, we employ (2.3) together with (2.4) in this paper.

In (2.3), the unknown function is $\varrho$, and the independent variables are $\theta$ and $t$, which stand for the phase and time, respectively. We regard $x$ and $\omega$ as parameters. We also use the notations

Note that we denote the $j$-th derivative of $\Gamma(\theta)$ by $\Gamma^{(j)} \; (j=1,2,\ldots)$, especially the first derivative by $\Gamma^\prime(\theta).$ Obviously, $\Gamma^{(0)}=\Gamma$.

As in the original Kuramoto model, the vanishing diffusion case is worth considering:

It is obvious that $\bar{\varrho}=1/2\pi$ is the trivial stationary solution to (1.6), (2.3) and (2.5). Moreover, on the basis of appropriate assumptions on $\varrho_0$, $g$ and $G$, the following properties of $\varrho$ will be derived (see, for instance, Lemma 2.1 in [11], Lemmas 1.1-1.2 in[24] and Lemma 4.3 of this paper):

for arbitrary $T>0$, which are natural as the properties of the probability distribution function. In this paper, we discuss the solvability and some stability properties of (2.3), as well as the convergence of the solution of (2.3) to that of (2.5) as $D$ tends to zero.

3.   Related work

Mathematical arguments concerning the solvability of the Kuramoto-Sakaguchi equation (1.6), which corresponds to the original Kuramoto model (1.3), was first presented by Lavrentiev et al.[24][25]. In their former work[24], they constructed the classical global-in-time solution when the support of $g(\omega)$ is compact.

Later, they removed this restriction[25] by applying the a-priori estimates derived from the energy method. They also studied the regularity of the unknown function with respect to $\omega$. Concenring the stability, a pioneering work was conducted by Strogatz and Mirollo[34], who were concerned with the linear stability of the trivial stationary solution $\bar{\varrho} = 1/2\pi$. Through the investigation of the spectrum of the linearized operator, they verified the existence of the critical coupling strength over which the coherent state becomes stable. Recently, Ha and Xiao[11] discussed the nonlinear stability of $\bar{\varrho}$ and convergence of the solution as $D$ tends to zero. However, their estimate remained in the space $L_\infty$ with respect to $\theta$, as we discuss in detail in Section 6. They also verified the instability of $\bar{\varrho}$ when the support of $g(\omega)$ is sufficiently narrow[12].

For the case of vanishing diffusion $D=0$ in (1.6), Chiba[5] argued the nonlinear stability of the trivial stationary solution under the assumption of unbounded support of $g(\omega)$.

Concerning the non-local coupling model, however, there are no mathematical arguments as far as we know.

4.   Global-in-time solvability

In this section, we discuss the global-in-time solvability of (2.3). We first prepare the definitions of function spaces.

4.1.   Function spaces

In this subsection, we define the function spaces used throughout this paper. Let $T>0$, and ${\mathcal G}$ be an open set in ${\bf{R}}$. Hereafter, $L_2({\mathcal G})$ stands for a set of square-integrable functions defined on ${\mathcal G}$, equipped with the norm

The inner product is defined by

where $\bar{z}$ denotes the complex conjugate of $z \in {\bf{C}}$.

For simplicity, we denote the $L_2$-norm of a function $f(\theta,t;x,\omega)$ with respect to $\theta$ merely by $\|f\|$ or $\|f(\cdot,t;x,\omega)\|.$

Hereafter, let us use the notation $\Omega \equiv (0,2\pi)$ for simplicity.

By $C({\mathcal G})$ and $C^k({\mathcal G}) \; \bigl( k \in {\bf{N}} \bigcup \{+\infty\} \bigr)$, we mean the spaces of real continuous and $k$-times continuously differentiable functions on $\Omega$, respectively. The notation $C_0^\infty({\mathcal G})$ denotes a set of $C^\infty({\mathcal G})$ functions with a compact support in ${\mathcal G}$.

For a Banach space $E$ with the norm $\|\cdot\|_E$, we denote the space of $E$-valued measurable functions $u(t)$ on the interval $(a,b)$ by $L_p(a,b;E)$, whose norm is defined by

Likewise, we denote by $C(a,b;E)$ (and by $C^k(a,b;E)$) the space of continuous functions (resp. $k$-times continuously differentiable functions) from $(a,b)$ into $E$.

Subject to the definition by Temam[35], we say, for a fixed parameter $\omega \in {\bf{R}}$, a $2\pi$-periodic function

which is expanded in the Fourier series, belongs to the Sobolev space ${\mathcal H}^m(\Omega) \; (m>0)$ when it satisfies

Due to the definition of the Fourier series, the Fourier coefficients $a_n \; (n=\pm1,\pm2,\ldots)$ of a function $u$ are defined by

Note that in case $m \in {\bf{N}}$, the norm above is equal to the usual Sobolev norm

We also introduce the notation $C_{x,\omega}^\infty({\bf{R}}^2)$ to denote a set of functions defined on ${\bf{R}}^2$, which are infinitely smoothly differentiable with respect to both $x$ and $\omega$. Let us introduce the following notations:

where $T>0$ is an arbitrary number. In addition, we use notations $|\!|\!|u |\!|\!|_m \equiv \sup_{x,\omega} \|u(\cdot;x,\omega)\|_m$ and $|\!|\!|u |\!|\!| \equiv |\!|\!|u |\!|\!|_0$ for brevity. Hereafter, $c$'s represent constants in the estimate of some quantities. When we denote $c(t)$ with suffixes, it depends on $t$. For simplicity, we hereafter use notations $f^{(j,k)} \equiv \Bigl( \frac{\partial}{\partial \theta}\Bigr)^j \Bigl( \frac{\partial}{\partial t}\Bigr)^k f \; (j,k=0,1,2,\ldots)$ for a function $f=f(\theta,t)$ in general.

4.2.   Local-in-time solution

First, we state the existence and uniqueness of the local-in-time solution to problem (2.3).

Theorem 4.1. Let us assume $m \in {\bf{N}}$, and the following issues:

(ⅰ) $\Gamma(\theta)$ is $2\pi$-periodic, and satisfies $|\Gamma^{(k)}(\theta)| \leq C_k \; (k=0,1,2,\ldots,2m)$;

(ⅱ) $\int_0^{2\pi} \Gamma^{(k)}(\theta) \; {\rm d}\theta = 0 (k=0,1,2,\ldots,2m)$;

(ⅲ)$\|\Gamma^{(k)}\| = \widetilde{C}_k < +\infty \; (k=0,1,2,\ldots,2m)$;

(ⅳ) $g(\omega) \in L_1^{(1)}$;

(ⅴ) $G(x) \in L_1({\bf{R}}) \bigcap L_2({\bf{R}}) \bigcap C^\infty({\bf{R}})$ satisfies $G \geq 0$, $\displaystyle \int_{\bf{R}} G(x) \; {\rm d}x \leq M_0$;

(ⅵ) $\varrho_0 \in \overline{\mathcal H}^{2m} \bigcap L_1^{(1)}$.

Then, there exists a certain $T_*>0$ and a solution $\varrho(\theta,t;x,\omega)$ to (2.3) on $(0,T_*)$ such that

where

Before proceeding to the proof of Theorem 4.1, we prepare the following lemmas.

Lemma 4.2. Let $f=f(\theta,t;x,\omega)$ be a function satisfying $|\!|\!|f(t)|\!|\!| < \infty \; \forall t > 0$ in general. Then, the following estimates hold:

(ⅰ) If $f(\theta,t;x,\omega) \geq 0$ and $\int_\Omega f(\theta,t;x,\omega) \; {\rm d}\theta=1$ for each $(t,x,\omega)$, then

(ⅱ)

(ⅲ)

Proof. Here we only show the proof of (ⅱ). By virtue of the Schwarz inequality, we have

The estimate (ⅲ) is obtained in a similar manner, and statement (ⅰ) is obtained easily.

Hereafter, for arbitrary $T>0$, we use the notation

Lemma 4.3. For an arbitrary $T>0$, if there exists a solution $\varrho$ to (2.3) that belongs to ${\mathcal V}_{(*)}^{2m}(T)$, then the following issues hold:

(ⅰ) $\varrho(\theta,t;x,\omega) \geq 0 \; {\rm a.e.} \; {\rm on} \; \theta\in \Omega, \; t \in (0,T), \; (x,\omega) \in {\bf{R}}^2$

(ⅱ) $\displaystyle \int_\Omega \varrho(\theta,t;x,\omega) \; {\rm d}\theta = 1 \ \ t \in (0,T), \; (x,\omega) \in {\bf{R}}^2.$

Remark 2. If $m$ in Theorem 4.1 is large enough ($m \geq 2$), then the first statement in Theorem 4.1 holds in the pointwise sense due to the Sobolev's embedding theorem and maximum principle (see[11]).

Proof. We verify the statement by using the Stampacchia's truncation method. Let us define

It is obvious that $\varrho_+$ and $\varrho_-$ stand for the positive and negative parts of $\varrho$, respectively, satisfying $\varrho = \varrho_+ - \varrho_-$. Then, by multiplying (2.3)$_1$ by $\varrho_-$ and integrating over $\Omega$, we obtain

Here we used the estimate:

which is derived by integration by parts and Lemma 4.2. Taking into account $\varrho_- \bigr|_{t=0}=0,$ we arrive at $\varrho_- =0$ on $t \in (0,T)$ by virtue of the Gronwall's inequality. This implies the first statement. The second one is proved by direct calculation (see, [24]), and we omit the proof.

4.3.   Proof of Theorem 4.1

We carry out the proof of Theorem 4.1 in three steps below.

(ⅰ) Existence of a solution $\varrho$ that belongs to ${\mathcal V}_{(*)}^{2m}(T_*)$ on a certain time interval $(0,T_*)$;

(ⅱ) proof of $\varrho \in L_1^{(1)}(T_*)$, and consequently $\varrho \in {\mathcal V}^{2m}(T_*)$;

(ⅲ) uniqueness of the solution in ${\mathcal V}^{2m}(T_*)$.

4.3.1.   Existence of solution

We apply the semi-discrete approximation used by Sjöberg[32] and Tsutsumi[36] for the study of the KdV equation. Let us take $N \in {\bf{N}}$, $h=2\pi/N$, $\theta_j = jh \; (j=1,2,\ldots,N)$ and denote as $\partial_+$ and $\partial_-$ the difference operators defined by

Now, instead of problem (2.3), we consider the following differential-difference equation:

where

Note that this function is defined on the continuous interval with respect to $\theta.$

By virtue of the estimate (4.9) we will verify later, it is clear that problem (4.1) above has a unique solution $\varrho_N$ for every $N \in {\bf{N}}$. Then we derive some bounds for $\varrho_N$ and its differences, which are uniform with respect to $N$. To do that, in the space of grid-functions we define the scalar product and the norm by

respectively. As Sjöberg[32] and Tsutsumi[36] did, we assume $N=2n+1$ with $n \in {\bf{N}}$, and then

forms an orthonormal basis with respect to the scalar products $(\cdot,\cdot)_h$ and $(\cdot,\cdot)$.

The following lemmas are due to the work by Sjöberg[32] and Tsutsumi[36]; thus, we omit the proof here.

Lemma 4.4. If $f_1$ is a real $N$-periodic grid-function, i.e., if $f_1(x_j)=f_1(x_{j+N}) \;(j=1,2,\ldots,N)$, and if $f_2(x)$ is another real $N$-periodic grid-function, then the following equalities hold:

Lemma 4.5.  Let $\tau_1$ and $\tau_2$ be non-negative integers with $\tau_1+\tau_2=\tau$, and $\psi$, a function of the form

Then,

holds with some constants $c_{4i} \; (i=2,3).$

For the proof of Lemma 4.5, see Lemma 2.2 in the work by Sjöberg[32]. In addition, the following lemma is useful.

Lemma 4.6. Let $f=f(\theta)$, and $\psi$ be the discrete Fourier series of $f$, that is,

with $b_k = \frac{1}{\sqrt{2\pi}} (e^{i k \theta}, f)_h.$ Then, the discrete version of the Parseval type equality

holds for non-negative integers $\tau_1$ and $\tau_2$.

Proof. We first verify the statement when $\tau_1=\tau_2=0.$ In fact, since

noting that

where $\delta_{kk^\prime}$ is the Kronecker's delta, we have

On the other hand, since

we have

Accordingly, by noting (4.4) again, we have

Finally, the statement holds in case $\tau_1+\tau_2>0$, since $\partial_+^{\tau_1} \partial_-^{\tau_2} \psi$ is the Fourier series of $\partial_+^{\tau_1} \partial_-^{\tau_2} f$. We verify this for $(\tau_1,\tau_2)=(1,0)$ for simplicity. In fact, by definition,

On the other hand,

However, it is easy to see that this equals the rightmost-hand side of (4.5); therefore, $\partial_+ \psi$ is the discrete Fourier series of $\partial_+ f.$ For other pairs of $(\tau_1,\tau_2)$, we are able to show the desired statement in the similar manners.

On the basis of Lemma 4.6, we derive some estimates of $\varrho_N$ and its differences.

Lemma 4.7. The following estimates hold:

where $c_{4(k)} \; (k=1,2)$ and $\tilde{c}_{4(j)}\; (j=1,2,\ldots,m)$ are positive constants independent on $t$, $x$ and $\omega$.

Proof. Let us multiply (4.1)$_1$ by $\varrho_N$. Then, by virtue of Lemma 4.2, we have

Making use of

for two real $N$-periodic grid-functions $f_1$ and $f_2$ in general, where $(f_1)_- (\theta_j) = f_1(\theta_{j-1})$, we have

where $\varepsilon$ is a certain positive constant, and $C_\varepsilon$, a constant dependent on $\varepsilon$ (hereafter we use these notations in the same meaning). Here we have used the estimates

and by means of the mean value theorem and Schwarz's inequality,

where $\theta_0 \in (\theta_{j-i},\theta_j)$. In addition, we have used the Young's inequality

where we take $\varepsilon < D.$ Thus, taking the supremum with respect to $x$ in (4.6), we arrive at

By virtue of the comparison theorem, $\sup_x \|\varrho_N(\cdot,t;x,\omega)\|_h$ is estimated from above by the solution of the ordinary differential equation

which has a solution on $t \in (0,T_*)$ with some $T_* > 0$, and satisfies the estimate of the form

with a constant $c_{44}$ independent on $\omega$. Next, we multiply (4.1)$_1$ by $\partial_+ \partial_- \varrho_N$. With the aid of (4.2)-(4.3), we show some of the elementary calculations below:

Combining these, we have

By taking $h$ sufficiently small, and due to estimate (4.8) we have already obtained, we derive

which yields

for $t \in (0,T_*)$ due to the preceding discussion.

Similarly, multiplying (4.1)$_1$ by $\partial_+^2 \partial_-^2 \varrho_N$ leads to

We expand the last term in the left-hand side as

Each term is estimated as follows.

Here we applied a similar estimate as (4.7). Thus, by using (4.8) we have the estimate of the form

which yields the boundedness of $\|\partial_+ \partial_- \varrho_N (\cdot,t;x,\omega) \|_h^2$ on $(0,T_*)$. In a similar manner, we obtain the estimates of higher difference terms with respect to $\theta$.

When $m \geq 2$, by virtue of Lemma 4.5, we easily obtain

under the assumptions of Theorem 4.1. Then, as Sjöberg[32] and Tsutsumi[36] did, we consider the discrete Fourier series of $\varrho_N$, which is denoted as $\{\phi_N\}$:

Estimate (4.9) and Lemmas 4.5-4.7 yield that the sequence of functions $\displaystyle \{ \phi_N \}$ is uniformly bounded and equicontinuous on $0 \leq \theta \leq 2\pi$, $0 \leq t \leq T$. With the aid of the Arzera-Ascoli theorem, we see that $\displaystyle \{ \phi_N \}$ contains a subsequence which converges to a certain function $\varrho$ as $N \rightarrow +\infty$. In addition, it is clear that

in $L_2(\Omega)$ for each $(t,x,\omega)$. Therefore, this $\varrho$ is the desired solution to (4.1).

Finally, we discuss the regularity of $\varrho$ with respect to $x$. Let us define

which clearly satisfy

Thus, under the assumptions of Theorem 4.1, we can show that

with respect to $x$ for each $(\theta,t)$ and $\omega$ in the same line with the arguments by Lavrentiev[24]. The same argument holds concerning the regularity with respect to $\omega$, and we finally arrive at the desired regularity of $\varrho$.

4.3.2.   Proof of $\varrho \in {\mathcal V}^{2m}(T_*)$

Before proceeding to the uniqueness part, we mention that the solution that was guaranteed to exist in the previous process also belongs to $L_1^{(1)}(T_*)$ thanks to Lemma 4.3, and consequently, to ${\mathcal V}^{2m}(T_*).$ This is directly obtained from Lemma 4.3.

4.3.3.   Uniqueness of solution

Finally, we discuss the uniqueness part of the statement. Assume that there exist two solutions $\varrho_i \; (i=1,2)$ to (2.3) on $(0,T^*)$ with the same initial data, and let us define $\tilde{\tilde{\varrho}} \equiv \varrho_1-\varrho_2.$

Then, it satisfies:

We multiply (4.10)$_1$ by $\tilde{\tilde{\varrho}}$. Then, integration by parts yields

These yield

By virtue of the Gronwall's inequality and the initial condition (4.10)$_3$, we have

which indicates the uniqueness of the solution in the desired function space.

4.4.   Global-in-time solution

Now we discuss the global-in-time solvability of (2.3). Let $\varrho_0$ be provided, which satisfies the assumption (ⅵ) in Theorem 4.1. In accordance with Theorem 4.1, we first construct the local-in-time solution $\varrho$ on $t \in (0,T_*)$, where $T_*$ is the time provided in that theorem. Then, we have the a-priori estimate.

Lemma 4.8. Let $T>0$ be an arbitrary number. If there exists a solution to 2.3 on $(0,T)$, estimates of the form

hold with certain constants $c_{4(k)}^\prime$ independent of $t$.

Proof. For the sake of simplicity, we introduce the notation $\tilde{\varrho} \equiv \varrho-\bar{\varrho}$ and derive the estimate of its norm which leads to the desired estimates. From (2.3), it is obvious that $\tilde{\varrho}$ satisfies

Multiply (4.12)$_1$ by $\tilde{\varrho}$, and making use of Lemma 4.2 and the periodicity of $F[\tilde{\varrho}+\bar{\varrho},\tilde{\varrho}+\bar{\varrho}]$ with respect to $\theta$ yields

On the other hand, in the same line as the arguments by Lavrentiev[24], we have

where we have applied the Young's inequality in the last inequality.

Thus, after taking the supremum with respect to $(x,\omega)$, we have the estimate of the form

Therefore, if we take $\varepsilon^\prime$ so small that $\varepsilon^\prime < D$ holds, then by virtue of the classical Gronwall's inequality, we have the estimate of the form (see, for instance, p.85 of[35])

Next, we show the estimate of $\tilde{\varrho}^{(1,0)}$, which satisfies

Then, due to the estimates

and the Young's inequality, and taking the supremum with respect to $x$ and $\omega$, we have the estimate of the form

with constants $\chi_1^{(i,0)} \; (i=0,1)$. Now we divide the second term in the left-hand side of (4.13) into two terms by using a small constant $\varepsilon>0$, and apply the Poincaré's inequality

to the first term:

Then, we obtain

Summing up this and (4.15) multiplied by a positive constant $m^{(1,0)}$, which will be specified later, we have

Therefore, we take $\varepsilon$, $\varepsilon^\prime$ and $m^{(1,0)}$ in the following manner:

(ⅰ) Take $\varepsilon$ and $\varepsilon^\prime$ so that $\varepsilon^\prime < \varepsilon < D$ holds;

(ⅱ) Then, take $m^{(1,0)}>0$ so small that

hold.

Then, in the same line with the deduction of (4.14), we have

Similarly, for $k=2,3,\ldots,(2m-2),$ we have the estimate of the form:

For estimates of $\tilde{\varrho}^{(l,0)} \; (l=2m-1,2m)$, we introduce the Friedrichs mollifier $\Phi \in C_0^\infty(\Omega)$ with respect to $\theta$[28], and define

with a constant $\delta>0.$ We also define

for functions $f_1, \; f_2 \in L_2({\bf{R}})$ in general. Note that when they are defined on $\Omega$, we extend them onto ${\bf{R}}$ preserving the regularity[23].

Now, by operating $\Phi_\delta *$ to (4.12)$_1$, we have

where

Operating $\displaystyle \Bigl( \frac{\partial}{\partial \theta} \Bigr)^l$ to (4.18) yields

Hereafter we use the notation ${}_l C_i$ to denote the binomial coefficient of l choose i.

Then, we multiply (4.19) by $\tilde{\varrho}_{(\delta)}^{(l,0)}$, and with the similar process as above, obtain the estimate

with constants $\{\chi_l^{(j,0)}\}_{j=0}^l.$ By noting $\displaystyle \lim_{\delta \rightarrow 0} |\!|\!| \widetilde{H}_{(\delta)}^{(l,0)} |\!|\!| =0$, we obtain the following estimate from (4.20) by letting $\delta$ tend to zero:

We now multiply each estimate for $\tilde{\varrho}^{(j,0)} \; (j=1,2,\ldots,2m)$ in (4.15), (4.17) and (4.21) by a positive constant $m^{(j,0)}$, which will be specified later, and sum up them. We also introduce $m^{(0,0)}=1$ for simplicity. These yield

It should be noted that a straightforward estimate, like that by Ha and Xiao[11] will require $D$ to be monotonically and increasingly dependent on $m$. Instead, we partially apply the Poincaré's inequality as before:

We take $\varepsilon$, $\varepsilon^\prime$ and $\{m^{(i,0)}\}_{i=1}^{2m}$ as follows.

(ⅰ) Take $\varepsilon$ and $\varepsilon^\prime$ so that $\varepsilon^\prime < \varepsilon < D$ hold;

(ⅱ) Take $m^{(1,0)}>0$ so small that

hold;

(ⅲ) Take $m^{(2,0)}>0$ so small that

hold;

(ⅳ) As for $\{m^{(i,0)}\}_{i=3}^{2m-1}$, take them so small inductively that:

hold;

(ⅴ) Finally, take $m^{(2m,0)}$ so small that

hold.

Thus, (4.22) becomes

where the coefficients of each term in the left-hand side are all positive. Thus, as we have obtained (4.14) and (4.16), the Gronwall's inequality again yields

Now, by virtue of Lemma 4.8, $\varrho\Bigr|_{t=T_*}$ satisfies the assumption (ⅵ) in Theorem 4.1 imposed on $\varrho_0$. Thereby, we are able to extend $\varrho$ over the time interval $(T_*,2T_*)$ as a solution to (2.3). This solution, which is now defined over $(0,2T_*)$, again satisfies the estimate (4.11), and we then extend it over the region $(2T_*,3T_*)$. Iterating this procedure sufficiently many times, we are able to obtain the solution over the time interval $(0,T)$ for arbitrary $T>0.$ We summarize these arguments as follows.

Theorem 4.9. Let $T$ be an arbitrary positive number. Then, under the same assumptions as in Theorem 4.1, there exists a solution $\varrho(\theta,t;x,\omega) \in {\mathcal V}^{2m}(T)$ to 2.3 on $(0,T)$.

Remark 3. From these considerations, it is obvious that Theorem 4.9 holds when $g$ is compactly supported. we can also extend these arguments when $g$ is the Dirac's delta function, i.e. $g(\omega)=\delta(0)$, since Lemma 4.2 holds in that case.

Corollary 4.10. Under the assumptions in Theorem 4.1 with (ⅳ) replaced by $g(\omega)=\delta(0)$, the same statement as in Theorem 4.9 holds.

5.   Stability of trivial solution

In this section, we discuss the nonlinear stability of the trivial stationary solution $\bar{\varrho}\equiv \frac{1}{2\pi}$. As in the proof of Lemma 4.8, we introduce the notation $\tilde{\varrho} \equiv \varrho-\bar{\varrho},$ and derive the estimate of its norm with respect to time.

The asymptotic stability of $\tilde{\varrho}$ reads

Theorem 5.1. In addition to the assumptions in Theorem 4.1, we assume

Then, $\tilde{\varrho}$ is asymptotically stable in $\overline{\mathcal H}^{2m}$ and satisfies the inequality

with certain positive constants $c_{5i} \; (i=1,2)$.

Proof. The line of the argument is similar to that of Lemma 4.8, but this time we have to confirm the non-positiveness of the left-hand side of the energy type inequalities. First, let us multiply (4.12)$_1$ by $\tilde{\varrho}$. By making use of Lemma 4.2, we have the estimate:

Thus, we have the energy estimate

For the estimate of $|\!|\!|\tilde{\varrho}(t) |\!|\!|$ only, if $D > 2\pi^2 M_0 \Bigl( C_1 + \widetilde{C}_1 \Bigr)$ holds, with the aid of the Poincaré's inequality, we have

which leads to the estimate of the form

by virtue of the Gronwall's inequality. Estimate (5.2) implies the asymptotic stability of $\bar{\varrho}$ in $\overline{\mathcal H}^0.$

Next, we show the estimate up to the first-order spatial derivative. First, by applying the Poincaré's inequality to the second term of the left-hand side in (5.1) partially, as in the previous section, we have

Then, let us sum up (5.3) and (4.15), and we have

As in the previous section, we take $\varepsilon$ and $m^{(1,0)}$ in the following manner:

(ⅰ) Take $\varepsilon$ so small that $\displaystyle \frac{D-\varepsilon}{4\pi^2} -\frac{M_0}{2}\Bigl( C_1 + \widetilde{C}_1 \Bigr) > 0$ holds;

(ⅱ) Then, take $m^{(1,0)}>0$ so small that

hold.

By applying the Gronwall's inequality to (5.4), these lead to the estimate of the form

that is, the asymptotic stability of $\bar{\varrho}$ holds in $\overline{\mathcal H}^1.$

Similarly, we make use of (4.17) and (4.21) to deduce

Summing up (5.3), (5.5) and (5.6) multiplied by constants $\{ m^{(i,0)} \}_{i=0}^{2m}$ with $m^{(0,0)}=1$, we arrive at

Now, we determine $\varepsilon$ and $\{m^{(i,0)}\}_{i=1}^{2m}$ as follows.

(ⅰ) Take $\varepsilon$ so small that $\displaystyle \frac{D-\varepsilon}{4\pi^2} - \frac{M_0}{2}\Bigl( C_1+\widetilde{C}_1\Bigr) >0$;

(ⅱ) Take $m^{(1,0)}>0$ so small that

hold;

(ⅲ) Take $m^{(2,0)}>0$ so small that

hold;

(ⅵ) As for $\{m^{(q,0)}\}_{q=3}^{2m-1}$, take them so small inductively that:

hold;

(ⅴ) Finally, take $m^{(2m,0)}$ so small that

hold.

These yield the estimate of the form

which directly leads to the desired statement.

Remark 4. For the original Kuramoto-Sakaguchi equation (1.6), Ha[11] deduced a similar result concerning the asymptotic stability of $\bar{\varrho}=1/2\pi$. On the basis of some regularity and decay of $\varrho_0$ with respect to $\omega$, they estimated the $H^3$ norm with weight with respect to $\omega$. Their arguments are based on the energy method, which is similar to the one presented here. The advantage of our method is to make the estimate sharper by dividing the terms of a higher derivative into two terms before applying the Poincaré's inequality. Indeed, as we have mentioned before, we need monotonically increasing $D$ with respect to $m$ with the procedure used by Ha and Xiao[11]. As we mentioned in Section 2, the original Kuramoto model (1.3) (or (1.6)) can be regarded as a specific case of the non-local coupling model (2.2) (resp. (2.3))[20]. Therefore, the stability of the incoherent state in (1.6) is verified through similar arguments and assumptions in Theorem 4.1 (see also[13]).

6.   Vanishing diffusion limit

Finally, we discuss the vanishing limit of the diffusion coefficient. In order to show the dependency of the solution on the diffusion coefficient clearly, we denote the solution of (2.3) as $\varrho_{(D)}$, whereas that of (2.5) is denoted as $\varrho_{(0)}.$ As we have stated, Ha and Xiao[11] held a similar discussion for the original Kuramoto-Sakaguchi equation (1.6). However, they estimated the norm of $\varrho_{(D)}$ by using the polynomial of $D$, which resulted in the convergence in $L_\infty(\Omega)$ with respect to $\theta$. In the discussion below, we apply the compactness argument for deriving higher order convergence than their result. This method is also applicable to the original Kuramoto-Sakaguchi equation[13]. We first prepare some lemmas below.

Lemma 6.1. Let $T>0$ be an arbitrary number. Then, the sequence $\{{\varrho}_{(D)}\}_{D>0}$ is uniformly bounded with respect to $D$ in ${\mathcal V}^{2m}(T)$.

Proof. What we have to verify are

with some constants $\{c_{l,k}(T)\}_{2k+l \leq 2m}$ and $\{c_{l,k}^\prime(T)\}_{2k+l \leq 2m}$ dependent on $T$. For $k=0$, we have already verified (6.1) in the previous section.

Now we provide the estimates of the temporal derivative of ${\varrho}$ inductively. Let us assume that the boundedness of $\{{\varrho}_{(D)}^{(0,i)}\}_{D>0}$ has been proven for $i=1,2,\ldots,k-1 \; (k \leq m-1).$ Then, by noting

and Lemma 4.2, we have the estimate of the form

where $\{ c_{k,j}^{\prime\prime} \}$ are the positive constants. By virtue of the assumption of the induction and the Gronwall's inequality, we then have the estimate of the form

As for $k=m$, we use the mollifier again. Recalling the notation defined in Section 4, we consider

where ${\varrho}_{(D)(\delta)} = \Phi_\delta * \varrho_{(D)},$ and

We operate the temporal derivative $\displaystyle \Bigl( \frac{\partial}{\partial t} \Bigr)^m$ to (6.5) and obtain

Then, after the energy type estimate, we make $\delta$ tend to zero. By making use of the fact

we obtain (6.4) with $k=m$ with the aid of the assumptions of the induction.

Next, we estimate the term including both temporal and spatial derivatives. We only show the case $l+2k+2 \geq 2m$, when the mollifier is necessary again. By applying the $l$-th order spatial derivative to (6.6) with $m$ replaced with $k$, we have

Now we show some examples of the energy type estimates. In case $(i,j)=(l,k)$, we have

For $(i,j)=(l-1,k)$, we have

Otherwise, we have

By combining these and (6.7), and applying the Young's and Schwarz's inequalities, we derive the energy estimate of the form

After making $\delta \rightarrow 0$ and applying the Gronwall's inequality, we have the estimate of the form

Thus, we have shown the first part of the statement. Estimate (6.2) is derived from (6.3), (6.8), and the estimates we have already obtained. These complete the proof.

By virtue of Lemma 6.1, we see that the sequence $\{\varrho_{(D)}\}_{D>0}$ includes a sub-sequence, denoted as $\{\varrho_{(D)}\}$ again, which is convergent in the weak-star sense as $D$ tends to zero:

Then, in the relationship

if we make $D$ tend to zero, we have

which means $\hat{\varrho}^\prime = \frac{\partial \hat{\varrho}}{\partial t}$.

The next lemma clarifies the space to which this sequence converges.

Lemma 6.2. Let $T>0$ be an arbitrary number. Then, the sequence $\{{\varrho}_{(D)}\}_{D>0}$ forms the Cauchy sequence in ${\mathcal V}^{2m}(T)$.

Proof. Let us define $\breve{\varrho} \equiv \varrho_{(D)} - \varrho_{(D^\prime)}$ for arbitrary $D$, $D^\prime>0$, which satisfies

Then, with the aid of the estimates

we have

Thus, by virtue of the Gronwall's inequality and the fact that $\breve{\varrho}|_{t=0}=0$, we have the estimate of the form

This implies that the sequence $\{\varrho_{(D)}\}_{D>0}$ is the Cauchy sequence in $\overline{\mathcal H}^0.$ The estimates of higher derivative terms are obtained in a similar manner. To show this, we subtract (6.7) with $D$ replaced with a certain $D^\prime > 0$ from itself and obtain

where $\breve{\varrho}_{(\delta)} \equiv \Phi_\delta * \breve{\varrho}$, and $\breve{H}_{(\delta)} = H_{(D)(\delta)}-H_{(D^\prime)(\delta)}.$

We show inductively that $\breve{\varrho}_{(\delta)}^{(l,k)} \rightarrow 0$ hold as $D$ and $D^\prime$ tend to zero based on the assumption that $\{\varrho_{(D)}^{(i,j)}\}_{D>0}$ form the Cauchy sequences with $0\leq i \leq l-1$ and $0 \leq j \leq k$. Let us show some examples of the estimates. In case $(i,j)=(l,k)$, we have

In case $(i,j)=(l-1,k)$, we have

Otherwise, we have

After applying the Schwarz's inequality, we make $\delta$ tend to zero. Then, on the basis of the assumption of the induction, we have the estimate of the form:

where $C(D,D^\prime) \rightarrow 0$ as $D$ and $D^\prime$ tend to zero.

From these considerations, the sequence $\{\varrho_{(D)}\}_{D>0}$ forms the Cauchy sequence in ${\mathcal V}^{2m}(T)$. This completes the proof of Lemma 6.2.

By Lemma 6.2, we see that $\hat{\varrho}$ belongs to ${\mathcal V}^{2m}(T)$. Now, we show that $\hat{\varrho}$ certainly satisfies (2.5). To do this, we take an arbitrary function $h(\theta,t) \in C^1(0,T;C_0^\infty(\Omega))$ satisfying $h(\theta,t)|_{t=T}=0$, $h(\theta,t)|_{t=0}\ne 0$, and consider

In virtue of (6.9)-(6.10), if we make $D$ tend to zero,

On the other hand, thanks to the Rellich's theorem[27], we have

strongly as $D \rightarrow 0$; therefore,

holds. Thus, we arrive at

which means that $\hat{\varrho}$ certainly satisfies (2.5)$_1$. Next, integrate (6.11) and (6.12) by part with respect to $t$, and the assumptions on $h$ yield

respectively. Comparing (6.13) and (6.14) with the aid of (6.9)-(6.10) implies $\hat{\varrho}|_{t=0} = \varrho_{(0)}$, and so the initial condition (2.5)$_3$ is satisfied. The periodicity (2.5)$_2$ obviously holds due to the function space to which $\hat{\varrho}$ belongs. Thus, $\hat{\varrho}=\varrho_{(0)}$. We summarize these arguments as follows.

Theorem 6.3. Let $T>0$ be an arbitrary number, and impose the same assumptions as in Theorem 4.1. Then, the solution $\varrho_{(D)}$ of (2.3) converges to that of (2.5) in ${\mathcal V}^{2m}(T)$, which is denoted as $\varrho_{(0)}$, as $D$ tends to zero.

7.   Conclusion

In this paper, we discussed the mathematical analysis of the nonlinear Fokker-Planck equation of Kuramoto's non-local coupling model of oscillators. We first showed the local and global-in-time solvability, and then the nonlinear asymptotic stability of the incoherent state. Finally, we verified the existence of the vanishing diffusion limit solution as the diffusion coefficient tends to zero.

Our future work will be concerned with the mathematical stability analysis of the chimera state of this model and the coupled oscillator model on the complex graph. We will also tackle the bifurcation problem.