Destabilization of synchronous periodic solutions for patch models: a criterion by period functions

Abstract.

In this paper, we study the destabilization of synchronous periodic solutions for patch models. By applying perturbation theory for matrices, we derive asymptotic expressions of the Floquet spectra and provide a destabilization criterion for synchronous periodic solutions arising from closed orbits or degenerate Hopf bifurcations in terms of period functions. Finally, we apply the main results to the well-known two-patch Holling-Tanner model.

Key words and phrases:

Patch model, destabilization, periodic solution, bifurcation.

2020 Mathematics Subject Classification:

Primary 34D20, 35B10; Secondary 34C23, 34C25.

†Email: schen@ccnu.edu.cn.

‡Corresponding author: hjc@mail.ccnu.edu.cn.

This work was partly supported by the National Natural Science Foundation of China (Grant No. 12101253, 12231008) and the Scientific Research Foundation of CCNU (Grant No. 31101222044).

Shuang Chen ${}^{\,{\dagger}}$ , Jicai Huang ${}^{\,{\ddagger}}$

School of Mathematics and Statistics, and Hubei Key Laboratory of Mathematical Sciences,

Central China Normal University, Wuhan, Hubei 430079, China

1. Introduction

Patch models have been extensively used to understand the spatial spread of infectious diseases and the effect of population dispersal on the total abundance and the total populations distribution (see, for instance, [1, 3, 12, 13, 14, 15, 25, 32, 34] and the references therein).

In this paper, we investigate a general $n$ -patch model with cross-diffusion-like couplings:

\displaystyle\ \ \begin{aligned} \frac{du_{1}^{(j)}}{dt}&=\delta\sum_{l=1}^{m}\sum_{i=1}^{n}\left(d_{1l}\left(u_{l}^{(i)}-u_{l}^{(j)}\right)\right)+f_{1}(u_{1}^{(j)},u_{2}^{(j)},\cdot\cdot\cdot,u_{m}^{(j)}),\ \ j=1,2,...,n,\\ \frac{du_{2}^{(j)}}{dt}&=\delta\sum_{l=1}^{m}\sum_{i=1}^{n}\left(d_{2l}\left(u_{l}^{(i)}-u_{l}^{(j)}\right)\right)+f_{2}(u_{1}^{(j)},u_{2}^{(j)},\cdot\cdot\cdot,u_{m}^{(j)}),\ \ j=1,2,...,n,\\ &\ \ \ \cdots\cdots\\ \frac{du_{m}^{(j)}}{dt}&=\delta\sum_{l=1}^{m}\sum_{i=1}^{n}\left(d_{ml}\left(u_{l}^{(i)}-u_{l}^{(j)}\right)\right)+f_{m}(u_{1}^{(j)},u_{2}^{(j)},\cdot\cdot\cdot,u_{m}^{(j)}),\ \ j=1,2,...,n.\end{aligned}

(1.1)

In the setting of population dynamics, the state variables $u_{i}^{(j)}$ are the population densities of $i$ -th species in the $j$ -th patch, $n\geq 2$ denotes the number of the patches, $d_{ij}$ are the diffusion coefficients, and $\delta>0$ indicates the coupling strength. Let

\mathcal{U}(t)=\left(\begin{array}[]{c}U^{1}\\ U^{2}\\ \vdots\\ U^{n}\end{array}\right),\ \mathcal{F}(\mathcal{U}(t))=\left(\begin{array}[]{c}F(U^{1})\\ F(U^{2})\\ \vdots\\ F(U^{n})\end{array}\right),\ \mathcal{E}=\left(\begin{array}[]{cccc}(n-1)E&-E&\cdots&-E\\ -E&(n-1)E&\cdots&-E\\ \vdots&\vdots&\ddots&\vdots\\ -E&-E&\cdots&(n-1)E\end{array}\right),

where $U^{j}=(u_{1}^{(j)},u_{2}^{(j)},\cdot\cdot\cdot,u_{m}^{(j)})^{T}$ , $F=(f_{1},f_{2},\cdot\cdot\cdot,f_{m})^{T}$ and $E=(d_{kl})_{m\times m}$ is a $m\times m$ matrix. Then we can rewrite patch model (1.1) in the compact form

\displaystyle\frac{d}{dt}\mathcal{U}(t)=-\delta\mathcal{E}\mathcal{U}(t)+\mathcal{F}(\mathcal{U}(t)).

(1.2)

Throughout this paper, we assume that $f_{i}$ ( $i=1,...,m$ ) are sufficiently smooth, and use ^T to denote the transpose of a matrix or a vector.

A solution $\mathcal{U}_{0}$ of patch model (1.2) is called a synchronous periodic solution if $\mathcal{U}_{0}=(\phi^{T},...,\phi^{T})^{T}$ and $\phi:\mathbb{R}\to\mathbb{R}^{m}$ is a periodic function with the minimum period $P>0$ . In this case, this periodic function $\phi$ is also a periodic solution of the underlying kinetic system

\displaystyle\frac{dU}{dt}=:\dot{U}=F(U),

(1.3)

where $U=(u_{1},u_{2},...,u_{m})^{T}\in\mathbb{R}^{m}$ and $F(U)=(f_{1}(U),f_{2}(U),\cdot\cdot\cdot,f_{m}(U))^{T}$ . Furthermore, if $\phi$ is a (Lyapunov) stable periodic solution of the kinetic system (1.3), then the corresponding synchronous periodic solution $\mathcal{U}_{0}$ is also stable in system (1.2) without the cross-diffusion-like couplings. A natural question arises:

$\bullet$

Can the synchronous periodic solution $\mathcal{U}_{0}(t)$ become unstable in system (1.2) with the cross-diffusion-like couplings?

The instability driven by the cross-diffusion-like couplings is called the destabilization of synchronous periodic solutions for patch model (1.1). It is also called the Turing instability of periodic solutions [35], in order to celebrate Turing’s discovery in [30], i.e., diffusion could destabilize stable equilibrium solutions of reaction-diffusion systems.

We are interested in the destabilization of synchronous periodic solutions for patch model (1.1). This is directly motivated by various phenomena and problems arising from real-world applications. For example, [4] recently shown that unstable states play a vital role in transient dynamics and the resilience of ecological systems to environmental change. [25] once found that the destabilization of periodic solutions in chemically reacting systems can lead to complicated oscillations and chaos. It is significant to understand the effect of the connectivity of subregions on infectious disease transmission [12, 15]. Along this direction, we also need to further investigate the impact of the cross-diffusion-like couplings on periodic oscillations.

The main obstacle to investigate the destabilization of synchronous periodic solutions is that it is difficult to analyze the Floquet spectra of the related linearizations about synchronous periodic solutions. The obstacle becomes evident after we present the linearization of patch model (1.1) about a synchronous periodic solution $\mathcal{U}_{0}$ , i.e., the following periodic system

\displaystyle\frac{d}{dt}Y(t)=(-\delta\mathcal{E}+\widehat{J}(t))Y(t),

(1.4)

where $Y(t)\in\mathbb{R}^{mn}$ , $\widehat{J}(t)={\rm diag}(J(t),J(t),\cdot\cdot\cdot,J(t))$ and $J(t)=F_{U}(\phi(t))$ . Note that general patch models always involve multiple patches. Then the related periodic system (1.4) is high-dimensional, and it is challenging to give the explicit expressions of the Floquet spectra for high-dimensional periodic systems, even for three-dimensional systems. See some reviews in such as [6, 7, 19, 35].

Bifurcations of invariant sets (e.g. equilibria, periodic solutions, homoclinic loops, heteroclinic loops, etc) in a parametrically perturbed system give rise to periodic solutions. The period of bifurcating periodic solutions arising from an invariant set can be well defined in terms of system parameters when parameters are near a bifurcation point. This gives the period function of bifurcating periodic solutions [5, 17]. In our recent work [7], we provided a criterion for the destabilization of synchronous periodic solutions bifurcating from double homoclinic loops. Based on the Lyapunov-Schmidt reduction, we obtained the characteristic function to determine the Floquet spectra associated with synchronous periodic solutions, while the period functions of bifurcating periodic solutions are not well-defined at bifurcation points [7, 17]. It is also interesting and challenging to deal with bifurcating periodic solutions whose period function is at least continuously differentiable. Some typical examples include periodic solutions arising from the parametric perturbations of equilibria and limit cycles. See [17] for instance.

Our goal is to further give criteria for the destabilization of synchronous periodic solutions in the term of period functions if the related period functions are continuously differentiable. Similar criteria were previously proposed for the diffusion-derived instability of spatially homogeneous periodic solutions in reaction-diffusion systems. For example, Maginu [27] in 1979 and Ruan [29] in 1998 once considered the diffusion-derived instability of spatially homogeneous periodic solutions for reaction-diffusion systems in the entire space. They established the relation between the Floquet spectra and the period functions of bifurcating periodic solutions for some certain perturbations of the kinetic systems. After that, they gave criteria for the instability of spatially homogeneous periodic solutions in terms of the dominant term of the related period functions (see Lemma 2 in [27] and Formula (4.16) in [29]).

In order to give a criterion for the destabilization, we consider a single parametric perturbation of the kinetic system (1.3) as follows:

\displaystyle(I_{m}+\epsilon E)\frac{dU}{dt}=F(U),

(1.5)

where $\epsilon\in\mathbb{R}$ is a small parameter, and $E$ and $F$ are defined as in (1.2). If the kinetic system (1.3) has a hyperbolic periodic solution $\phi$ , then by the classical bifurcation theory [2, 5, 36], there exists a sufficiently small $\epsilon_{0}>0$ such that the perturbed system (1.5) admits a family of bifurcating periodic solutions $\phi_{\epsilon}$ for $\epsilon\in(-\epsilon_{0},\epsilon_{0})$ that bifurcate from the hyperbolic periodic solution $\phi$ . Let $P(\epsilon)$ denote the minimum period of $\phi_{\epsilon}$ , and call $P(\cdot)$ the period function of the family of bifurcating periodic solutions $\phi_{\epsilon}$ . Furthermore, the period function $P(\epsilon)$ is continuously differentiable in the open set $(-\epsilon_{0},\epsilon_{0})$ . By applying the dominant term of the period function $P(\epsilon)$ , we give a criterion for the destabilization of synchronous periodic solutions (see Theorem 2.3). The argument is mainly based on the perturbation theory for matrices that was developed in our recent work [7]. Actually, we present the asymptotic expression of the related Floquet spectra for synchronous periodic solutions.

It is worth mentioning that our result actually improves Yi’s work [35]. Yi recently introduced a single parametric perturbation of patch model (1.1) as follows:

\displaystyle(I_{mn}+\epsilon\mathcal{E})\frac{d}{dt}\mathcal{U}(t)=\mathcal{F}(\mathcal{U}(t)),

(1.6)

where $I_{mn}$ is the $mn\times mn$ identity matrix, and $\mathcal{E}$ and $\mathcal{F}$ are defined as in (1.2). Under the assumption that the perturbed patch model (1.6) possesses a family of bifurcating periodic solutions

((U_{p}^{1}(\cdot,\epsilon))^{T},(U_{p}^{2}(\cdot,\epsilon))^{T},...,(U_{p}^{n}(\cdot,\epsilon))^{T})^{T}

with period functions $\widehat{P}(\epsilon)$ that vanish asymptotically in one patch and persist in the other $(n-1)$ patches, i.e.,

\displaystyle((U_{p}^{1}(t,\epsilon))^{T},(U_{p}^{2}(t,\epsilon))^{T},...,(U_{p}^{n}(t,\epsilon))^{T})^{T}\to((\phi(t))^{T},(\phi(t))^{T},\cdot\cdot\cdot,(\phi(t))^{T},0)^{T}\ \ \ \mbox{ as }\ \epsilon\to 0,

uniformly in $t\in\mathbb{R}$ . To determine the destabilization of synchronous periodic solutions for patch model (1.1), following the idea of Maginu [27] and Ruan [29], Yi [35] presented a criterion for the destabilization of the synchronous periodic solution $\mathcal{U}_{0}=(\phi^{T},...,\phi^{T})^{T}$ in terms of the first-order derivative of $\widehat{P}(\epsilon)$ at $\epsilon=0$ . Here, we give a simpler criterion that only requires $P^{\prime}(0)<0$ , where $P^{\prime}(0)$ denotes the first-order derivative of $P(\epsilon)$ at $\epsilon=0$ , without any additional conditions. See Theorem 2.5.

The present paper is devoted to the destabilization of synchronous periodic solutions for patch models and appears to be our second paper on this topic. We refer to our first paper [7], where the related period functions are not continuously differentiable at bifurcations points. Instead, here we explore the case of continuously differentiable period functions. The theory we developed is applied to give criteria for the destabilization of synchronous periodic solutions, bifurcating from closed orbits [31, 36] and degenerate Hopf bifurcation [9, 11], in $n$ -patch models with two-dimensional kinetic systems. It was once found that diffusion-driven instability of periodic solutions for reaction-diffusion systems can not be induced by the identical diffusion rates [21]. As an easy by-product, we prove that patch model with the identical diffusion rates never undergoes the destabilization of synchronous periodic solutions.

This paper is organized as follows. The criterion for the destabilization of synchronous periodic solutions is given in section 2, and then we apply our results to general $n$ -patch models with two-dimensional kinetic systems. Sections 3.1 and 3.2 are devoted to synchronous periodic solutions arising from closed orbits and degenerate Hopf bifurcation, respectively. Finally, the well-known two-patch Holling-Tanner model is provided in Section 3.3 to illustrate the main results.

2. Destabilization of synchronous periodic solutions

In this section, we give the main result on the destabilization of synchronous periodic solutions. The proof is based on two fundamental lemmas on the perturbation theory for matrices that were developed in our recent work [7]. We present them in Appendix A for convenience.

Let $\phi$ be a periodic solution with the minimum period $P>0$ for the kinetic system (1.3). Then the linearization of the kinetic system (1.3) about this periodic solution $\phi$ is governed by

\displaystyle\frac{dU}{dt}=F_{U}(\phi(t))U=J(t)U,\ \ \ \ \ \ U\in\mathbb{R}^{m}.

(2.1)

Throughout this section, we make the following assumption:

(H)

All Floquet multipliers $\gamma_{1},...,\gamma_{m}$ of the linearized system (2.1) satisfy

$\displaystyle\gamma_{1}=1,\ \ |\gamma_{2}|<1,...,|\gamma_{m}|<1.$ (2.2)

Under this assumption, the periodic solution $\phi$ is stable with respect to the kinetic system (1.3) and the synchronous periodic solution $\mathcal{U}_{0}(t)=(\phi(t)^{T},...,\phi(t)^{T})^{T}$ is also stable with respect to patch model (1.6) if $\delta=0$ . By the discussion in the Introduction, there exists a sufficiently small $\epsilon_{0}>0$ such that the perturbed system (1.5) with $\epsilon\in(-\epsilon_{0},\epsilon_{0})$ admits a family of bifurcating periodic solutions $\phi_{\epsilon}$ for $\epsilon\in(-\epsilon_{0},\epsilon_{0})$ that bifurcate from this stable solution $\phi$ . We use $P(\epsilon)$ to denote the related period function.

We are interested in the effect of the cross-diffusion-like couplings on the stability of this synchronous periodic solution. Consider the linearized system (1.4) of patch model (1.2) about the periodic solution $\mathcal{U}_{0}$ . Note that all Floquet multipliers of the linearized system (2.1) satisfy (2.2). Then $n$ Floquet multipliers of system (1.4) with $\delta=0$ are one, and all other Floquet multipliers have modulii less than one. Let $\Psi(t,\delta)$ denote the principal fundamental matrix solution of system (1.4). Then we can compute that the kernel $\ker(\Psi(P,0)-I_{mn})$ is spanned by the following vectors in $\mathbb{R}^{mn}$ :

\xi_{1}=\left(\begin{array}[]{c}\phi_{1}(0)\\ 0\\ \vdots\\ 0\end{array}\right),\ \ \xi_{2}=\left(\begin{array}[]{c}0\\ \phi_{1}(0)\\ \vdots\\ 0\end{array}\right),\cdots,\ \xi_{n}=\left(\begin{array}[]{c}0\\ 0\\ \vdots\\ \phi_{1}(0)\end{array}\right),

where $0$ is the zero vector in $\mathbb{R}^{m}$ , $\xi_{j}\in\mathbb{R}^{mn}$ and $\phi_{1}$ is given by

\displaystyle\phi_{1}(t):=\frac{d}{dt}\phi(t),\ \ \ \ \ t\in\mathbb{R}.

Here we use the fact that $\phi_{1}$ is a solution of the linearized system (2.1). Next we give an important property on the monodromy operator $\Psi(P,\delta)$ .

Lemma 2.1.

Suppose that $\phi(\cdot,\epsilon)$ $(\epsilon\in(-\epsilon_{0},\epsilon_{0}))$ are periodic solutions bifurcating from $\phi(\cdot)$ in the perturbed system (1.5). Define

\displaystyle\tilde{\xi}_{1}:=(\xi_{1}+\xi_{2}+\cdot\cdot\cdot+\xi_{n}),\ \ \ \ \ \tilde{\xi}_{j}(\delta):=\frac{1}{n}(\xi_{j-1}-\xi_{n})+\delta(\eta_{j-1}-\eta_{n}),\ \ \ \ \ j=2,...,n,

where $\delta\geq 0$ and $\eta_{j}$ for $j=1,...,n$ are given by

\eta_{1}=\left(\begin{array}[]{c}\frac{\partial}{\partial\epsilon}\phi(0,0)\\ 0\\ \vdots\\ 0\end{array}\right),\ \ \eta_{2}=\left(\begin{array}[]{c}0\\ \frac{\partial}{\partial\epsilon}\phi(0,0)\\ \vdots\\ 0\end{array}\right),\cdots,\ \eta_{n}=\left(\begin{array}[]{c}0\\ 0\\ \vdots\\ \frac{\partial}{\partial\epsilon}\phi(0,0)\end{array}\right).

Then

\displaystyle\Psi(P,\delta)\tilde{\xi}_{1}=\tilde{\xi}_{1}

(2.3)

and for sufficiently small $\delta\geq 0$ ,

\displaystyle\Psi(P,\delta)\tilde{\xi}_{j}(\delta)=(1-\delta nP^{\prime}(0))\tilde{\xi}_{j}(\delta)+O(\delta^{2}),\ \ \ \ \ j=2,...,n.

(2.4)

Proof.

Note that $\widetilde{\mathcal{U}}(t):=\left(\left(\phi_{1}(t)\right)^{T},\left(\phi_{1}(t)\right)^{T},...,\left(\phi_{1}(t)\right)^{T}\right)^{T}$ satisfies (1.4). Then the monodromy operator $\Psi(P,\delta)$ has the Floquet multiplier $\gamma_{1}(\delta)\equiv 1$ for each $\delta\geq 0$ and satisfies (2.3).

In order to prove (2.4), we first consider the case $j=2$ . Let $\Psi_{1}(t,\delta)$ be the unique solution of the initial value problem:

\displaystyle\begin{aligned} \frac{d}{dt}\Psi_{1}(t,\delta)=&\ (-\delta\mathcal{E}+\widehat{J}(t))\Psi_{1}(t,\delta),\\ \Psi_{1}(0,\delta)=&\ \tilde{\xi}_{2}(\delta).\end{aligned}

(2.5)

Differentiating with respect to $\delta$ and then setting $\delta=0$ , we have

\displaystyle\begin{aligned} \frac{\partial}{\partial t}\frac{\partial}{\partial\delta}\Psi_{1}(t,0)=&\ -\mathcal{E}\Psi_{1}(t,0)+\widehat{J}(t)\frac{\partial}{\partial\delta}\Psi_{1}(t,0),\\ \frac{\partial}{\partial\delta}\Psi_{1}(0,0)=&\ \eta_{1}-\eta_{n}.\end{aligned}

(2.6)

Note that $\Psi_{1}(t,0)$ satisfies

\displaystyle\frac{d}{dt}\Psi_{1}(t,0)=\widehat{J}(t)\Psi_{1}(t,0),\ \ \ \ \ \Psi_{1}(0,0)=\frac{1}{n}(\xi_{1}-\xi_{n}).

Then we can compute

\Psi_{1}(t,0)=\left(\begin{array}[]{c}\frac{1}{n}\phi_{1}(t)\\ 0\\ \vdots\\ 0\\ -\frac{1}{n}\phi_{1}(t)\end{array}\right),\ \ \ \ \ \mathcal{E}\Psi_{1}(t,0)=\left(\begin{array}[]{c}E\phi_{1}(t)\\ 0\\ \vdots\\ 0\\ -E\phi_{1}(t)\end{array}\right),

where $0$ is the zero vector in $\mathbb{R}^{m}$ .

Recall that $\phi(t,\epsilon)$ is a periodic solution bifurcating from $\phi(t)$ in the perturbed system (1.5). Then substituting $U(t)=\phi(t,\epsilon)$ into (1.5) and differentiating system (1.5) with respect to $\epsilon$ , we have

\displaystyle\frac{\partial}{\partial t}\left(\frac{\partial}{\partial\epsilon}\phi(t,0)\right)=-E\phi_{1}(t)+J(t)\left(\frac{\partial}{\partial\epsilon}\phi(t,0)\right).

(2.7)

By the above equation, we can check that the initial value problem (2.6) has the unique solution

\displaystyle\frac{\partial}{\partial\delta}\Psi_{1}(t,0)=\left(\begin{array}[]{c}\frac{\partial}{\partial\epsilon}\phi(t,0)\\ 0\\ \vdots\\ 0\\ -\frac{\partial}{\partial\epsilon}\phi(t,0)\end{array}\right).

Note that $\Psi_{1}(P,\delta)$ is analytic in the parameter $\delta$ . We have the following expansion:

\displaystyle\Psi_{1}(P,\delta)=\Psi_{1}(P,0)+\delta\frac{\partial}{\partial\delta}\Psi_{1}(P,0)+O(\delta^{2}).

This together with (2.5) and the fact that

\Psi(P,0)\frac{\xi_{1}-\xi_{n}}{n}=\frac{\xi_{1}-\xi_{n}}{n}

yields

\displaystyle\Psi(P,\delta)\tilde{\xi}_{2}(\delta)=\frac{1}{n}(\xi_{1}-\xi_{n})+\delta\left(\begin{array}[]{c}\frac{\partial}{\partial\epsilon}\phi(P,0)\\ 0\\ \vdots\\ 0\\ -\frac{\partial}{\partial\epsilon}\phi(P,0)\end{array}\right)+O(\delta^{2}).

(2.14)

Since $\phi(t+P(\epsilon),\epsilon)=\phi(t,\epsilon)$ , differentiating with respect to $\epsilon$ , we have

\displaystyle\frac{\partial}{\partial\epsilon}\phi(P,0)=-P^{\prime}(0)\phi_{1}(0)+\frac{\partial}{\partial\epsilon}\phi(0,0).

Substituting this into (2.14) yields

\displaystyle\Psi(P,\delta)\tilde{\xi}_{2}(\delta)=(1-\delta nP^{\prime}(0))\tilde{\xi}_{2}(\delta)+O(\delta^{2}).

Similarly, we can prove that (2.4) holds for $j=3,...,n$ . This finishes the proof. ∎

In order to study the stability of $\mathcal{U}_{0}(t)$ with respect to system (1.2), we give the following lemma on the Floquet multipliers of (1.5).

Lemma 2.2.

For sufficiently small $|\epsilon|\geq 0$ , let $P(\epsilon)$ be the period function of periodic solutions $\phi(t,\epsilon)$ bifurcating from $\phi(t)$ in the perturbed system (1.5). Then for sufficiently small $\delta>0$ , system (1.5) has $n$ Floquet multipliers $\gamma_{j}(\delta)$ perturbed from one, and $\gamma_{j}(\delta)$ are in the form

\displaystyle\gamma_{1}(\delta)\equiv 1,\ \ \ \ \gamma_{j}(\delta)=1-nP^{\prime}(0)\delta+O(\delta^{2}),\ \ \ j=2,...,n.

Proof.

Note that the monodromy operator $\Psi(p,\delta)$ is analytic in $\delta$ . Then for sufficiently small $\delta>0$ , the monodromy operator $\Psi(p,\delta)$ has $n$ eigenvalues $\gamma_{j}(\delta)$ ( $j=1,2,...,n$ ) arising from one. This finishes the proof for the first statement.

By the proof of Lemma 2.1, we have $\gamma_{1}(\delta)\equiv 1$ for $\delta\geq 0$ . To prove the expressions for $\gamma_{j}(\delta)=1-nP^{\prime}(0)\delta+O(\delta^{2})$ , $j=2,...,n$ , we define $W_{+}(\delta)$ and $\Lambda_{+}(\delta)$ by $W_{+}(\delta)=\left(\tilde{\xi}_{1},\tilde{\xi}_{2}(\delta),...,\tilde{\xi}_{n}(\delta)\right)$ and

\displaystyle\Lambda_{+}(\delta)=\left(\begin{array}[]{cccc}1&0&\cdot\cdot\cdot&0\\ 0&1-\delta nP^{\prime}(0)&\cdot\cdot\cdot&0\\ \cdot\cdot\cdot&\cdot\cdot\cdot&\cdot\cdot\cdot&\cdot\cdot\cdot\\ 0&0&\cdot\cdot\cdot&1-\delta nP^{\prime}(0)\end{array}\right),

where $\tilde{\xi}_{1}$ , $\tilde{\xi}_{2}(\delta)$ ,…, $\tilde{\xi}_{n}(\delta)$ are defined in Lemma 2.1. Then by (2.3) and (2.4),

\displaystyle\Psi(P,\delta)W_{+}(\delta)=W_{+}(\delta)\Lambda_{+}(\delta)+O(\delta^{2})

(2.16)

for sufficiently small $\delta\geq 0$ .

By (2.2), Lemma A.1 and the definition of $W_{+}(\delta)$ , there exists a small $\delta_{0}>$ and continuous functions $\tilde{\xi}_{j}(\delta)\in\mathbb{R}^{mn}$ for $j=n+1,...,mn$ and $|\delta|<\delta_{0}$ such that

W_{-}(\delta):=\left(\tilde{\xi}_{n+1}(\delta),...,\tilde{\xi}_{mn}(\delta)\right)

satisfies

\displaystyle\Psi(P,\delta)W_{-}(\delta)=W_{-}(\delta)\Lambda_{-}(\delta),\ \ \ |{\rm det}(W_{+}(\delta),W_{-}(\delta))|>0,\ \ \ |\delta|<\delta_{0},

(2.17)

where $\Lambda_{-}(\delta)\in\mathbb{R}^{(m-1)n\times(m-1)n}$ is non-singular and continuous in $\delta$ , and all eigenvalue of $\Lambda_{-}(\delta)$ is bounded away from $\gamma_{j}(\delta)$ in the complex plane for $|\delta|<\delta_{0}$ . By (2.16) and (2.17), we get

\displaystyle\Psi(P,\delta)(W_{+}(\delta),W_{-}(\delta))=(W_{+}(\delta),W_{-}(\delta))\left\{\left(\begin{array}[]{cc}\Lambda_{+}(\delta)&0\\ 0&\Lambda_{-}(\delta)\end{array}\right)+W_{0}(\delta)\right\},

where $W_{0}(\delta)$ is continuous in $\delta$ and $W_{0}(\delta)=O(|\delta|^{2})$ for sufficiently small $|\delta|$ . Therefore, the proof is finished by Lemma A.2. ∎

Now we state the main result in the following.

Theorem 2.3.

Suppose that $\phi$ is a periodic solution with the minimum period $P>0$ for the kinetic system (1.3) that satisfies assumption (H). Let $P(\epsilon)$ denote the period function of bifurcating periodic solutions arising from $\phi$ in the perturbed system (1.5). Then for sufficiently small $\delta>0$ , the synchronous periodic solution $\mathcal{U}_{0}(t)=(\phi(t)^{T},...,\phi(t)^{T})^{T}$ is unstable with respect to patch model (1.1) if $P^{\prime}(0)<0$ .

Proof.

If $P^{\prime}(0)<0$ , then by Lemma 2.2, the Floquet multipliers $\gamma_{j}(\delta)$ $(j=2,...,n)$ satisfy

\displaystyle\gamma^{\prime}_{j}(0)=-nP^{\prime}(0)>0.

Recall that $\gamma_{j}(0)=1$ for $j=2,...,n$ . Then for sufficiently small $\delta>0$ , the linearized system (1.4) has at least $(n-1)$ Floquet multipliers that have modulii greater than one. This proves that $\mathcal{U}_{0}(t)$ is unstable. Thus, the proof is now complete. ∎

Remark 2.4.

More recently, Yi [35] considered patch model (1.1) and studied the destabilization of the synchronous periodic solution $\mathcal{U}_{0}$ using the perturbed patch model (1.6). In order to prove the instability of $\mathcal{U}_{0}$ , [35] required that the perturbed patch model (1.6) has a periodic solution $((U_{p}^{1}(\cdot,\epsilon))^{T},(U_{p}^{2}(\cdot,\epsilon))^{T},...,(U_{p}^{n}(\cdot,\epsilon))^{T})^{T}$ that asymptotically vanishes in one patch and persists in the other $(n-1)$ patches, i.e.,

\displaystyle((U_{p}^{1}(t,\epsilon))^{T},(U_{p}^{2}(t,\epsilon))^{T},...,(U_{p}^{n}(t,\epsilon))^{T})^{T}\to((\phi(t))^{T},(\phi(t))^{T},\cdot\cdot\cdot,(\phi(t))^{T},0)^{T},\ \ \ \mbox{ as }\ \epsilon\to 0.

Here without this condition, we rigorously prove the instability of $\mathcal{U}_{0}$ under the condition that $P^{\prime}(0)<0$ . This improves the result in [35].

3. Application to $n$ -patch models with two-dimensional kinetic system

As an application of Theorem 2.3, we consider an $n$ -patch model with two-dimensional kinetic system

		$\displaystyle\frac{du_{j}}{dt}=\delta\sum_{i\in\Omega}\left(d_{11}(u_{i}-u_{j})+d_{12}(v_{i}-v_{j})\right)+f(u_{j},v_{j},\alpha),\ \ \ \ \$		$\displaystyle\ \ j\in\Omega:=\{1,2,...,n\},$		(3.1)
		$\displaystyle\frac{dv_{j}}{dt}=\delta\sum_{i\in\Omega}(d_{21}(u_{i}-u_{j})+d_{22}(v_{i}-v_{j}))+g(u_{j},v_{j},\alpha),\ \ \ \ \$		$\displaystyle\ \ j\in\Omega:=\{1,2,...,n\},$		(3.1)

where $(u_{j}(t),v_{j}(t))^{T}\in\mathbb{R}^{2}$ represent the population densities of two species in the $j$ -th patch, $n$ is an integer greater or equal to $2$ , and $\delta>0$ indicates the coupling strength. Here the parameters $d_{ij}$ in (3.1) indicate the diffusion coefficients. In particular, when $i\neq j$ , the parameters $d_{ij}$ are the cross-diffusion rates.

The underlying kinetic system of the above patch model reads as the following system of ordinary differential equations

		$\displaystyle\dot{u}=f(u,v,\alpha),$		(3.2)
		$\displaystyle\dot{v}=g(u,v,\alpha),$		(3.2)

where $\alpha$ in $\mathbb{R}$ is a system parameter, and the functions $f$ and $g$ are sufficiently smooth. If the kinetic system (3.2) has a stable periodic solution $\phi(t)\in\mathbb{R}^{2}$ , then patch model (3.1) has a synchronous periodic solution $((\phi(t))^{T},...,(\phi(t))^{T})^{T}\in\mathbb{R}^{2n}$ that are also stable in the absence of diffusion. Our aim is to discuss whether this stable synchronous periodic solution could become unstable in the presence of diffusion. Specially, we focus on periodic solutions arising from closed orbits and degenerate Hopf bifurcation, which are called large- and small-amplitude bifurcating periodic solutions, respectively.

3.1. Application to large-amplitude bifurcating periodic solutions

Consider the kinetic system with $\alpha=\alpha_{0}$ which is in the form

		$\displaystyle\dot{u}=f(u,v,\alpha_{0})=:f(u,v),$		(3.3)
		$\displaystyle\dot{v}=g(u,v,\alpha_{0})=:g(v,v).$		(3.3)

Let $\psi(t)=(u_{0}(t),v_{0}(t))^{T}\in\mathbb{R}^{2}$ be a periodic solution of (3.3) that satisfies the following hypothesis:

•

The periodic solution $\psi(t)$ is stable and has minimum period $p>0$ . The related Floquet multipliers $\gamma$ and $\tilde{\gamma}$ satisfy

$\displaystyle\gamma=1,\ \ \ 0<\tilde{\gamma}<1.$ (3.4)

In the view of bifurcation theory, this periodic solution $\psi(t)$ can bifurcate from a perturbation of a closed orbit in the kinetic system (3.2) with $\alpha$ near $\alpha_{0}$ .

It is clear that the linearization of system (3.3) about $\psi(t)$ is in the form

\displaystyle\frac{d}{dt}\left(\begin{array}[]{c}\tilde{u}(t)\\ \tilde{v}(t)\end{array}\right)=\left(\begin{array}[]{cc}f_{u}(\psi(t))&f_{v}(\psi(t))\\ g_{u}(\psi(t))&g_{v}(\psi(t))\end{array}\right)\left(\begin{array}[]{c}\tilde{u}(t)\\ \tilde{v}(t)\end{array}\right).

Then by Lemma 7.3 in [19, p.120], we have

\displaystyle\tilde{\gamma}=\gamma\tilde{\gamma}=\exp\left(\int_{0}^{p}\left(f_{u}(\psi(t))+g_{v}(\psi(t))\right)dt\right)<1.

Following the discussion in Section 2, we consider an auxiliary planar system of the form

\displaystyle(I_{2}+\epsilon D)\left(\begin{array}[]{c}\dot{u}\\ \dot{v}\end{array}\right)=\left(\begin{array}[]{c}f(u,v,\alpha)\\ g(u,v,\alpha)\end{array}\right),

(3.10)

where $I_{2}$ is the $2\times 2$ identity matrix, and the matrix $D$ is in the form

\displaystyle D=\left(\begin{array}[]{cc}d_{11}&d_{12}\\ d_{21}&d_{22}\end{array}\right).

Consider the perturbation system (3.10) with $\alpha=\alpha_{0}$ . It has the expansion with respect to $\epsilon$ as follows:

\displaystyle\ \ \ \begin{aligned} \dot{u}&=f(u,v)-\epsilon(d_{11}f(u,v)+d_{12}g(u,v))+O(\epsilon^{2})=:f(u,v)+\epsilon f_{1}(u,v)+O(\epsilon^{2}),\\ \dot{v}&=g(u,v)-\epsilon(d_{21}f(u,v)+d_{22}g(u,v))+O(\epsilon^{2})=:g(u,v)+\epsilon g_{1}(u,v)+O(\epsilon^{2}).\end{aligned}

(3.12)

We summarize some results on periodic solutions bifurcating from $\psi(t)$ in system (3.12).

Lemma 3.1.

Suppose that the kinetic system (3.3) has a stable periodic solution satisfying the conditions in (3.4). Then there exists a small constant $\epsilon_{0}>0$ such that for each $\epsilon$ with $0\leq|\epsilon|<\epsilon_{0}$ , system (3.12) has exactly one limit cycle $\psi(t,\epsilon):=(u_{0}(t,\epsilon),v_{0}(t,\epsilon))\in\mathbb{R}^{2}$ with period $P(\epsilon)$ bifurcating from $\psi(t)$ . Moreover, the period function $P(\epsilon)$ has the expansion of the form

\displaystyle\ \ \ \ \ \ \ \begin{aligned} P(\epsilon)=&\ p+\epsilon\int_{0}^{p}\left(\left(\frac{2fg(f_{u}-g_{v})+(g^{2}-f^{2})(f_{u}+g_{v})}{(f^{2}+g^{2})^{3/2}}|_{(u,v)^{T}=\psi(t)}\right)\right.\\ &\ \times\frac{1}{\sqrt{f^{2}(\psi(t))+g^{2}(\psi(t))}}\left(I(t)+\frac{1}{1-\tilde{\gamma}}I(p)e^{h(t)}\right)\\ &\left.-\left(\frac{ff_{1}+gg_{1}}{f^{2}+g^{2}}|_{(u,v)^{T}=\psi(t)}\right)\right)dt+O(\epsilon^{2})\\ =:&\ p+P_{1}\epsilon+O(\epsilon^{2}),\ \ \ \ \ 0\leq|\epsilon|<\epsilon_{0},\end{aligned}

(3.13)

where $f_{1}$ and $g_{1}$ are defined as in (3.12), and $I(t)$ and $h(t)$ are in the form

	$\displaystyle I(t)\!\!\!$	$\displaystyle=$	$\displaystyle\!\!\!e^{h(t)}\int_{0}^{t}\left(e^{-h(s)}(fg_{1}-gf_{1})\|_{(u,v)^{T}=\psi(s)}\right)ds,$		(3.14)
	$\displaystyle h(t)\!\!\!$	$\displaystyle=$	$\displaystyle\!\!\!\int_{0}^{t}\left(f_{u}(\psi(s))+g_{v}(\psi(s))\right)ds.$

Proof.

We can prove the existence of perturbed periodic solutions using [5, Theorem 2.1, p.352]. Now we give the expansion of the period function $P(\epsilon)$ . For sufficiently small $|\epsilon|$ , by the formulas (2.16), (2.18) and (4.5) in [31], we have

	$\displaystyle P(\epsilon)\!\!\!$	$\displaystyle=$	$\displaystyle\!\!\!p+\int_{0}^{p}\left(\left(\frac{2fg(f_{u}-g_{v})+(g^{2}-f^{2})(f_{u}+g_{v})}{(f^{2}+g^{2})^{3/2}}\|_{(u,v)^{T}=\psi(t)}\right)\rho_{1}(t,\epsilon)\right.$
			$\displaystyle\!\!\!\left.-\epsilon\left(\frac{ff_{1}+gg_{1}}{f^{2}+g^{2}}\|_{(u,v)^{T}=\psi(t)}\right)\right)dt+O(\epsilon^{2}),$

where $\rho_{1}(t,\epsilon)$ is in the form

	$\displaystyle\rho_{1}(t,\epsilon)\!\!\!$	$\displaystyle=$	$\displaystyle\!\!\!\frac{\epsilon}{\sqrt{f^{2}(\psi(t))+g^{2}(\psi(t))}}I(t)+\frac{\sqrt{f^{2}(\psi(0))+g^{2}(\psi(0))}}{\sqrt{f^{2}(\psi(t))+g^{2}(\psi(t))}}e^{h(t)}c(\epsilon),$
	$\displaystyle c(\epsilon)\!\!\!$	$\displaystyle=$	$\displaystyle\!\!\!\epsilon\frac{1}{(1-\tilde{\gamma})\sqrt{f^{2}(\psi(0))+g^{2}(\psi(0))}}I(p)+O(\epsilon^{2}).$

Thus, we can compute (3.13). This finishes the proof. ∎

Now we state the results on destabilization of large-amplitude periodic solutions.

Proposition 3.2.

Suppose that the kinetic system (3.2) with $\alpha=\alpha_{0}$ has a stable periodic solution $\psi(t)=(u_{0}(t),v_{0}(t))\in\mathbb{R}^{2}$ that satisfies (3.4). If the constant $P_{1}$ defined in (3.13) satisfies $P_{1}<0$ , then for sufficiently small $\delta>0$ , the corresponding synchronous periodic solution $(\psi(t),...,\psi(t))^{T}$ is unstable with respect to patch model (3.1).

Proof.

If $P_{1}<0$ , then we have

\displaystyle P^{\prime}(0)<0.

Therefore, the proof is finished by Theorem 2.3. ∎

We remark that the formula in Lemma 3.1 gives an analytic formula to determine the sign of $P_{1}$ , although the expression of $P_{1}$ is complicated. This formula also provides a possibility to numerically give criteria for the destabilization of synchronous periodic solutions. As an easy by-product of Proposition 3.2, we have that the destabilization of synchronous periodic solutions can not be induced by the identical diffusion rates. More precisely, we have the following result.

Proposition 3.3.

Suppose that the kinetic system (3.2) with $\alpha=\alpha_{0}$ has a stable periodic solution $\psi(t)=(u_{0}(t),v_{0}(t))\in\mathbb{R}^{2}$ that satisfies (3.4). Then for sufficiently small $\delta>0$ , the corresponding synchronous periodic solution $(\psi(t),...,\psi(t))^{T}$ is stable with respect to patch model (3.1) with the identical diffusion rates.

Proof.

We first prove that the period function $P(\epsilon)$ satisfies $P^{\prime}(0)>0$ . If $d_{11}=d_{22}=d_{0}>0$ and $d_{12}=d_{21}=0$ , then the functions $f_{1}$ and $g_{1}$ in (3.12) are in the form

\displaystyle f_{1}(u,v)=-d_{0}f(u,v),\ \ \ g_{1}(u,v)=-d_{0}g(u,v),\ \ \ (u,v)\in\mathbb{R}^{2}.

This yields

\displaystyle f(u,v)g_{1}(u,v)-g(u,v)f_{1}(u,v)=0.

Then $I(t)$ in (3.14) satisfies $I(t)\equiv 0$ for $t\in\mathbb{R}$ . So we can compute

\displaystyle\begin{aligned} P(\epsilon)=p+\epsilon d_{0}p+O(\epsilon^{2}).\end{aligned}

(3.15)

This implies $P^{\prime}(0)>0$ .

Now we prove the stability of $(\psi(t),...,\psi(t))^{T}$ with respect to patch model (3.1). Let the notations be defined as in Section 2. Since $P^{\prime}(0)>0$ , by Lemma 2.2 we have that the Floquet multipliers $\gamma_{j}(\delta)$ satisfy $\gamma_{1}(\delta)\equiv 1$ for $\delta\geq 0$ and $\gamma^{\prime}_{j}(0)<0$ for $j=2,3,...,n$ . Consequently, $\gamma_{j}(\delta)<1$ for sufficiently small $\delta>0$ and $j=2,...,n$ . This together with the fact that $\Psi(p,\delta)$ is continuous with respect to $\delta$ and the condition (3.4) yields that all of the Floquet multipliers of (1.4) have modulii less than one except $\gamma_{1}(\delta)$ . Then $(\psi(t),...,\psi(t))^{T}$ is stable with respect to patch model (3.1) with the identical diffusion rates. This finishes the proof. ∎

3.2. Application to small-amplitude bifurcating periodic solutions

In this section, we study the destabilization of periodic solutions arising from Hopf bifurcation. Based on the normal form theory and the formal series method, we give the conditions under which the destabilization of Hopf bifurcating periodic solutions appears.

Without loss of generality, throughout this section we make the following hypotheses:

•

The functions $f$ and $g$ in the kinetic system (3.2) are $C^{\infty}$ in $(u,v,\alpha)$ , and $f(0,0,\alpha)=g(0,0,\alpha)=0$ for all $\alpha\in\mathbb{R}$ .
•

The kinetic system (3.2) has a center-type equilibrium at the origin $O:=(0,0)^{T}$ for $\alpha=0$ , that is, the Jacobian matrix $J(O)$ of system (3.2) with $\alpha=0$ at the origin has a pair of purely imaginary eigenvalues $\lambda_{1,2}=\pm{\bf i}\mu_{0}$ for $\mu_{0}>0$ .

Let $J(\alpha)$ denote the Jacobian matrix of the kinetic system (3.2) at the origin, that is,

J(\alpha)=\left(\begin{array}[]{cc}f_{u}(0,0,\alpha)&f_{v}(0,0,\alpha)\\ g_{u}(0,0,\alpha)&g_{v}(0,0,\alpha)\end{array}\right)=:\left(\begin{array}[]{cc}J_{11}(\alpha)&J_{12}(\alpha)\\ J_{21}(\alpha)&J_{22}(\alpha)\end{array}\right).

Then we can compute

\displaystyle J_{11}(0)+J_{22}(0)=0,\ \ \ J_{11}(0)J_{22}(0)-J_{12}(0)J_{21}(0)=\mu_{0}^{2}>0.

(3.16)

By [9, Lemma 1.1, p.384], the kinetic system (3.2) with $\alpha=0$ in complex coordinates has the following Poincaré-Birkhoff normal form

\displaystyle\dot{z}={\bf i}\mu_{0}z+C_{1}z^{2}\bar{z}+C_{2}z^{3}\bar{z}^{2}+\cdot\cdot\cdot+C_{k}z^{k+1}\bar{z}^{k}+O(|z|^{2k+3}).

The constants $C_{j}$ are called the $j$ th Lyapunov coefficients of system (3.2) with $\alpha=0$ at the center-type equilibrium $O$ . If the Lyapunov coefficients $C_{j}$ satisfy

\displaystyle{\rm Re}(C_{1})=\cdot\cdot\cdot={\rm Re}(C_{k-1})=0,\ \ \ {\rm Re}(C_{k})\neq 0,

(3.17)

then we say that the kinetic system (3.2) could undergo a Hopf bifurcation of order $k$ for $k\geq 1$ at the origin, and the origin is a weak focus of order $k$ .

Recall that the auxiliary system is defined by (3.10). When $\epsilon=0$ , system (3.10) is reduced to the kinetic system (3.2). It is clear that for sufficiently small $|\epsilon|$ , one can transform system (3.10) into the following system

\displaystyle\left(\begin{array}[]{c}\dot{u}\\ \dot{v}\end{array}\right)=(\mathcal{M}(\epsilon))^{-1}\left(\begin{array}[]{cc}1+\epsilon d_{22}&-\epsilon d_{12}\\ -\epsilon d_{21}&1+\epsilon d_{11}\end{array}\right)\left(\begin{array}[]{c}f(u,v,\alpha)\\ g(u,v,\alpha)\end{array}\right),

(3.24)

where

\mathcal{M}(\epsilon)=\det(I_{2}+\epsilon D)=(d_{11}d_{22}-d_{12}d_{21})\epsilon^{2}+(d_{11}+d_{22})\epsilon+1.

Note that system (3.24) always has an equilibrium at the origin for all $\alpha\in\mathbb{R}$ . A direct computation yields that the Jacobian matrix $J(\alpha,\epsilon)$ of system (3.24) at the origin is in the form

J(\alpha,\epsilon)=(\mathcal{M}(\epsilon))^{-1}\left(\begin{array}[]{cc}1+\epsilon d_{22}&-\epsilon d_{12}\\ -\epsilon d_{21}&1+\epsilon d_{11}\end{array}\right)\left(\begin{array}[]{cc}J_{11}(\alpha)&J_{12}(\alpha)\\ J_{21}(\alpha)&J_{22}(\alpha)\end{array}\right),

and the trace $\mathcal{T}(\alpha,\epsilon)$ of $J(\alpha,\epsilon)$ and the determinant $\mathcal{D}(\alpha,\epsilon)$ of $J(\alpha,\epsilon)$ are given by

\displaystyle\begin{aligned} \mathcal{T}(\alpha,\epsilon)=&\ (\mathcal{M}(\epsilon))^{-1}\left\{J_{11}(\alpha)+J_{22}(\alpha)\right.\\ &\left.+\epsilon\left(d_{22}J_{11}(\alpha)+d_{11}J_{22}(\alpha)-d_{12}J_{21}(\alpha)-d_{21}J_{12}(\alpha)\right)\right\},\\ \mathcal{D}(\alpha,\epsilon)=&\ (\mathcal{M}(\epsilon))^{-1}\det(J(\alpha))=(\mathcal{M}(\epsilon))^{-1}(J_{11}(\alpha)J_{22}(\alpha)-J_{12}(\alpha)J_{21}(\alpha)).\end{aligned}

(3.25)

Next we state the results on the periodic solutions bifurcating from the origin.

Lemma 3.4.

Let $\lambda_{1,2}(\alpha)=A(\alpha)\pm{\bf i}B(\alpha)$ denote the eigenvalues of the matrix $J(\alpha)$ . Suppose that the kinetic system (3.2) with $\alpha=0$ has a weak focus of order $k$ at the origin, and $A^{\prime}(0)\neq 0$ . Then the following statements hold:

(i)

There exist two small constants $\epsilon_{0}>0$ and $\tilde{r}_{0}>0$ , and a smooth function $\alpha(r_{0},\epsilon)$ for $0<r_{0}<\tilde{r}_{0}$ and $|\epsilon|\leq\epsilon_{0}$ such that system (3.24) with $\alpha=\alpha(r_{0},\epsilon)$ has exactly one limit cycle $(u(t,r_{0},\epsilon),v(t,r_{0},\epsilon))$ with period $P(r_{0},\alpha(r_{0},\epsilon))$ near the origin.
(ii)

If the $k$ th Lyapunov coefficient $C_{k}$ satisfies ${\rm Re}(C_{k})<0$ (resp. ${\rm Re}(C_{k})>0$ ), then the perturbed limit cycle is stable (resp. unstable), and each perturbed limit cycle $(u(t,r_{0},\epsilon),v(t,r_{0},\epsilon))$ passes through $(r_{0},0)$ .

Proof.

Let $\lambda_{1,2}(\alpha,\epsilon)=A(\alpha,\epsilon)\pm{\bf i}B(\alpha,\epsilon)$ denote the eigenvalues of the Jacobian matrix $J(\alpha,\epsilon)$ . System (3.24) can be transformed into

\displaystyle\begin{aligned} \dot{u}&=A(\alpha,\epsilon)u-B(\alpha,\epsilon)v+f_{1}(u,v,\alpha,\epsilon),\\ \dot{v}&=B(\alpha,\epsilon)u+A(\alpha,\epsilon)v+g_{1}(u,v,\alpha,\epsilon).\end{aligned}

(3.26)

See the detailed proof in [35, Appendix A]. Set $(u,v)=(r\cos\theta,r\sin\theta)$ . Then

\displaystyle\begin{aligned} \dot{r}&=A(\alpha,\epsilon)r+\cos\theta\cdot f_{1}(r\cos\theta,r\sin\theta,\epsilon)+\sin\theta\cdot g_{1}(r\cos\theta,r\sin\theta,\epsilon),\\ r\dot{\theta}&=B(\alpha,\epsilon)r+\cos\theta\cdot g_{1}(r\cos\theta,r\sin\theta,\epsilon)-\sin\theta\cdot f_{1}(r\cos\theta,r\sin\theta,\epsilon).\end{aligned}

(3.27)

Let $(r(t,r_{0},\alpha,\epsilon),\theta(t,r_{0},\alpha,\epsilon))$ denote the solution of (3.27) with $(r(0),\theta(0))=(r_{0},0)$ . Then for sufficiently small $|\alpha|+|\epsilon|$ , we can define the displacement map by

\displaystyle H(r_{0},\alpha,\epsilon)=r(2\pi,r_{0},\alpha,\epsilon)-r_{0}.

Since $H(0,\alpha,\epsilon)=0$ for all $\alpha$ and $\epsilon$ , we write $H(r_{0},\alpha,\epsilon)$ as

\displaystyle H(r_{0},\alpha,\epsilon)=r_{0}\tilde{H}(r_{0},\alpha,\epsilon),

where $\tilde{H}(r_{0},\alpha,\epsilon)$ has the expansion as the form

\displaystyle\tilde{H}(r_{0},\alpha,\epsilon)=\left(\exp\left(2\pi\frac{A(\alpha,\epsilon)}{B(\alpha,\epsilon)}\right)-1\right)+H_{1}(\alpha,\epsilon)r_{0}+O(r_{0}^{2}).

By (3.17) and (3.27), we further have

\displaystyle\tilde{H}(r_{0},0,0)=\frac{2\pi{\rm Re}(C_{k})}{\mu_{0}}r_{0}^{2k}+O(r_{0}^{2k+1}).

(3.28)

A direct computation yields

\displaystyle\begin{aligned} \frac{\partial\tilde{H}}{\partial\alpha}(0,0,0)=&\ \frac{\partial\tilde{H}}{\partial\alpha}(0,\alpha,0)|_{\alpha=0}\\ =&\ \frac{\partial}{\partial\alpha}\left(\exp\left(2\pi\frac{A(\alpha,0)}{B(\alpha,0)}\right)-1\right)|_{\alpha=0}\\ =&\ \frac{2\pi A^{\prime}(0)}{\mu_{0}}\neq 0.\end{aligned}

(3.29)

This together with the implicit function theorem yields that there exist two small constants $\epsilon_{0}>0$ and $\tilde{r}_{0}>0$ , and a smooth function $\alpha(r_{0},\epsilon)$ for $0<r_{0}<\tilde{r}_{0}$ and $|\epsilon|\leq\epsilon_{0}$ such that $\tilde{H}(r_{0},\alpha(r_{0},\epsilon),\epsilon)=0$ for sufficiently small $|r_{0}|+|\epsilon|$ . Thus, we obtain (i). The statements in (ii) can be proved by (3.27). Therefore, the proof is now complete. ∎

By the above lemma, we can further verify that there exists a constant $\alpha_{0}>0$ such that for each $(\alpha,\epsilon)$ with $|\alpha|<\alpha_{0}$ and $|\epsilon|<\epsilon_{0}$ , system (3.24) has a small-amplitude periodic orbit bifurcating from the origin. Let $T(\alpha,\epsilon)$ denote the corresponding period. Then we call $T(\alpha,\epsilon)$ the period function for this family of periodic solutions arising from Hopf bifurcation.

Note that the period function $T(\alpha,\epsilon)$ depends not only on $\alpha$ but also on $\epsilon$ . Then the results in [8, 17, 20], where the period function depends on a single parameter, is not applicable to $T(\alpha,\epsilon)$ . Now we give the formula of $T(\alpha,\epsilon)$ in the next lemma.

Lemma 3.5.

Let $T(\alpha,\epsilon)$ denote the period function of the Hopf bifurcating periodic solutions. Then $T(\alpha,\epsilon)$ satisfies the formula

\displaystyle\frac{2\pi}{T(\alpha,\epsilon)}=B(\alpha,\epsilon)+\sum_{j=1}^{k}{\rm Im}(C_{j}(\alpha,\epsilon))r_{0}^{2j}+O(r_{0}^{2k+1}),

(3.30)

where $r_{0}$ satisfies

\displaystyle\begin{aligned} r_{0}^{2k}=R_{0}(\alpha,\epsilon)=&-\frac{\epsilon}{2{\rm Re}(C_{k})}\left(d_{22}J_{11}(0)+d_{11}J_{22}(0)-d_{12}J_{21}(0)-d_{21}J_{12}(0)\right)\\ &-\frac{A^{\prime}(0)}{{\rm Re}(C_{k})}\alpha+(|(\alpha,\epsilon)|^{2}).\end{aligned}

(3.31)

for sufficiently small $|\alpha|+|\epsilon|$ .

Proof.

By [9, Lemma 1.1, p.384], system (3.26) can be normalized into

\displaystyle\dot{z}=(A(\alpha,\epsilon)+{\bf i}B(\alpha,\epsilon))z+\sum_{j=1}^{k}C_{j}(\alpha,\epsilon)z^{j+1}\bar{z}^{j}+O(|z|^{2k+3}).

Let $z(t,r_{0},\epsilon)$ denote the bifurcating periodic solution with $z(0,r_{0},\epsilon)=r_{0}$ . Define

\tau=t/T(\alpha,\epsilon),\ \ \ \ \ z(t,r_{0},\epsilon)=r_{0}e^{2\pi{\bf i}\tau}\eta(\tau,r_{0},\epsilon),

where $\eta(\tau+1,r_{0},\epsilon)=\eta(\tau,r_{0},\epsilon)$ for each $\tau\in\mathbb{R}$ . Then

\displaystyle(2\pi{\bf i})\eta+\frac{d\eta}{d\tau}=T(\alpha,\epsilon)\eta\left\{\lambda(\alpha,\epsilon)+\sum_{j=1}^{k}C_{j}(\alpha,\epsilon)r_{0}^{2j}(\eta\bar{\eta})^{j}\right\}+O(r_{0}^{2k+1}).

(3.32)

We expand $\eta(\tau,r_{0},\epsilon)$ as the form

\displaystyle\eta(\tau,r_{0},\epsilon)=\eta_{0}(\tau)+\eta_{1}(\tau,r_{0},\epsilon)+\eta_{2}(\tau,r_{0},\epsilon)+\cdot\cdot\cdot+\eta_{k}(\tau,r_{0},\epsilon)+\cdot\cdot\cdot,

where $\eta_{j}(\tau,\alpha,\epsilon)$ are periodic functions with period one, the homogeneous polynomials of $j$ -th degree with respect to $r_{0}$ and $\epsilon$ , and satisfy

\eta_{0}(0)=1,\ \ \ \ \eta_{j}(0,r_{0},\epsilon)=0,\ \ \ j\geq 1.

Substituting the expansion of $\eta(\tau,r_{0},\epsilon)$ into (3.32) and then comparing the term of the zeroth degree with respect to $r_{0}$ and $\epsilon$ , we have

(2\pi{\bf i})\eta_{0}+\frac{d\eta_{0}}{d\tau}=(2\pi{\bf i})\eta_{0},\ \ \ \ \eta_{0}(0)=1.

This yields $\eta_{0}(\tau)\equiv 1$ for $\tau\in\mathbb{R}$ . Comparing the term of the first degree yields

(2\pi{\bf i})\eta_{1}+\frac{d\eta_{1}}{d\tau}=(2\pi{\bf i})\eta_{1}+d_{1}(\alpha,\epsilon),\ \ \ \ \eta_{1}(0)=0,

where $d_{1}(\alpha,\epsilon)$ is a constant term. Recall that $\eta_{1}(\tau,r_{0},\epsilon)$ is periodic with respect to $\tau$ . Then $\eta_{1}(\tau,r_{0},\epsilon)=0$ for each $\tau\in\mathbb{R}$ . Similarly, we can prove that $\eta_{j}(\tau,r_{0},\epsilon)=0$ for $2\leq j\leq k$ . This together with (3.32) yields (3.30). By applying the implicit function theorem, (3.28), (3.29) and

$\displaystyle\frac{\partial\tilde{H}}{\partial\epsilon}(0,0,0)\!\!\!$	$\displaystyle=$	$\displaystyle\!\!\!\frac{\partial\tilde{H}}{\partial\epsilon}(0,0,\epsilon)\|_{\epsilon=0}$
	$\displaystyle=$	$\displaystyle\!\!\!\frac{\partial}{\partial\epsilon}\left(\exp\left(2\pi\frac{A(0,\epsilon)}{B(0,\epsilon)}\right)-1\right)\|_{\epsilon=0}$
	$\displaystyle=$	$\displaystyle\!\!\!\frac{\pi}{\mu_{0}}\left(d_{22}J_{11}(0)+d_{11}J_{22}(0)-d_{12}J_{21}(0)-d_{21}J_{12}(0)\right),$

we obtain (3.31). This finishes the proof. ∎

Finally, we establish the existence and destabilization of Hopf bifurcating periodic solutions for patch model (3.1).

Proposition 3.6.

Suppose that the kinetic system (3.2) with $\alpha=0$ has a weak focus of order $k$ at the origin, $A^{\prime}(0)\neq 0$ , and the related Lyapunov coefficients $C_{j}$ satisfy

\displaystyle{\rm Re}(C_{1})=\cdot\cdot\cdot={\rm Re}(C_{k-1})=0,\ \ {\rm Re}(C_{k})<0.

Then there exists a constant $\alpha_{0}>0$ such that

(i)

if $|\alpha|<\alpha_{0}$ and $\alpha A^{\prime}(0)>0$ , then the kinetic system (3.2) has a stable periodic solution $\psi(t,\alpha):=(u(t,\alpha),v(t,\alpha))$ bifurcating from the origin, and $\psi(t,\alpha)$ tends to the origin as $\alpha\to 0$ .
(ii)

if $|\alpha|<\alpha_{0}$ and $\alpha A^{\prime}(0)>0$ , then patch model (3.1) has a synchronous periodic solution $\mathcal{U}(t,\alpha):=(\psi(t,\alpha),...,\psi(t,\alpha))^{T}\in\mathbb{R}^{2n}$ .

Proof.

Let $\epsilon=0$ in system (3.24). Since $A^{\prime}(0)\neq 0$ and the $k$ th Lyapunov coefficient $C_{k}$ satisfies ${\rm Re}(C_{k})<0$ , by the formulas (3.28) and (3.29) we have

\displaystyle\alpha(r_{0},0)=-\frac{{\rm Re}(C_{k})}{A^{\prime}(0)}r_{0}^{2k}+O(r_{0}^{2k+1}).

(3.33)

By the proof for Lemma 3.4, there exists a sufficiently small $\alpha_{0}>0$ such that the kinetic system (3.2) with $|\alpha|<\alpha_{0}$ and $\alpha A^{\prime}(0)>0$ has a stable periodic solution $(u(t,\alpha),v(t,\alpha))$ bifurcating from the origin. This implies the existence of perturbed periodic solutions for (3.1). Therefore, the proof is now complete. ∎

Proposition 3.7.

Suppose that the conditions in Theorem 3.6 hold. If $|\alpha|<\alpha_{0}$ , $\alpha A^{\prime}(0)>0$ ,

\displaystyle{\rm Im}(C_{1})=\cdot\cdot\cdot={\rm Im}(C_{k_{1}-1})=0,

\displaystyle\ \ \ {\rm Im}(C_{k_{1}})\neq 0,

(3.34)

where $k_{1}$ satisfies $1\leq k_{1}\leq k$ , and one of the following two conditions holds:

(C1)

$k_{1}<k$ and

${\rm Im}(C_{k_{1}})\left(d_{22}J_{11}(0)+d_{11}J_{22}(0)-d_{12}J_{21}(0)-d_{21}J_{12}(0)\right)>0.$

(C2)

$k_{1}=k$ and

\mu_{0}(d_{11}+d_{22})+\frac{{\rm Im}(C_{k}(0,0))}{{\rm Re}(C_{k}(0,0))}\left(d_{22}J_{11}(0)+d_{11}J_{22}(0)-d_{12}J_{21}(0)-d_{21}J_{12}(0)\right)<0.

then there exists a small constant $\hat{\alpha}_{0}$ with $0<\hat{\alpha}_{0}<\alpha_{0}$ such that for each $\alpha$ with $0<|\alpha|<\hat{\alpha}_{0}$ and sufficiently small $\delta>0$ , the synchronous periodic solution $\mathcal{U}(t,\alpha)$ is unstable with respect to patch model (3.1).

Proof.

We first compute $\frac{\partial B}{\partial\epsilon}(0,0)$ . It is clear that

B(\alpha,\epsilon)=\frac{\sqrt{4\mathcal{D}(\alpha,\epsilon)-(\mathcal{T}(\alpha,\epsilon))^{2}}}{2},

where $\mathcal{D}(\alpha,\epsilon)$ and $\mathcal{T}(\alpha,\epsilon)$ are defined by (3.25). Then

	$\displaystyle\frac{\partial B}{\partial\epsilon}(0,0)\!\!\!$	$\displaystyle=$	$\displaystyle\!\!\!\frac{4\frac{\partial\mathcal{D}}{\partial\epsilon}(\alpha,\epsilon)-2\mathcal{T}(\alpha,\epsilon)\frac{\partial\mathcal{T}}{\partial\epsilon}(\alpha,\epsilon)}{4\sqrt{4\mathcal{D}(\alpha,\epsilon)-(\mathcal{T}(\alpha,\epsilon))^{2}}}\|_{(\alpha,\epsilon)=(0,0)}$
		$\displaystyle=$	$\displaystyle\!\!\!-\frac{\mu_{0}}{2}(d_{11}+d_{22}).$

Now we compute $\frac{\partial T}{\partial\epsilon}(0,0)$ . When $1\leq k_{1}<k$ , by Lemma 3.5 and (3.34) we have

	$\displaystyle-\frac{2\pi}{T^{2}(\alpha,0)}\frac{\partial T}{\partial\epsilon}(\alpha,0)$
	$\displaystyle=\frac{\partial B}{\partial\epsilon}(\alpha,0)+\frac{k_{1}}{k}{\rm Im}(C_{k_{1}}(0,0))(R_{0}(\alpha,0))^{\frac{k_{1}}{k}-1}\frac{\partial R_{0}}{\partial\epsilon}(\alpha,0)+O(\|\alpha\|)$
	$\displaystyle=\left\{-\frac{k_{1}{\rm Im}(C_{k_{1}}(0,0))}{2k{\rm Re}(C_{k}(0,0))}\left(d_{22}J_{11}(0)+d_{11}J_{22}(0)-d_{12}J_{21}(0)-d_{21}J_{12}(0)\right)+O(\alpha)\right\}$
	$\displaystyle\ \ \ \times(R_{0}(\alpha,0))^{\frac{k_{1}}{k}-1}-\frac{\mu_{0}}{2}(d_{11}+d_{22})+O(\alpha),$

where $R_{0}(\alpha,0)$ is defined as in (3.31). Note that $R_{0}(\alpha,0)\to 0^{+}$ as $\alpha\to 0$ and $k_{1}<k$ . Then under the conditions that ${\rm Re}(C_{k})<0$ and ${\bf(C1)}$ , we have that $\frac{\partial T}{\partial\epsilon}(\alpha,0)<0$ for sufficiently small $|\alpha|$ .

When $k_{1}=k$ , by Lemma 3.5 and (3.34) we have

	$\displaystyle-\frac{2\pi}{T^{2}(0,0)}\frac{\partial T}{\partial\epsilon}(0,0)$
	$\displaystyle=\frac{\partial B}{\partial\epsilon}(0,0)+{\rm Im}(C_{k}(0,0))\frac{\partial R_{0}}{\partial\epsilon}(0,0)$
	$\displaystyle=-\frac{\mu_{0}}{2}(d_{11}+d_{22})-\frac{{\rm Im}(C_{k}(0,0))}{2{\rm Re}(C_{k}(0,0))}\left(d_{22}J_{11}(0)+d_{11}J_{22}(0)-d_{12}J_{21}(0)-d_{21}J_{12}(0)\right).$

Then under the condition (C2), we have that $\frac{\partial T}{\partial\epsilon}(\alpha,0)<0$ for sufficiently small $|\alpha|$ . Therefore, the proof is finished by Theorem 2.3. ∎

Remark 3.8.

Yi [35] recently applied the results in [20] to give the period function for small-amplitude periodic solutions bifurcating from a weak focus of order one, whose first Lyapunov coefficient has nonzero real part. Following that, Yi [35] gave a criterion for the destabilization of the synchronous periodic solutions arising from a weak focus of order one. Here we consider the destabilization of the synchronous periodic solutions arising from a higher-order weak focus. This phenomenon is called degenerate Hopf bifurcation [11]. It is worth mentioning that [8, 17, 20] considered the period function of perturbed periodic solutions appearing in one-parameter Hopf bifurcation. However, there are two parameters involved in determining the related period function for a higher-order weak focus, the results in [8, 17, 20] are not applicable to this case.

3.3. Application to the two-patch Holling-Tanner model

In this section, we consider the two-patch Holling-Tanner model as illustration for our results. A similar argument can be also applied to explore the destabilization of synchronous periodic solutions for various patch models with two-dimensional kinetic systems. It is worth mentioning that Theorem 2.3 are applicable to patch models not only with two-dimensional kinetic systems but also with high-dimensional kinetic systems, e.g. epidemic models [26, 33], ecological systems [4, 23], chemical reaction models [10, 28], etc.

Consider a two-dimensional kinetic system

\displaystyle\begin{aligned} \frac{du}{dt}&=\dot{u}=ru\left(1-\frac{u}{K}\right)-q(u)v=:P(u,v),\\ \frac{dv}{dt}&=\dot{v}=sv\left(1-h\frac{v}{u}\right)=:Q(u,v),\end{aligned}

(3.35)

where $u$ and $v$ are the population densities of a prey and a predator, respectively. Here $r$ and $s$ indicate the growth rates of the prey and the predator respectively, $K$ measures the prey environmental carrying capacity in the absence of predation, and $h$ presents a measure of food quality. The functional response $q(u)$ is of Holling type II and has the form

q(u)=\frac{mu}{u+a},

where we require $a>0$ and $m>0$ . This is the classical Holling-Tanner model, which exhibits interesting oscillatory behaviors (see, for instance, [16, 22]).

To study spatial aspects, we consider a two-patch model with the kinetic system (3.35) on each patch and cross-diffusion-like couplings between the two patches. The dynamics are governed by

		$\displaystyle\frac{du_{j}}{dt}=\delta\sum_{i\in\Omega}\left(d_{11}(u_{i}-u_{j})+d_{12}(v_{i}-v_{j})\right)+P(u_{j},v_{j}),\ \ \ \ \$		$\displaystyle\ \ \ \ j\in\Omega:=\{1,2\},$		(3.36)
		$\displaystyle\frac{dv_{j}}{dt}=\delta\sum_{i\in\Omega}(d_{21}(u_{i}-u_{j})+d_{22}(v_{i}-v_{j}))+Q(u_{j},v_{j}),\ \ \ \ \$		$\displaystyle\ \ \ \ j\in\Omega:=\{1,2\},$		(3.36)

where $d_{ij}$ are the diffusion coefficients and $\delta$ is the coupling strength. Then the linearization of patch model (3.36) about a synchronous periodic solution $\mathcal{U}_{0}(t)=(\phi(t),\phi(t))$ is

		$\displaystyle\frac{dX_{j}}{dt}=\delta\sum_{i\in\Omega}\left(d_{11}(X_{i}-X_{j})+d_{12}(Y_{i}-Y_{j})\right)+P_{u}(\phi(t))X_{j}+P_{v}(\phi(t))Y_{j},$		$\displaystyle\ \ j\in\Omega:=\{1,2\},$		(3.37)
		$\displaystyle\frac{dY_{j}}{dt}=\delta\sum_{i\in\Omega}(d_{21}(X_{i}-X_{j})+d_{22}(Y_{i}-Y_{j}))+Q_{u}(\phi(t))X_{j}+Q_{v}(\phi(t))Y_{j},$		$\displaystyle\ \ j\in\Omega:=\{1,2\},$		(3.37)

In the following, we use two concrete examples to demonstrate the destabilization of synchronous periodic solutions for patch model (3.36). Examples 1 and 2 illustrate the cases of large- and small-amplitude bifurcating periodic solutions, respectively. Here we shall use Lyapunov exponents (see [18]) of the linearized systems to describe the destabilization. If the linearized systems have a positive Lyapunov exponent, then the corresponding synchronous periodic solutions become unstable.

Example 1. Consider the Holling-Tanner model (3.35) and fix the parameters as follows:

\displaystyle a=1,\ \ \ h=0.5,\ \ \ K=5,\ \ \ m=1,\ \ \ r=1,\ \ \ s=0.1.

Numerical simulation with MATLAB shows that the kinetic system (3.35) has a stable periodic solution $\varphi_{1}(t)$ for $t\in\mathbb{R}$ that surrounds a unstable focus. See Figure 1.

Refer to caption — Figure 1. A stable periodic solution (blue cycle) surrounds a unstable focus (red dot) in model (3.35), where $a=1$ , $h=0.5$ , $K=5$ , $m=1$ , $r=1$ and $s=0.1$ .

Set $\mathcal{U}_{1}(t):=(\varphi_{1}(t),\varphi_{1}(t))$ . When $\delta=0$ , patch model (3.36) is decoupled and the synchronous periodic solution $\mathcal{U}_{1}$ is stable with respect to patch model (3.36). With the aid of MATLAB, we obtain the results as follows:

(i)

Set $\delta=0.01$ , $d_{11}=d_{22}=1$ and $d_{12}=d_{21}=0$ . Numerical simulation shows that the largest Lyapunov exponent is zero. See Figure 2. This coincides with Proposition 3.3, which shows that the synchronous periodic solution $\mathcal{U}_{1}$ is still stable with respect to patch model (3.36) with the identical diffusion rates.
(ii)

Set $\delta=0.1$ , $d_{11}=d_{21}=d_{22}=1$ and $d_{12}=10$ . Numerical simulation shows that the largest Lyapunov exponent is about 0.0031. See Figure 2. This implies the destabilization of the synchronous periodic solution $\mathcal{U}_{1}$ can be induced by cross-diffusion-like couplings.

Example 2. Consider the Holling-Tanner model (3.35) and fix $K=m=r=1$ . By a direct computation, the positive equilibria of model (3.35) are determined by the roots of equation

\displaystyle\beta u^{2}+(a\beta-\beta+s)u-a\beta=0,\ \ \ \beta:=hs.

(3.38)

Clearly, the above equation has a unique positive root that is denoted by $u_{*}$ . So $E_{*}:=(u_{*},u_{*}/h)$ is the unique positive equilibrium of model (3.35). Furthermore, if the following equations hold:

	$\displaystyle 2u_{}^{2}+(a+s-1)u_{}+as\!\!\!$	$\displaystyle=$	$\displaystyle\!\!\!0,$		(3.39)
	$\displaystyle-(3+a)u_{}^{3}-6au_{}(a+u_{*})+a^{2}(1-a)\!\!\!$	$\displaystyle=$	$\displaystyle\!\!\!0,$		(3.40)

then by [16, Theorem 3.2] and statements in [16, p.161], this equilibrium $(u_{*},u_{*}/h)$ is a weak focus of multiplicity two and asymptotically stable. Consequently, the corresponding first- and second-order Lyapunov coefficients $C_{1}$ and $C_{2}$ satisfy

{\rm Re}(C_{1})=0,\ \ \ {\rm Re}(C_{2})<0.

Set $s=0.1$ and let $u$ in (3.38) be replaced by $u_{*}$ . Solving (3.38), (3.39) and (3.40) yields a unique positive solution $a\approx 0.336238$ , $h\approx 0.222132$ and $u_{*}\approx 0.085693$ . Now we fix these $a$ and $h$ , i.e., $a\approx 0.336238$ and $h\approx 0.222132$ , and vary $s$ . At this equilibrium $E_{*}$ , we can obtain the Jacobian matrix $J(s)=(J_{ij}(s))$ of model (3.35) satisfies

\displaystyle J_{11}(0.1)\approx 0.1,\ \ \ J_{12}(0.1)\approx-0.203097,\ \ \ J_{21}(0.1)\approx 0.450183,\ \ \ J_{22}(0.1)\approx-0.1,

and has a pair of purely imaginary eigenvalue $\lambda_{\pm}\approx\pm 0.285361\sqrt{-1}$ . By the formula of the first-order Lyapunov coefficient $C_{1}$ in page 90 of [20] (see also Appendix A in [35]), we can compute that the imaginary part ${\rm Im}(C_{1})$ of the first-order Lyapunov coefficient $C_{1}$ is

\displaystyle{\rm Im}(C_{1})\approx-1.872272.

Note that the trace $\mathcal{T}(s)$ of $J(s)$ is $(0.1-s)$ . Then by Proposition 3.6, there exists a sufficiently small $s_{0}>0$ such that for $0<0.1-s<s_{0}$ , model (3.35) has a stable periodic solution $\varphi_{2}(t)$ arising from $E_{*}$ . See Figure 3.

Set $\mathcal{U}_{2}(t):=(\varphi_{2}(t),\varphi_{2}(t))$ . When $\delta=0$ , patch model (3.36) is decoupled and the synchronous periodic solution $\mathcal{U}_{2}$ is stable with respect to patch model (3.36). If the diffusion rates $d_{ij}$ satisfy

\displaystyle d_{22}J_{11}(0.1)+d_{11}J_{22}(0.1)-d_{12}J_{21}(0.1)-d_{21}J_{12}(0.1)<0,

(3.41)

then by Proposition 3.7, the synchronous periodic solution $\mathcal{U}_{2}$ becomes unstable. With the aid of MATLAB, we obtain the results as follows:

(i)

Set $\delta=0.01$ , $d_{11}=1$ , $d_{12}=1$ , $d_{21}=-100$ and $d_{22}=5$ . Numerical simulation shows that the largest Lyapunov exponent is zero. See Figure 4. Then the synchronous periodic solution $\mathcal{U}_{2}$ is still stable.
(ii)

Set $\delta=0.01$ , $d_{11}=1$ , $d_{12}=1$ , $d_{21}=100$ and $d_{22}=5$ . Numerical simulation shows that the largest Lyapunov exponent is about 0.0079. See Figure 4. Then the synchronous periodic solution $\mathcal{U}_{2}$ becomes unstable. This coincides with Proposition 3.7, which shows that stable synchronous periodic solution $\mathcal{U}_{1}$ becomes unstable with respect to patch model (3.36) if the condition (3.41) holds.

Appendix A. perturbation of eigenvalues for matrices

In this appendix, we present the perturbation theory for matrices developed in our recent work [7]. For each $\zeta=(\zeta_{1},...,\zeta_{N})^{T}$ in $\mathbb{C}^{N}$ or $\mathbb{R}^{N}$ , set

\displaystyle|\zeta|:=\max\{|\zeta_{1}|,|\zeta_{2}|,...,|\zeta_{N}|\}.

Let $\|\cdot\|$ denote the norm of a matrix, i.e., the maximum row sum of the absolute values of the entries.

In order to give our criterion for the destabilization of synchronous periodic solutions, it is necessary to give the asymptotic expressions of the Floquet spectra. The argument is based on two perturbation results which were proved in our recent work [7].

Consider a matrix $A\in\mathbb{R}^{N\times N}$ ( $N\geq 2$ ) as follows:

\displaystyle A=\left(\begin{array}[]{cc}A_{1}&0\\ 0&A_{2}\end{array}\right),

where $A_{1}={\rm diag}(a_{1},a_{2},...,a_{N_{1}})\in\mathbb{R}^{N_{1}\times N_{1}}$ and $A_{2}\in\mathbb{R}^{N_{2}\times N_{2}}$ for integers $N_{1}$ and $N_{2}$ with $0<N_{1}=N-N_{2}<N$ . Additionally, the spectra $\sigma(A_{1})$ and $\sigma(A_{2})$ of $A_{1}$ and $A_{2}$ satisfy the following assumption:

•

The spectra $\sigma(A_{1})$ and $\sigma(A_{2})$ are separated by a simple closed positively oriented cycle $\Gamma$ in the complex plane, and $\sigma(A_{1})$ lies in the interior of the closed cycle $\Gamma$ .

For a sufficiently small constant $\delta_{0}>0$ , let $B:(-\delta_{0},\delta_{0})\to\mathbb{R}^{N\times N}$ denote a matrix function that is analytic in $\delta\in(-\delta_{0},\delta_{0})$ and satisfies that $\|B(\delta)\|=O(|\delta|)$ for sufficiently small $|\delta|$ . Consider a perturbation of $A$ as follows:

\displaystyle A(\delta):=A+B(\delta),\ \ \ \ \ \delta\in(-\delta_{0},\delta_{0}).

(A.1)

By the classical perturbation theory of eigenvalues for matrices (see, for instance, [24, pp. 63-64]), the eigenvalues of $A(\delta)$ are continuous in $\delta$ . By continuity, we can choose sufficiently small $\delta_{0}>0$ such that all eigenvalues of $A(\delta)$ perturbed from $\sigma(A_{1})$ lie in the interior of the closed cycle $\Gamma$ and all others lie outside the domain surrounded by $\Gamma$ in the complex plane. Following that, we can define a family of parametric projections $\mathcal{P}(\delta)$ by

\mathcal{P}(\delta)=\frac{1}{2\pi{\bf i}}\oint_{\Gamma}(\lambda I_{N}-A(\delta))\,d\lambda,\ \ \ \ \ \delta\in(-\delta_{0},\delta_{0}).

Concerning these parametric projections $\mathcal{P}(\delta)$ , we have the next lemma.

Lemma A.1.

[7, Lemma A.1] There exists a sufficiently small constant $\hat{\delta}_{0}$ with $0<\hat{\delta}_{0}\leq\delta_{0}$ such that the following statements hold:

(i)

$\mathcal{P}(\delta)$ is analytic in the interval $(-\hat{\delta}_{0},\hat{\delta}_{0})$ .
(ii)

Let the ranges of $\mathcal{P}(0)$ and $I-\mathcal{P}(0)$ be spanned by $\{\xi_{j}:j=1,...,N_{1}\}$ and $\{\xi_{j}:j=N_{1}+1,...,N\}$ , respectively. Then there exists an operator-valued function $U(\cdot):(-\hat{\delta}_{0},\hat{\delta}_{0})\to\mathbb{R}^{N\times N}$ which is analytic and invertible, such that for each $\delta\in(-\hat{\delta}_{0},\hat{\delta}_{0})$ , the sets $\{U(\delta)\xi_{j}:j=1,...,N_{1}\}$ and $\{U(\delta)\xi_{j}:j=N_{1}+1,...,N\}$ form the bases of the ranges of $\mathcal{P}(\delta)$ and $I-\mathcal{P}(\delta)$ , respectively.

We continue to study a special case of the perturbation (A.1), i.e., the perturbation $A(\delta)$ is in the form

\displaystyle A(\delta)=\left(\begin{array}[]{cc}A_{1}(\delta)&0\\ 0&A_{2}(\delta)\end{array}\right)+A_{3}(\delta),\ \ \ \ \ \delta\in(-\delta_{0},\delta_{0}),

(A.2)

where $A_{1}(\delta)\in\mathbb{R}^{N_{1}\times N_{1}}$ , $A_{2}(\delta)\in\mathbb{R}^{N_{2}\times N_{2}}$ and $A_{3}(\delta)\in\mathbb{R}^{N\times N}$ are continuous and satisfy

\displaystyle A_{1}(\delta)=A_{1}+O(|\delta|),\ \ \ A_{2}(\delta)=A_{2}+O(|\delta|),\ \ A_{3}(\delta)=O(|\delta|^{2})

for sufficiently small $|\delta|$ . Then we have the following result.

Lemma A.2.

[7, Lemma A.2] Consider the perturbation $A(\delta)$ of the form (A.2). Suppose that

\displaystyle A_{1}(\delta)=\left(\begin{array}[]{cccc}a_{1}+b_{1}\delta&0&\cdot\cdot\cdot&0\\ 0&a_{2}+b_{2}\delta&\cdot\cdot\cdot&0\\ \cdot\cdot\cdot&\cdot\cdot\cdot&\cdot\cdot\cdot&\cdot\cdot\cdot\\ 0&0&\cdot\cdot\cdot&a_{N_{1}}+b_{N_{1}}\delta\end{array}\right),\ \ \ \delta\in(-\delta_{0},\delta_{0}),

where $b_{j}$ ( $j=1,2,...,N_{1}$ ) are real constants. Then for each $j=1,2,...,N_{1}$ , the eigenvalue $\lambda_{j}(\delta)$ perturbed from $\lambda_{j}=a_{j}$ has the asymptotic expression as follows:

\displaystyle\lambda_{j}(\delta)=a_{j}+b_{j}\delta+O(\delta^{2})

for sufficiently small $|\delta|\geq 0$ .

Declarations

Ethical Approval

Not applicable.

Competing interests

There are no financial or non-financial competing interests.

Authors’ contributions

S. Chen and J. Huang wrote the paper. All authors read and approved the manuscript.

Funding

This work was partly supported by the National Natural Science Foundation of China (Grant No. 12101253, 12231008) and the Scientific Research Foundation of CCNU (Grant No. 31101222044).

Availability of data and materials

Data sharing not applicable to this article as no datasets were generated or analysed during the current study.

References

[1] L. Allen, B. Bolker, Y. Lou, A. Nevai, Asymptotic profiles of the steady states for an SIS epidemic patch model, SIAM J. Appl. Math. 67 (2007), 1283–1309.
[2] A. Andronov, E. Leontovich, I. Gordon, A. Maier, Theory of Bifurcation of Dynamic Systems on a Plane, Israel Program for Sci. Transl., Wiley, New York, 1973.
[3] R. Arumugam, F. Guichard, F. Lutscher, Persistence and extinction dynamics driven by the rate of environmental change in a predator-prey metacommunity, Theoretical Ecology 13 (2020), 629–643.
[4] R. Arumugam, F. Lutscher, F. Guichard, Tracking unstable states: ecosystem dynamics in a changing world, Oikos 130 (2021), 525–540.
[5] E. Coddington, N. Levinsion, Theory of Ordinary Differential Equations, McGraw-Hill, New-York, 1955.
[6] S. Chen, J. Duan, Instability of small-amplitude periodic waves from fold-Hopf bifurcation, J. Math. Phys. 63 (2022), 112702, 18 pp.
[7] S. Chen, J. Huang, Destabilization of synchronous periodic solutions for patch models, J. Differential Equations 364 (2023), 378–411.
[8] S. Chen, J. Huang, Periodic traveling waves with large speed, Z. Angew. Math. Phys. 74 (2023), No. 102, 28 pp.
[9] S.-N. Chow, C. Li, D. Wang, Normal forms and bifurcation of planar vector fields, Cambridge University Press, Cambridge, 1994.
[10] M. Dolnik, I. Epstein, A coupled chemical burster: The chlorine dioxide-iodide reaction in two flow reactors, J. Chem. Phys. 98 (1993), 1149–1155.
[11] W. Farr, C. Li, I. Labouriau, W. Langford, Degenerate Hopf bifurcation formulas and Hilbert’s 16th problem, SIAM J. Math. Anal. 20 (1989), 13–30.
[12] D. Gao, Transmission Dynamics of Some Epidemiological Patch Models, PhD Thesis, University of Miami, 2012.
[13] D. Gao, Travel frequency and infectious diseases, SIAM J. Appl. Math. 79 (2019), 1581–1606.
[14] D. Gao, How does dispersal affect the infection size? SIAM J. Appl. Math. 80 (2020), 2144–2169.
[15] D. Gao, Y. Lou, Impact of state-dependent dispersal on disease prevalence, J. Nonlinear Sci. 31 (2021), 41 pp.
[16] A. Gasull, R.E. Kooij, J. Torregrosa, Limit cycles in the Holling-Tanner model, Publ. Mat. 41 (1997), 149–167.
[17] A. Gasull, V. Mañosa, J. Villadelprat, On the period of the limit cycles appearing in one-parameter bifurcations, J. Differential Equations 213 (2005), 255–288.
[18] J. Guckenheimer, P. Holmes, Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields, Appl. Math. Sci., Vol. 42, Springer, New York, 1983.
[19] Jack K. Hale, Ordinary Differential Equations, Dover Publ., New York, 1980.
[20] B. Hassard, N. Kazarinoff, Y. Wan, Theory and Application of Hopf Bifurcation, Cambridge University Press, Cambridge, 1981.
[21] D. Henry, Geometric Theory of Semilinear Parabolic Equations, Lecture Notes in Mathematics, Vol. 840, Springer, New York, 1981.
[22] S. Hsu, T. Huang, Global stability for a class of predator-prey systems, SIAM J. Appl. Math. 55 (1995), 763–783.
[23] H. Jiang, K.-Y. Lam, Y. Lou, Three-patch models for the evolution of dispersal in advective environments: varying drift and network topology, Bull. Math. Biol. 83 (2021), 1–46.
[24] T. Kato, Perturbation Theory for Linear Operators, Springer, 1980.
[25] I. Lengyel, I. Epstein, Diffusion-induced instability in chemically reacting systems: steady state multiplicity, oscillation, and chaos, Chaos 1(1) (1991), 69–76.
[26] M. Lu, D. Gao, J. Huang, H. Wang, Relative prevalence-based dispersal in an epidemic patch model, J. Math. Biol. 86 (2023), 35 pp.
[27] K. Maginu, Stability of spatially homogeneous periodic solutions of reaction-diffusion equations, J. Differential Equations 31 (1979), 130–138.
[28] P. Moore, W. Horsthemke, Localized patterns in homogeneous networks of diffusively coupled reactors, Phys. D 206 (2005), 121–144.
[29] S. Ruan, Diffusion-driven instability in the Gierer-Meinhardt model of morphogenesis, Nat. Resour. Model. 11 (1998) 131–142.
[30] A. M. Turing, The chemical basis of morphogenesis, Philos. Trans. R. Soc. Lond. B 237 (1952), 37–72.
[31] M. Urabe, Infinitesimal deformation of cycles, J. Sci. Hiroshima Univ. Ser. A 18 (1954), 37–53.
[32] P. van den Driessche, Spatial structure: patch models, in: Fred Brauer, Pauline van den Driessche, Jianhong Wu (Eds.), Mathematical Epidemiology, Springer, Berlin, 2018, pp. 170-189.
[33] W. Wang, G. Mulone, Threshold of disease transmission in a patch environment, J. Math. Anal. Appl. 285 (2003), 321–335.
[34] W. Wang, X.-Q. Zhao, An epidemic model in a patchy environment, Math. Biosci. 190 (2004), 97–112.
[35] F. Yi, Turing instability of the periodic solutions for reaction-diffusion systems with cross-diffusion and the patch model with cross-diffusion-like coupling, J. Differential Equations 281 (2021), 379–410.
[36] Y. Q. Ye et al., Theory of Limit Cycles, Transl. of Math. Monogr., Vol. 66, Amer. Math. Soc., Providence, RI, 1986.

Destabilization of synchronous periodic solutions for patch models: a criterion by period functions

Abstract.

Key words and phrases:

2020 Mathematics Subject Classification:

1. Introduction

2. Destabilization of synchronous periodic solutions

Lemma 2.1.

Proof.

Lemma 2.2.

Proof.

Theorem 2.3.

Proof.

Remark 2.4.

3. Application to nn-patch models with two-dimensional kinetic system

3.1. Application to large-amplitude bifurcating periodic solutions

Lemma 3.1.

Proof.

Proposition 3.2.

Proof.

Proposition 3.3.

Proof.

3.2. Application to small-amplitude bifurcating periodic solutions

Lemma 3.4.

Proof.

Lemma 3.5.

Proof.

Proposition 3.6.

Proof.

Proposition 3.7.

Proof.

Remark 3.8.

3.3. Application to the two-patch Holling-Tanner model

Appendix A. perturbation of eigenvalues for matrices

Lemma A.1.

Lemma A.2.

Declarations

References

3. Application to $n$ -patch models with two-dimensional kinetic system