# Dynamical Analysis of a Delayed Reaction-Diffusion Predator-Prey System

Authors: Yanuo Zhu; Yongli Cai; Shuling Yan; Weiming Wang
Journal: Abstract and Applied Analysis
 (2012)

---

## Abstract

This work deals with the analysis of a delayed diffusive predator-prey system under Neumann boundary conditions. The dynamics are investigated in terms of the stability of the nonnegative equilibria and the existence of Hopf bifurcation by analyzing the characteristic equations. The direction of Hopf bifurcation and the stability of bifurcating periodic solution are also discussed by employing the normal form theory and the center manifold reduction. Furthermore, we prove that the positive equilibrium is asymptotically stable when the delay is less than a certain critical value and unstable when the delay is greater than the critical value.


---

## Body

## 1. Introduction

The study on the dynamics of predator-prey systems is one of the dominant subjects in ecology and mathematical ecology due to its universal existence and importance [1]. A prototypical predator-prey interaction model is of the following form:
(1.1)dudt=a(u)-f(u)g(v),dvdt=σf(u)g(v)-z(v),
where u(t) and v(t) are the densities of the prey and predator at time t>0, respectively.Furthermore, the functiona(u) is growth rate of the prey in the absence of predation, which is given by
(1.2)a(u)=αumin{1,K-uK-ε},α>0,ε≥0,K>0.
If ε=0, this reduces to the traditional logistic form a(u)=αu(1-u/K), see [2] and the references therein. Here, the parameter α stands for the specific growth rate of the prey u, and K for carrying capacity of the prey in the absence of predators.The productf(u)g(v) gives the rate at which prey is consumed, and f(u)g(v)/v is termed as the functional response [3]. In particular, these functions can be defined by
(1.3)f(u)=cu,g(v)=vmv+1,c>0,m≥0,
where c denotes the capture rate, and m the half capturing saturation constant.The proportionality constantσ is the rate of prey consumption. And the function z(v) is given by
(1.4)z(v)=γv+lv2,γ>0,l≥0,
where γ denotes the natural death rate of the predators, and l>0 can be used to model predator in traspecific competition that is not direct competition for food, such as some type of territoriality, see [2]. In this paper, we discuss the case of l=0, which is used in a much more traditional case. Based on the above discussions, we can obtain the following model:
(1.5)dudt=αu(1-uK)-cuvmv+1,dvdt=v(-γ+cσumv+1).Settingβ=α/K,  b=cσ, model (1.5) leads to the following dimensionless equation:
(1.6)dudt=u(α-βu)-cuvmv+1,dvdt=v(-γ+bumv+1),
where b>0 denotes conversion rate.In recent years, the models involving time delay and spatial diffusion have been extensively studied by many authors and many interesting results have been obtained, including the stability, the existence of Hopf bifurcation, and direction of bifurcating periodic solutions, see [1, 4–18]. In this paper, we mainly focus on the effects of both spatial diffusion and time delay factors on system (1.6) with Neumann boundary conditions as follows:
(1.7)∂u∂t=u(α-βu)-cuvmv+1+d1Δu,x∈Ω,t>0,∂v∂t=v(-γ+bu(t-τ)mv(t-τ)+1)+d2Δv,x∈Ω,t>0,∂u∂ν=∂v∂ν=0,x∈∂Ω,t>0,u(x,0)=u0(x)≥0,v(x,0)=v0(x)≥0,x∈Ω,
where Ω is a bounded open domain in ℝ with a smooth boundary ∂Ω, and Δ=∂2/∂x2 denotes the Laplacian operator in ℝ. d1>0 and d2>0 denote the diffusion coefficients of the prey u and predator v, respectively. ν is the outward unit normal vector on ∂Ω. τ>0 can be regarded as the gestation of the predator. System (1.7) includes not only the dispersal processes, but also some of the past states of the system.Throughout this paper, we restrict ourselves to the one-dimensional spatial domainΩ=(0,π) for the sake of convenience.The remaining parts of the paper are structured in the following way. In Section2, we analyze the distribution of the roots of the characteristic equation and give various conditions on the stability of a unique positive equilibrium and the existence of Hopf bifurcation with time delay. In Section 3, applying the normal form theory [19, 20] and the center manifold reduction of partial functional differential equations [21], we derive the explicit algorithm in order to determine the direction of the Hopf bifurcation, the stability, and other properties on bifurcating periodic solutions. Finally, a brief discussion is given.


## 2. Stability of Positive Equilibrium and Existence of Hopf Bifurcation

In this section, by analyzing the associated characteristic equation of system (1.7) at the positive equilibrium, we investigate the stability of the positive equilibria of system (1.7).It is straightforward to see that system (1.7) has the following two boundary equilibria:(i)
E1=(0,0) (total extinct) which is saddle point, hence it is unstable;(ii)
E2=(α/β,0) (extinct of the predator) which is saddle point if bα>βγ, or stable if bα<βγ.To find the positive equilibrium, we set(2.1)u(α-βu)-cuvmv+1=0,v(-γ+bumv+1)=0,
which yields
(2.2)mbβu2+(c-mα)bu-cγ=0.Obviously, system (1.7) has a unique positive equilibrium E*=(u*,v*) with bα>βγ, where
(2.3)u*=mbα-bc+4mbcβγ+b2(mα-c)22mbβ,v*=bu*-γmγ.Setu¯=u-u*, v¯=v-v* and drop the bars for simplicity of notations, then system (1.7) can be transformed into the following equivalent system:
(2.4)∂u∂t=(u+u*)(α-β(u+u*))-c(u+u*)vm(v+v*)+1+d1Δu,x∈Ω,t>0,∂v∂t=(v+v*)(-γ+b(u(t-τ)+u*)m(v(t-τ)+v*)+1)+d2Δv,x∈Ω,t>0,∂u∂ν=∂v∂ν=0,x∈∂Ω,t>0,u(x,0)=u0(x)-u*≥0,v(x,0)=v0(x)-v*≥0,x∈Ω.Assume thatu0(x),v0(x)∈C([-τ,0];X) and X is defined by
(2.5)X={(u,v):u,v∈W2,2(Ω):∂u∂ν=∂v∂ν=0,x∈∂Ω}
with the inner product 〈·,·〉.Denote(u(t),v(t))=(u(t,x),v(t,x)) and U(t)=(u(t),v(t))T. Then system (2.4) can be rewritten as an abstract differential equation in the phase space C([-τ,0];X) as follows:
(2.6)∂U(t)∂t=dΔU(t)+L(Ut)+F(Ut),
where d=diag(d1,d2), Ut(θ)=U(t+θ), -τ≤θ≤0, and L:C([-τ,0];X)→X, F:C([-τ,0];X)→X are given by
(2.7)L(ϕ)=(-βu*ϕ1(0)-cu*(mv*+1)2ϕ2(0)bv*(mv*+1)ϕ1(-τ)-mbu*v*(mv*+1)2ϕ2(-τ)),(2.8)F(ϕ)=(-βϕ12(0)-c(mv*+1)2ϕ1(0)ϕ2(0)+cmu*(cm+1)3ϕ22(0)(ϕ2(0)+mv*mv*+1ϕ2(-τ))(bmv*+1ϕ1(-τ)-bmu*(mv*+1)2ϕ2(-τ))),
respectively, where ϕ(θ)=Ut(θ),  ϕ=(ϕ1,ϕ2)T∈C([-τ,0];X).The linearization of (2.6) is given by
(2.9)∂U(t)∂t=dΔU(t)+L(Ut),
and its characteristic equation is
(2.10)λy-dΔy-L(eλ·y)=0,
where y∈dom(Δ) and y≠0, dom(Δ)⊂X.It is well known that the eigenvalue problem(2.11)-Δψ=μψ,x∈(0,π),∂ψ∂x|x=0=∂ψ∂x|x=π=0,
has eigenvalues 0=μ0≤μ1≤μ2≤⋯≤μn≤μn+1≤⋯, with the corresponding eigenfunctions ψn(x). Substituting
(2.12)y=∑n=0∞‍ψn(x)(y1ny2n)
into characteristic equation (2.10), we obtain
(2.13)(-βu*-d1μn-cu*(mv*+1)2bv*mv*+1e-λτ-bmu*v*(mv*+1)2e-λτ-d2μn)(y1ny2n)=λ(y1ny2n).Hence, we can conclude that the characteristic equation (2.10) is equivalent to the sequence of the following characteristic equations:
(2.14)λ2+(An+De-λτ)λ+Bn+Cne-λτ=0(n=0,1,2,…),
where
(2.15)An=βu*+(d1+d2)μn,Bn=βu*d2μn+d1d2μn2,Cn=bd1mμnβu*v*(mv*+1)2+bmβ(u*)2v*(mv*+1)2+bcu*v*(mv*+1)3,D=bmu*v*(mv*+1)2.
The stability of the positive equilibrium E*=(u*,v*) can be determined by the distribution of the roots of (2.14) (n=0,1,2,…), that is, the equilibrium E*=(u*,v*) is locally asymptotically stable if all the roots of (2.14) (n=0,1,2,…) have negative real parts. Note that λ=0 is not a root of (2.14) for any n=0,1,2,…. Next, we analyze the behaviour of system (1.7) in two situations:  with/without delay effect.


### 2.1. Caseτ=0

Equation (2.14) with τ=0 is equivalent to the following quadratic equation:
(2.16)λ2+(An+D)λ+Bn+Cn=0,
where An,  Bn,  Cn,  and  D are defined as (2.15).Letλ1 and λ2 be two roots of (2.16), then for any n=0,1,2,…, we have
(2.17)λ1+λ2=-(An+D)<0,λ1λ2=Bn+Cn>0.
Then we can get the following theorem.Theorem 2.1.
Ifbα>βγ holds, the positive equilibrium E*=(u*,v*) of system (1.7) with τ=0 is asymptotically stable.In the following, we prove thatE*=(u*,v*) of system (1.7) is globally stable with τ=0.Theorem 2.2.
Ifbα>βγ holds, the positive equilibrium E*=(u*,v*) of system (1.7) with τ=0 is globally asymptotically stable.Proof.
To prove our statement, we need to construct a Lyapunov function. To this end, we define(2.18)V(t)=∫Ω‍V1(u,v)dx,
where
(2.19)V1(u,v)=∫u*u‍ξ-u*ξdξ+cb(1+mv*)∫v*v‍η-v*ηdη.
We claim that V1(u,v) is positive definite. In fact, set
(2.20)∂V1∂u=1-u*u=0,∂V1∂v=cb(1+mv*)(1-v*v)=0,
we can obtain (u,v)=(u*,v*). And the Hessian Matrix at (u*,v*) is given by
(2.21)H(E)|(u*,v*)=(1u*00cbv*(1+mv*)),
Hence H(E)|(u*,v*) is positive definite, which follows that
(2.22)min(V1(u,v))=V1(u*,v*)=0.
The time derivative ofV along the solutions of system (1.7) with τ=0, we have
(2.23)dVdt=∫Ω‍(u-u*u∂u∂t+cb(1+mv*)(v-v*)v∂v∂t)dx=∫Ω‍dV1dtdx+∫Ω‍{(d1Δu)∂V∂u+(d2Δv)∂V∂v}dx=∫Ω‍(-β(u-u*)2-cmu*(mv+1)(mv*+1)2(v-v*)2)dx-d1u*u2∫Ω‍(∂u∂x)2dx-d2v*v2∫Ω‍(∂v∂x)2dx,
then, we obtain dV/dt<0.
It is enough to see thatdV/dt satisfies Lyapunov’s asymptotic stability theorem, hence the positive equilibrium E*=(u*,v*) of system (1.7) with τ=0 is globally asymptotically stable.


### 2.2. Caseτ≠0

In the following, we prove the stability of the positive equilibriumE*=(u*,v*) of system (1.7) and the existence of Hopf bifurcation at the positive equilibrium E*=(u*,v*).Theorem 2.3.
Assume that>βγ, An2-D2-2Bn2<0 and Bn2-Cn2>0 hold, then the positive equilibrium E*=(u*,v*) is asymptotically stable for all τ≥0.Proof.
Letλ=μ(τ)+iω(τ) be a root of the characteristic equation (2.14), then we have
(2.24)μ2-ω2+Anμ+Bn+(Dωsinωτ+(Cn+Dμ)cosωτ)e-μτ=0,2μω+Anω-((Cn+D)sinωτ+Dωcosωτ)e-μτ=0,
where λ, μ, ω are functions of τ. A necessary condition for the stability of E*=(u*,v*) is that the characteristic equation has a purely imaginary solution λ=iω. Let μ(τ)=0 and ω(τ)≠0, then we can reduce (2.24) to
(2.25)-ω2+Dωsinωτ+Cncosωτ+Bn=0,Anω+Dωcosωτ-Cnsinωτ=0,
which lead to
(2.26)ω4+(An2-D2-2Bn)ω2+(Bn2-Cn2)=0.
Since Bn2-Cn2>0 and (An2-D2-2Bn)<0, these imply that (2.26) has no positive roots, that is, all roots of (2.14) have negative real parts.Theorem 2.4.
Assume thatbα>βγ, Bn2-Cn2<0, then there exists a sequence
(2.27)τj0=τ00+2jπω0(j=0,1,2,…),
where  τ00 is defined as (2.32), such that, for system (1.7), the following statements are true. (i)
Ifτ∈[0,τ00), then the positive equilibrium E*=(u*,v*) is asymptotically stable.(ii)
Ifτ>τ00, then the positive equilibrium E*=(u*,v*) is unstable.(iii)
τ=τ0j(n=0,1,2,…) are Hopf bifurcation values of system (1.7) and these Hopf bifurcations are all spatially homogeneous.Proof.
Letλ=iω(τ) be a root of the characteristic equation (2.14). By the same way in Theorem 2.3, then ω satisfies the following equation:
(2.28)ω4+(An2-D2-2Bn)ω2+(Bn2-Cn2)=0.
Since Bn2-Cn2<0, (2.28) has a unique positive root ω0 satisfying
(2.29)ω0=-(An2-D2-2Bn)+(An2-D2-2Bn)2-4(Bn2-Cn2)2
and from (2.25) we obtain
(2.30)sinω0τ=Dω03+AnCnω0-BnDω0D2ω02+Cn2≜G(ω0),cosω0τ=Cnω02-BnCn-AnDω02D2ω02+Cn2≜E(ω0),
then (2.16) has one imaginary root iω0 when
(2.31)τj0=τ00+2jπω0,(j=0,1,2,…),
where τ00 satisfies
(2.32)τ00={arc cos(E(ω0))ωn,ifG(ω0)>0,arc cos(2π-E(ω0))ω0,ifG(ω0)<0,Letλ(τ)=α(τ)+iω(τ) be the root of (2.14) satisfying α(τj0)=0 and ω(τj0)=ω0 when τ is close to τj0. Now by some simple calculations we obtain
(2.33)[dαdτ|τ=τj0]-1=Re[(2λ+An+De-λτ-(Cn+Dλ)τe-λτ(Cn+Dλ)λe-λτ)|τ=τj0]=2Cncosω0τ+Dω0sinω0τCn2+D2ω02-D2ωn2Cn2ω02+D2ω04+Anω0(-Dωncosω0τ+Cnsinω0τ)Cn2ω02+D2ω04,
since ω02=Dωnsinω0τ+Cncosω0τ and Anω0=-Dω0cosω0τ+Cnsinω0τ, and from the expression of ω02 in (2.29), we immediately see that
(2.34)[dαdτ|τ=τj0]-1=2ω02Cn2+D2ω02+An2-D2Cn2+D2ω02=(D2-An2)2+4Cn2Cn2+D2ω02.
Therefore, sign[dα/dτ|τ=τj0]-1=1, that is, dα/dτ|τ=τj0>0. This implies that all the roots that cross the imaginary axis at iω cross from left to right as τ increases.Hence the transversality condition holds and accordingly Hopf bifurcation occurs atτ=τ00, and τ=τj0(j=0,1,2,…) are Hopf bifurcation values of system (1.7) and these Hopf bifurcations are all spatially homogeneous. This completes the proof.


## 2.1. Caseτ=0

Equation (2.14) with τ=0 is equivalent to the following quadratic equation:
(2.16)λ2+(An+D)λ+Bn+Cn=0,
where An,  Bn,  Cn,  and  D are defined as (2.15).Letλ1 and λ2 be two roots of (2.16), then for any n=0,1,2,…, we have
(2.17)λ1+λ2=-(An+D)<0,λ1λ2=Bn+Cn>0.
Then we can get the following theorem.Theorem 2.1.
Ifbα>βγ holds, the positive equilibrium E*=(u*,v*) of system (1.7) with τ=0 is asymptotically stable.In the following, we prove thatE*=(u*,v*) of system (1.7) is globally stable with τ=0.Theorem 2.2.
Ifbα>βγ holds, the positive equilibrium E*=(u*,v*) of system (1.7) with τ=0 is globally asymptotically stable.Proof.
To prove our statement, we need to construct a Lyapunov function. To this end, we define(2.18)V(t)=∫Ω‍V1(u,v)dx,
where
(2.19)V1(u,v)=∫u*u‍ξ-u*ξdξ+cb(1+mv*)∫v*v‍η-v*ηdη.
We claim that V1(u,v) is positive definite. In fact, set
(2.20)∂V1∂u=1-u*u=0,∂V1∂v=cb(1+mv*)(1-v*v)=0,
we can obtain (u,v)=(u*,v*). And the Hessian Matrix at (u*,v*) is given by
(2.21)H(E)|(u*,v*)=(1u*00cbv*(1+mv*)),
Hence H(E)|(u*,v*) is positive definite, which follows that
(2.22)min(V1(u,v))=V1(u*,v*)=0.
The time derivative ofV along the solutions of system (1.7) with τ=0, we have
(2.23)dVdt=∫Ω‍(u-u*u∂u∂t+cb(1+mv*)(v-v*)v∂v∂t)dx=∫Ω‍dV1dtdx+∫Ω‍{(d1Δu)∂V∂u+(d2Δv)∂V∂v}dx=∫Ω‍(-β(u-u*)2-cmu*(mv+1)(mv*+1)2(v-v*)2)dx-d1u*u2∫Ω‍(∂u∂x)2dx-d2v*v2∫Ω‍(∂v∂x)2dx,
then, we obtain dV/dt<0.
It is enough to see thatdV/dt satisfies Lyapunov’s asymptotic stability theorem, hence the positive equilibrium E*=(u*,v*) of system (1.7) with τ=0 is globally asymptotically stable.


## 2.2. Caseτ≠0

In the following, we prove the stability of the positive equilibriumE*=(u*,v*) of system (1.7) and the existence of Hopf bifurcation at the positive equilibrium E*=(u*,v*).Theorem 2.3.
Assume that>βγ, An2-D2-2Bn2<0 and Bn2-Cn2>0 hold, then the positive equilibrium E*=(u*,v*) is asymptotically stable for all τ≥0.Proof.
Letλ=μ(τ)+iω(τ) be a root of the characteristic equation (2.14), then we have
(2.24)μ2-ω2+Anμ+Bn+(Dωsinωτ+(Cn+Dμ)cosωτ)e-μτ=0,2μω+Anω-((Cn+D)sinωτ+Dωcosωτ)e-μτ=0,
where λ, μ, ω are functions of τ. A necessary condition for the stability of E*=(u*,v*) is that the characteristic equation has a purely imaginary solution λ=iω. Let μ(τ)=0 and ω(τ)≠0, then we can reduce (2.24) to
(2.25)-ω2+Dωsinωτ+Cncosωτ+Bn=0,Anω+Dωcosωτ-Cnsinωτ=0,
which lead to
(2.26)ω4+(An2-D2-2Bn)ω2+(Bn2-Cn2)=0.
Since Bn2-Cn2>0 and (An2-D2-2Bn)<0, these imply that (2.26) has no positive roots, that is, all roots of (2.14) have negative real parts.Theorem 2.4.
Assume thatbα>βγ, Bn2-Cn2<0, then there exists a sequence
(2.27)τj0=τ00+2jπω0(j=0,1,2,…),
where  τ00 is defined as (2.32), such that, for system (1.7), the following statements are true. (i)
Ifτ∈[0,τ00), then the positive equilibrium E*=(u*,v*) is asymptotically stable.(ii)
Ifτ>τ00, then the positive equilibrium E*=(u*,v*) is unstable.(iii)
τ=τ0j(n=0,1,2,…) are Hopf bifurcation values of system (1.7) and these Hopf bifurcations are all spatially homogeneous.Proof.
Letλ=iω(τ) be a root of the characteristic equation (2.14). By the same way in Theorem 2.3, then ω satisfies the following equation:
(2.28)ω4+(An2-D2-2Bn)ω2+(Bn2-Cn2)=0.
Since Bn2-Cn2<0, (2.28) has a unique positive root ω0 satisfying
(2.29)ω0=-(An2-D2-2Bn)+(An2-D2-2Bn)2-4(Bn2-Cn2)2
and from (2.25) we obtain
(2.30)sinω0τ=Dω03+AnCnω0-BnDω0D2ω02+Cn2≜G(ω0),cosω0τ=Cnω02-BnCn-AnDω02D2ω02+Cn2≜E(ω0),
then (2.16) has one imaginary root iω0 when
(2.31)τj0=τ00+2jπω0,(j=0,1,2,…),
where τ00 satisfies
(2.32)τ00={arc cos(E(ω0))ωn,ifG(ω0)>0,arc cos(2π-E(ω0))ω0,ifG(ω0)<0,Letλ(τ)=α(τ)+iω(τ) be the root of (2.14) satisfying α(τj0)=0 and ω(τj0)=ω0 when τ is close to τj0. Now by some simple calculations we obtain
(2.33)[dαdτ|τ=τj0]-1=Re[(2λ+An+De-λτ-(Cn+Dλ)τe-λτ(Cn+Dλ)λe-λτ)|τ=τj0]=2Cncosω0τ+Dω0sinω0τCn2+D2ω02-D2ωn2Cn2ω02+D2ω04+Anω0(-Dωncosω0τ+Cnsinω0τ)Cn2ω02+D2ω04,
since ω02=Dωnsinω0τ+Cncosω0τ and Anω0=-Dω0cosω0τ+Cnsinω0τ, and from the expression of ω02 in (2.29), we immediately see that
(2.34)[dαdτ|τ=τj0]-1=2ω02Cn2+D2ω02+An2-D2Cn2+D2ω02=(D2-An2)2+4Cn2Cn2+D2ω02.
Therefore, sign[dα/dτ|τ=τj0]-1=1, that is, dα/dτ|τ=τj0>0. This implies that all the roots that cross the imaginary axis at iω cross from left to right as τ increases.Hence the transversality condition holds and accordingly Hopf bifurcation occurs atτ=τ00, and τ=τj0(j=0,1,2,…) are Hopf bifurcation values of system (1.7) and these Hopf bifurcations are all spatially homogeneous. This completes the proof.


## 3. Direction and Stability of Hopf Bifurcation

In the previous section, we have already obtained that system (1.7) undergoes Hopf bifurcation at the positive equilibrium E*=(u*,v*) when τ crosses through the critical value τj0(j=0,1,2,…). In this section, we will study the direction of the Hopf bifurcation and the stability of the bifurcating periodic solutions by employing the normal form method [19, 20] as well as center manifold theorem [21] for partial differential equations with delay. Then we compute the direction and stability of the Hopf bifurcation when τ0=τj0 for fixed j∈{0,1,2,…}.Without loss of generality, we denote the critical value ofτ by τ0 and set τ=τ0+μ, then μ=0 is the Hopf bifurcation value of system (2.6). Rescaling the time by t→t/τ to normalize the delay, system (2.6) can be written in the following form:
(3.1)∂U(t)∂t=τ0dΔU(t)+τ0L(Ut)+G(Ut,μ),
where
(3.2)L(ϕ)=(-βu*ϕ1(0)-cu*(mv*+1)2ϕ1(0)bv*(mv*+1)ϕ2(-1)-bmu*v*(mv*+1)2ϕ2(-1)),G(ϕ,μ)=μdΔϕ(0)+μL(ϕ)+(τ0+μ)F(ϕ,μ),F(ϕ,μ)=(-βϕ12(0)-c(mv*+1)2ϕ1(0)ϕ2(0)+cmu*(cm+1)3ϕ22(0)(ϕ2(0)-mv*mv*+1ϕ2(-1))(bmv*+1ϕ1(-1)-bmu*(mv*+1)2ϕ2(-1)))
for ϕ∈C([-τ,0];X).From Section2, we know that ±iω0τ0 are a pair of simple purely imaginary eigenvalues of the liner system
(3.3)∂U(t)∂t=τ0dΔU(t)+τ0L(Ut),
and the following liner functional differential equation:
(3.4)U˙(t)=τ0L(Ut).By the Riesz representation theorem, there exists a2×2 matrix function η(θ,τ0) of bounded variation for θ∈[-1,0], such that
(3.5)(τ0+μ)L(τ0)(ϕ)=∫-10‍dη(θ,τ0)ϕ(θ)forϕ∈C([-1,0],ℝ2).
In fact, we can choose
(3.6)η(θ,τ0)=(τ0+μ)H1δ(θ)+(τ0+μ)H2δ(θ+1),
where
(3.7)H1=(-βu*-cu*(mv*+1)200),H2=(00bv*mv*+1-bmu*v*(mv*+1)2),
and δ is the Dirac delta function.Forϕ∈C1([-1,0],ℝ2), we define A(0) as
(3.8)A(0)ϕ(θ)={ddθϕ(θ),θ∈[-1,0),∫-10ddθdη(θ,τ0)ϕ(θ),θ=0,R(0)ϕ(θ)={0,θ∈[-1,0),F(ϕ,μ),θ=0.
For ψ=(ψ1,ψ2)∈C1([-1,0],(ℝ2)*), we define
(3.9)A*(ψ(s))={-ddsψ(s),s∈(0,1],∫-10ddθψ(-ξ)dη(θ,0),s=0.
Then A(0) and A* are adjoint operators under the following bilinear form:
(3.10)〈ψ(s),ϕ(θ)〉=ψ-(0)ϕ(0)-∫-10‍∫ξ=0θ‍ψ-(ξ-θ)dη(0,θ)ϕ(ξ)dξ,
where η(θ)=η(0,θ).We note that±iω0τ0 are the eigenvalues of A(0). Since A(0) and A* are two adjoint operators, ±iω0τ0 are also eigenvalues of A*. We will first try to obtain eigenvector of A(0) and A* corresponding to the eigenvalue iω0τ0 and -iω0τ0, respectively.Letq(θ)=(1,ρ)Teiω0τ0θ,(θ∈[-1,0]) be the eigenvector of A(0) corresponding to the eigenvalue iω0τ0. Then we have A(0)q(θ)=iω0τ0q(θ) by the definition of eigenvector. Therefore, from (3.6), (3.10), and the definition of A(0), we can get
(3.11)τ0(-βu*-iω0-cu*(mv*+1)2bv*mv*+1e-iω0τ0-bmu*v*(mv*+1)2e-iω0τ0-iω0)(1ρ)=(00),
here,
(3.12)ρ=-(mv*+1)2(βu*+iω0)cu*.On the other hand, suppose thatq*(S)=D(1,r)eiω0τ0S is the eigenvector of A* corresponding to the eigenvalue -iω0τ0. By the definition of A*, we have
(3.13)τ0(-βu*+iω0bv*mv*+1e-iω0τ0-cu*(mv*+1)2-bmu*v*(mv*+1)2e-iω0τ0+iω0)(DDr)=(00),
where
(3.14)r=-mv*mv*+1-mβu*v*(mv*+1)+cu*ω0(mv*+1)2i,
and we also assume that 〈q*(S),q(θ)〉=1. To obtain the value of D, from (3.10) we have
(3.15)〈q*(S),q(θ)〉=D-{(1,r-)(1,ρ)T-∫-10‍∫ξ=0θ‍(1,r-)e-i(ξ-θ)ω0τ0dη(θ)(1,ρ)Teiξω0τ0dξ}=D-{1+ρr--∫-10‍(1,r-)θeiθω0τ0dη(θ)(1,ρ)T}=D-{1+ρr-+br-τ0v*mv*+1(mρu*mv*+1-1)e-iω0τ0}.
Thus, we can choose
(3.16)D-=(1+ρr-+br-τ0v*mv*+1(mρu*mv*+1-1))-1eiω0τ0,
such that 〈q*(S),q(θ)〉=0 and 〈q*(S),q(θ)〉=1, that is to say that let ϕ=(q(θ),q-(θ)), ψ=(q*(S),q-*(S))T, then (ψ,ϕ)=I, where I is the unit matrix. Then the center subspace of system (3.4) is P=span{q(θ),q-(θ)}, and the adjoint subspace P*=span{q*(S),q-*(S)}. Let f0=(f01,f02), where
(3.17)f01=(10),f02=(01),
by using the notation from [20], we also define
(3.18)cf0=c1f01+c2f02,
for c=(c1c2)∈𝒞=C([-1,0],X).And the center subspace of linear system (3.4) is given by PCN𝒞, where
(3.19)PCNϕ=ϕ(ψ,〈ϕ,f0〉)0·f0,ϕ∈𝒞,PCN𝒞=(q(θ)z+q-(θ)z-)·f0,z∈𝒞,
and 𝒞=PCN𝒞⊕PS𝒞, where PS𝒞 is the stable subspace.Following Wu [20], we know that the infinitesimal generator AU of linear system (3.4) satisfies
(3.20)U˙(t)=AUU
moreover, U∈dom(AU) if and only if
(3.21)U˙(t)∈𝒞,U(0)∈dom(Δ),U˙(0)=τ0ΔU(0)+τ0L(U(0)).As the formulas to be developed for the bifurcation direction and stability are all relative toμ=0 only, we set μ=0 in system (3.4) and can obtain the center manifold
(3.22)W(z,z-)=W20(θ)z22+W11(θ)zz-+W02(θ)z-22+⋯
with the range in PS𝒞. The flow of system (3.4) in the center manifold can be written as follows:
(3.23)Ut=ϕ(z(t),z-(t))T·f0+W(z(t),z-(t)),
where
(3.24)z˙(t)=iω0τ0z(t)+q*(0)〈G(ϕ(z(t),z-(t)))T·f0+W(z,z-,0),f0〉.We rewrite (3.24) as
(3.25)z˙(t)=iω0τ0z(t)+g(z(t),z-(t)),
with
(3.26)g(z(t),z-(t))=q*(0)〈G(ϕ(z(t),z-(t))T·f0+ω(z,z-,0),f0〉=g20z22+g11zz-+g02z-22+g21z-2z-2.DenoteG(ϕ,0)=τ0(G1,G2)T, and let
(3.27)f1(u,v)=(u+u*)(α-β(u+u*))-c(u+u*)vm(v+v*)+1,f2(u,v)=(v+v*)(-r+b(u(t-τ)+u*)m(v(t-τ)+v*)+1).
Afterwards, from Taylor formula, we have
(3.28)G1=-βϕ12(0)-c(mv*+1)2ϕ1(0)ϕ2(0)+cmu*(mv*+1)3ϕ22(0)+cm(mv*+1)3ϕ1(0)ϕ22(0)-cm2u*(mv*+1)4ϕ23(0),G2=bmv*+1(ϕ2(0)-mv*mv*+1ϕ2(-1))(ϕ1(-1)-mu*mv*+1ϕ2(-1))(1-mmv*+1ϕ2(-1)).From (3.26) and (3.28), we have
(3.29)g20=2D-τ0[-β-cρ(mv*+1)2+cmρ2u*(mv*+1)3]+2D-r-τ0[bρmv*+1e-iω0τ0+bmρv*(mv*+1)2(mρu*mv*+1-1)e-2iω0τ0],g11=D-τ0[-2β+(2mρ-ρu*(mv*+1)3-ρ-+ρ(mv*+1)2)(c+bmv*)]+D-τ0[br-mv*+1(ρ-e-iω0τ0+ρeiω0τ0)],g02=2D-τ0[-β-cρ-(mv*+1)2+cmρ-2u*(mv*+1)3]+2D-r-τ0[bm2ρ-2u*v*(mv*+1)3e2iω0τ0-mbρ-v*(mv*+1)2e2iω0τ0+bρ-mv*+1eiω0τ0],g21=2D-τ0[-β(2W111(0)+W202(0))]+2D-τ0[-c(mv*+1)2(ρW111(0)+12ρ-W201(0)+12W202(0)+W112(0))]+2D-τ0[cmu*(mv*+1)3(2ρW112(0)+ρ-W202(0))]+2D-τ0[br-mv*+1(ρW111(-1)+12ρ-W201(-1))]+2D-τ0[br-mv*+1(12W202(0)eiω0τ0+W112(0)e-iω0τ0)]+2D-τ0[-mbv*r-(mv*+1)2e-iω0τ0(W112(-1)+ρW111(-1))]+2D-τ0[-mbv*r-2(mv*+1)2eiω0τ0(W202(-1)+ρ-W201(-1))]+2D-τ0[m2bu*v*r-(mv*+1)3(2W112(-1)ρe-iω0τ0+W202(-1)ρ-eiω0τ0)].SinceW11(θ) and W20(θ) for θ∈[-1,0] appear in g21, we need to compute them. It follows from (3.26) that
(3.30)W˙(z,z-)=W20zz˙+W11z˙z-+W11zz-˙+W20z-z-˙+,…,AUW=AUW20z22+AUW11zz-+AUW20z-22+⋯.
In addition, W(z,z-) satisfies
(3.31)W˙=AUW+H(z,z-),
where
(3.32)H(z,z-)=H20z22+H11zz-+H02z-22+⋯=X0G(Ut,0)-ϕ(ψ,〈X0G(Ut,0),f0〉)·f0.
Thus, from (3.24) and (2.18), we can get
(3.33)[2iω0τ0-AU]W20=H20,-AUW11=H11,[-2iω0τ0-AU]W02=H02.
Note that AU has only two eigenvalues ±iω0τ0, therefore, (3.33) has unique solution Wij in Q given by
(3.34)W20=[2iω0τ0-AU]-1H20,W11=-AU-1H11,W02=[-2iω0τ0-AU]-1H02.
From (2.18), we know that for θ∈[-1,0)(3.35)H(z,z-)=-ϕ(θ)ψ(0)〈G(Ut,f0)〉·f0=-[q(θ)q*(0)+q-(θ)q-*(0)]〈G(Ut,f0)〉·f0=-[q(θ)g(z,z-)+q-(θ)g-(z,z-)]·f0.
Therefore, for θ∈[-1,0)(3.36)H20(θ)=-[g20q(θ)+g-02q-(θ)]·f0,(3.37)H11(θ)=-[g11q(θ)+g-11q-(θ)]·f0.
From (3.34), (3.36), and the definition of AU, it follows that
(3.38)W˙20(θ)=2iω0τ0W20(θ)+[g20q(θ)+g-02q-(θ)]·f0.
Noting that q(θ)=q(0)eiω0τ0θ,θ∈[-1,0), we have
(3.39)W20(θ)=ig20ω0τ0q(0)eiω0τ0θ·f0+ig-203ω0τ0q-(0)e-iω0τ0θ·f0+E1′e2iω0τ0θ.
Similarly, from (3.34) and (3.37), we obtain
(3.40)W11(θ)=-ig11ω0τ0q(0)eiω0τ0θ·f0+ig-11ω0τ0q-(0)e-iω0τ0θ·f0+E2′.In what follows, we will seek appropriateE1′ and E2′ in (3.39) and (3.40). From the definition of AU and (3.33), we have
(3.41)∫-10‍η(θ)W20(θ)=2iω0τ0W20(θ)-H20(θ),(3.42)∫-10‍η(θ)W11(θ)=-H11(0),
where η(θ)=η(0,θ), then
(3.43)H20=-[g20q(0)+g-02q-(0)]·f0+2τ0(-β-cρ(mv*+1)2+cmρ2u*(mv*+1)3bρmv*+1e-iω0τ0+bmρv*(mv*+1)2(mρu*mv*+1-1)e-2iω0τ0).
Substituting (3.43) into (3.41), note that
(3.44)(iω0τ0-∫-10‍eiω0τ0dη(θ))q-(0)=0,(-iω0τ0-∫-10‍eiω0τ0dη(θ))q(0)=0,
then we deduce
(3.45)2τ0(-β-cρ(mv*+1)2+cmρ2u*(mv*+1)3bρmv*+1e-iω0τ0+bmρv*(mv*+1)2(mρu*mv*+1-1)e-2iω0τ0)=g20q(θ)+g-02q-(θ)+2iω0τ0W20(θ)-∫-10‍dη(θ)W20(θ)=-g20q(0)+g-02q-(0)3+2iω0τ0E-∫-10‍dη(θ)[ig20ω0τ0q(0)eiω0τ0θ+ig-203ω0τ0q-(0)eiω0τ0θ+Ee2iω0τ0]=(2iω0τ0I-∫-10‍dη(θ)e2iω0τ0θ)E=[2iω0τ0(1001)-τ0(-βu*-cu*(mv*+1)2bv*(mv*+1)e-2iω0τ0-bmu*v*(mv*+1)2e-2iω0τ0)]E=τ0(2iω0+βu*cu*(mv*+1)2-bv*(mv*+1)e-2iω0τ02iω0+bmu*v*(mv*+1)2e-2iω0τ0)E.
Finally, we arrive at
(3.46)E=2E1′(-β-cρ(mv*+1)2+cmρ2u*(mv*+1)3bρmv*+1e-iω0τ0+bmρv*(mv*+1)2(mρu*mv*+1-1)e-2iω0τ0),
where
(3.47)E1′=(2iω0+βu*cu*(mv*+1)2-bv*(mv*+1)e-2iω0τ02iω0+bmu*v*(mv*+1)2e-2iω0τ0)-1,
and from the above equation we can find the value of E1 and E2.Following the similar steps, we also get(3.48)H11=-[g20q(0)+g-02q-(0)]·f0+τ0(-2β+(2mρ-ρu*(mv*+1)3-ρ-+ρ(mv*+1)2)(c+bmv*)br-mv*+1(ρ-e-iω0τ0+ρeiω0τ0)),
then
(3.49)τ0(-2β+(2mρ-ρu*(mv*+1)3-ρ-+ρ(mv*+1)2)(c+bmv*)br-mv*+1(ρ-e-iω0τ0+ρeiω0τ0))=g11q(0)+g-11q-(0)-∫-10‍dη(θ)W11(θ)=g11q(0)+g-11q-(0)-∫-10‍dη(θ)[ig11ω0τ0q(0)eiω0τ0θ+ig-11ω0τ0q-(0)eiω0τ0θ+E′]=∫-10‍dη(θ)E′=τ0(βu*cu*(mv*+1)2-bv*mv*+1bmu*v*(mv*+1)2)E′.
Finally, we arrive at
(3.50)E′=E2′(-2β+(2mρ-ρu*(mv*+1)3-ρ-+ρ(mv*+1)2)(c+bmv*)br-mv*+1(ρ-e-iω0τ0+ρeiω0τ0)),
where
(3.51)E2′=(βu*cu*(mv*+1)2-bv*mv*+1bmu*v*(mv*+1)2)-1.In a similar manner, we can calculateE'1 and E'2. Then g21 can be expressed. Based on the above analysis, it is enough to see that each gij is determined by the parameters. Thus, we can compute the following values which determine the direction and stability of the following bifurcating periodic orbits:
(3.52)C1(0)=i2ω0τ0(g11g22-2|g11|2-|g02|23)+g212,μ2=-Re(C1(0))Re(λ′(τ0)),β2=2Re(C1(0)),T2=-Im(C1(0))+μ2Im(λ′(τ0))ω0τ0.Theorem 3.1.
For system (1.7),(i)
μ2 determines the directions of the Hopf bifurcation: if μ2>0(<0), the direction of the Hopf bifurcation is forward (backward), that is, the bifurcating periodic solutions exist for τ>τ0=τj0(τ<τ0=τj0);(ii)
β2 determines the stability of the bifurcating periodic solutions on the center manifold: if β2<0(>0), the bifurcating periodic solutions are orbitally asymptotically stable (unstable);(iii)
T2 determines the period of the bifurcating periodic solutions: the period increases (decreases) if T2>0(<0).


## 4. Conclusions and Remarks

In this paper, under homogeneous Neumann boundary conditions, we have analyzed dynamical behaviors of the diffusion predator-prey system (1.7) with and without delay. The value of this study lies in two aspects. First, it presents the stability of positive equilibrium E*=(u*,v*) of system (1.7) with and without delay, and the existence of Hopf bifurcation, which indicates that the dynamical behaviors become rich and complex with delay. Second, it shows the analysis of stability of Hopf bifurcation, from which one can find that small sufficiently delays cannot change the stability of the positive equilibrium and large delays cannot only destabilize the positive equilibrium but also induce an oscillation near the positive equilibrium.Next, numerical simulations are performed to illustrate results with respect to the theoretical facts under the special example. Whenτ and di(i=1,2) are all equal 0, we demonstrate that the positive equilibrium E* is globally asymptotically stable (see, Figure 1), which means that if the intraspecific competitions of the prey and the predator dominate the inter-specific interaction between the prey and the predator, then both the prey and the predator populations are permanent [14].Figure 1
The phase portrait of model (1.7) without the diffusion and delay effects, that is, τ=0, d1=d2=0. The other parameters are taken as α=2, β=0.2, γ=0.5, b=2, andc=3, m=1. In this case, E1=(0,0) and E2=(10,0) are saddle points; the unique positive equilibrium E*=(0.6623,1.6491) is globally asymptotically stable. The dashed curve is the u-nullcline u(α-βu)-cuv/(mv+1)=0, and the dotted curve is the v-nullcline v(-r+bu/(mv+1))=0.In addition, we consider the dynamics of system (1.7) affecting by both spatial diffusion and time-delay factors with fixed parameters d1=0.01 and d2=0.1. In this case, τ00=1.0289, Re(C1(0))=-0.6689<0. If τ<τ00, the positive equilibrium E* is remain stability (see Figure 2), which indicates that the predators and preys can coexist in stable conditions.The solution of system (1.7) tends to the positive equilibrium E*. The parameters are taken as α=2, β=0.2, γ=0.5, b=2, c=3, m=1, d1=0.01, d2=0.1, τ=0.9, and the initial values u(x,t)=0.662+0.01tcosx, v(x,t)=1.649+0.01tsinx, t∈[-0.9,0], x∈[0,π]. In this case, τ<τ00=1.0289.
(a)
(b)While ifτ>τ00, the positive equilibrium E* losses its stability and Hopf bifurcation occurs, which means that a family of stable periodic solutions bifurcate from E* since μ2=0.5515>0, β2=-0.1338<0 (see Figure 3). The numerical result indicates that the predator coexists with the prey with oscillatory behaviors.The solution of system (1.7) tends to aperiodic orbit. The parameters are taken as α=2, β=0.2, γ=0.5, b=2, c=3, m=1, d1=0.01, d2=0.1, τ=1.3, and the initial values u(x,t)=0.662+0.01tcosx, v(x,t)=1.649+0.01tsinx, and t∈[-1.3,0], x∈[0,π]. In this case, τ>τ00=1.0289.
(a)
(b)In the present paper, we incorporate time delay into biological system due to the gestation of the predator, which causes stable equilibrium to become unstable and causes the populations to oscillate via Hopf bifurcation. That is to say that the effect of delay for the population dynamics is tremendous. From a biological perspective, the time delay of species may be related to the gestation of the predator, mature stage from juvenile to adult, the interaction time between prey and predator and others. In these cases, the methods and results in the present paper may provide a favorable value on controlling ecological population. It would be more accurate to describe the growth rate of population.



---

*Source: 323186-2012-10-30.xml*