# Dynamics and Patterns of a Diffusive Prey-Predator System with a Group Defense for Prey

Authors: Honglan Zhu; Xuebing Zhang
Journal: Discrete Dynamics in Nature and Society
 (2018)

---

## Abstract

We study a diffusive prey-predator system with a group defense for prey. Under Neumann boundary condition, we analyze local and stability of nonnegative constant steady states and the existence and nonexistence of nonconstant steady states. These results also exhibit the critical role of the system parameters leading to the formation of spatiotemporal patterns.


---

## Body

## 1. Introduction

The predator-prey system first proposed by [1, 2] is one of the fundamental ecological systems in both ecology and mathematical ecology. Based on different settings, various types of predator-prey models described by differential systems have been proposed and the dynamics of these systems are studied [3–6]. The basic form of these models is as follows:(1)dxdt=rx1-xK-Px,y,dydt=-sy+cPx,y,where r is the intrinsic growth rate and K is the environmental carrying capacity of prey population, and the function P(x) is the functional response; the constant c(>0) is the ratio of biomass conversion and s is the natural death rate of predator species. The simplest functional response is Lotka-Volterra function which is described as (2)Px,y=axy,0≤x≤ka,ky,x≥ka,which is also called Holling type I function. However, the curve defined by the Lotka-Volterra response function is a straight line through the origin and is unbounded. Thus, more reasonable response functions should be nonlinear and bounded. In 1913, Michaelis and Menten proposed the response function(3)Px,y=mxya+x,where m>0 denotes the maximal growth rate of the species and a>0 is the half-saturation constant. It is now referred to as a Michaelis-Menten function or a Holling type II function. Another class of response function is (4)Px,y=mx2ya+bx+x2which is called a sigmoidal response function, while the simplification(5)Px,y=mx2ya+x2is known as a Holling type III function. Some authors [7, 8] considered system (1) with following response function:(6)Px,y=mxya+x2which is called Holling type IV function. Besides, Beddington-DeAngelis type P(x,y)=βxy/α+x+my and more complicated functional response P(x,y)=βx2y/x2+my2 are also considered by some researchers [9, 10].Recently, some works consider the case when animals join together in herds in order to provide a self-defense from predators. In [11], the authors argued that it is more appropriate to model the response functions of prey that exhibit herd behavior in terms of the square root of the prey population. Inspired by this thought, the authors in [12] choose response function P(x)=x to reflect this fact. When motion is allowed, [13] considered the spatiotemporal behavior of a prey-predator system with a group defense for prey by means of extensive computer simulations. The proposed model is as follows:(7)∂u∂t=d1Δu+r1u1-uK-ρuαv,x,t∈Ω×0,+∞,∂v∂t=d2Δv+βuαv-r2v,x,t∈Ω×0,+∞ux,0=u0x,vx,0=v0x,x∈Ω,where u and v denote, respectively, the densities of prey and predator species. r is the growth rate of prey species, K is its carrying capacity, r2 is the mortality rate of predator species, ρ is the search efficiency of predator for prey, β is the biomass conversion coefficient, and α∈(0,1) represents a kind of aggregation efficiency. The local dynamics for nonspatial model was studied, such as Hopf bifurcation and existence of extinction domain. For model (7), the authors only give some numerical simulations to find some spatiotemporal features. Reference [14] considers the direction and the stability of the bifurcating periodic solutions for model (7) with α=1/2 under Neumann boundary conditions. Reference [15] investigated the global dynamics of nonspatial model including the nonexistence of periodic orbits and the existence and uniqueness of limit cycles. We refer readers to [16–21] as some other related works on predator-prey model with herd behavior.It is noted that up to now no one has studied the existence and nonexistence of positive steady state solutions of (7). Therefore, the main aim of this article is to study the existence and nonexistence of nonconstant positive solutions of the following elliptic system:(8)-d1Δu=r1u1-uK-ρuαv,inΩ,-d2Δv=βuαv-r2v,inΩ,∂u∂ν=∂v∂ν=0,on∂Ω,where ν is the outward unit normal vector on ∂Ω, and we impose a homogeneous Neumann type boundary condition, which implies that (8) is a closed system and has no flux across the boundary ∂Ω.The structure of this paper is arranged as follows. In Section2, we estimate the a priori bounds of positive solutions of (7). In Section 3, the local and global stabilities of nonnegative constant steady states of (7) are discussed. In Section 4, we give a priori estimate for the positive solutions of (8) by using maximum principle and Harnack inequality. In Section 5, we give a nonexistence result of nonconstant solutions of (8). In Section 6, we consider the existence of nonconstant positive solutions of (8). Finally, to support our theoretical predictions, some numerical simulations are given.


## 2. Basic Dynamics and a Priori Bound

Theorem 1.
For system (7), one has the following.(a)
Ifu0(x)≥0, v0(x)≥0, then system (7) has a unique solution (u(t,x),v(t,x)) such that u(t,x)>0, v(t,x)>0 for t∈(0,+∞) and x∈Ω¯.(b)
Any solution(u,v) of (7) satisfies (9)limsupt→+∞⁡ut,x≤K,∫Ωvx,tdx≤βρ1+r1r2KΩ.Proof.
(a) Define(10)fu,v=r1u1-uK-ρuαv,gu,v=βuαv-r2v. Then fv=-ρuα≤0 and gu=βαuα-1v≥0 in R+2¯={u≥0,v≥0}. Hence, (7) is a mixed qusi-monotone system. Consider following system:(11)dudt=r1u1-uK-ρuαv,dvdt=βuαv-r2v,u0=u0,v0=v0.Assume ut;u0,v0, vt;u0,v0 are the unique solution to system (11). Let (12)maxΩ¯⁡u0x=ϕM,maxΩ¯⁡v0x=ψM. Obviously, u_t,x,v_t,x=(0,0) and (u¯(t),v¯(t))=(u(t;ϕM,ψM),v(t;ϕM,ψM)) are a pair of lower-solution and upper-solution to system (7). Therefore, according to the Theorem  8.3.3 in [22] or Theorem5.3.2 in [23], system (7) has a unique globally defined solution (u(x,t),v(x,t)) which satisfies (13)0≤ux,t≤ut;ϕM,ψM,0≤vx,t≤vt;ϕM,ψM. The strong maximum principle implies that u(x,t),v(x,t)>0 when t>0 for all x∈Ω¯.
(b) By the first equation of (7), we easily obtain the fact that r1u1-u/K-ρuαv≤r1u1-u/K in 0,+∞×Ω; the first result follows easily from the simple comparison argument for parabolic problems, and thus there exists T∈(0,+∞) such that u(t,x)≤K+ϵ in [T,+∞)×Ω for an arbitrary constant ϵ>0.
For the estimate ofv(x,t), let U(t)=∫Ωu(x,t)dx, V(t)=∫Ωv(x,t)dx; then(14)dUdt=∫Ωutdx=d1∫ΩΔudx+∫Ωr1uu-uK-ρuαvdx=∫Ωr1uu-uK-ρuαvdx,(15)dVdt=∫Ωvtdx=d2∫ΩΔvdx+∫Ωβuαv-r2vdx=∫Ωβuαv-r2vdx.Multiplying (14) by β/ρ and adding it to (15), we have (16)βρU+Vt=-r2V+βr1ρ∫Ωuu-uKdx≤-r2βρU+V+βr2ρ+βr1ρρK+ϵΩ. Integration of the inequality leads to(17)∫Ωvx,tdx=Vt<βρUt+Vt≤βρ1+r1r2K+ϵΩ.


## 3. Stability of the Nonnegative Constant Steady States of (7)

In this section, we will analyze the stability of nonnegative constant steady states of (7). By the direct computation, we see that the possible nonnegative constant steady states of (7) are (18)E0=0,0,E1=K,0,E∗=u∗,v∗, where u∗=r2/β1/α,v∗=r1u∗(1-α)(1-u∗/K)/ρ. Obviously, the positive constant steady state E∗ exists if r2<βKα holds.Notation 1.
Let0=μ0<μ1<μ2<⋯<μn<⋯→∞ be the eigenvalues of -Δ on Ω under homogeneous Neumann boundary condition. We define the following space decomposition:(i)
S(μn) is the space of eigenfunctions corresponding to μi for n=0,1,2,….(ii)
Xij≔{c·ϕij:c∈R2}, where {ϕij} are orthonormal basis of S(μn) for j=1,2,…,dim[S(μn)].(iii)
X≔u=u,v∈C1Ω¯2:∂u/∂n=∂v/∂n=0, and so X=⊕i=1∞Xi, where Xi=⊕j=1dim[S(μj)]Xij.LetE~ be a nonnegative constant steady state of (7); then the linearization of (7) at a constant solution E~ can be expressed by(19)ut=DΔ+Ju,where D=diag(d1,d2), u=ux,t,vx,tT, and (20)J=r1-2r1u∗K-ραu∗α-1v∗-ρu∗α-βv∗u∗α-1βu∗α-r2. In view of Notation 1, we can induce the eigenvalues of system (19) confined on the subspace Xi. If λ is an eigenvalue of (19) on Xi, it must be an eigenvalue of the matrix -μnD+J for each n≥0. It is easy to see that λ satisfies the characteristic equation(21)λ2+d1μn+d2μn-r1+2r1u∗K+ραu∗α-1v∗λ+d2μnd1μn-r1+2r1u∗K+ραu∗α-1v∗+ρβu∗2α-1v∗=0.Theorem 2.
(i) The trivial equilibriumE0=(0,0) is unstable.
(ii) IfβKα<r2, then E1=(K,0) is globally asymptotically stable.
(iii) If1-α-2-α/Ku∗<0, then E∗ is locally asymptotically stable.Proof.
(i) ForE0=(0,0), the eigenvalues are (22)λ1n=r1-d1μn,λ2n=-d2μn. Obviously, E0 is unstable.
(ii) ForE1=(K,0), the eigenvalues are (23)λ1n=-r1-d1μn,λ2n=βKα-r2-d2μn. If βKα<r2, then λ1n and λ2n are all negative. Therefore E1 is locally asymptotically stable. Indeed, E1 is globally asymptotically stable.
On account of Theorem1, we have limsupt→+∞⁡maxΩ¯⁡u(·,t)≤K, and thus there exists T1∈(0,+∞) such that, for an arbitrary constant 0<ϵ<r2/β1/α-K, (24)u·,t≤K+ϵ,t≥T1. It follows from the second equation of (7) that (25)vt-d2Δv≤vβK+ϵ⁡α-r2,t≥T1. Therefore, limsupt→+∞⁡maxΩ¯⁡v(·,t)≤0, and there exists T2>T1 such that (26)v·,t≤ϵ,t≥T2. It follows from the first equation of (7) that (27)ut-d1Δu≥r1u1-uK-ρuαϵ,t>T2,x∈Ω. On account of and the arbitrariness of ϵ>0, we have liminft→+∞⁡minΩ¯⁡u(·,t)≥K. This combined with limsupt→+∞⁡maxΩ¯⁡u(·,t)≤K allows us to derive(28)limt→+∞⁡maxΩ¯⁡u·,t-K=0. Hence, E1 is globally asymptotically stable when r2>Kαβ.
(iii) WhenE∗=(u∗,v∗) exists, the corresponding characteristic equation is as follows:(29)λ2+d1μn+d2μn-r11-α-2-αKu∗λ+d2μnd1μn-r11-α-2-αKu∗+ρβu∗2α-1v∗=0.Obviously, we have (30)λ1n+λ2n=r11-α-2-αKu∗-d1μn-d2μn,λ1nλ2n=d2μnd1μn-r11-α-2-αKu∗+ρβu∗2α-1v∗. If 1-α-2-α/Ku∗<0, then λ1n+λ2n<0 and λ1nλ2n>0. Hence, all the roots of (29) have negative real part which means that E∗ is locally asymptotically stable when 1-α-2-α/Ku∗<0.


## 4. The Prior Estimate

In this section, we will give some a priori estimates of positive solutions to (8). Firstly, we give two known lemmas.Lemma 3 (Harnack inequality (cf. [24])).
Letω∈C2(Ω)∩C1(Ω¯) be a positive classical solution to (31)Δωx+cxωx=0inΩ,∂ω∂n=0on∂Ω. Then there exists a positive constant C such that (32)maxΩ¯⁡ω≥CminΩ¯⁡ω.Lemma 4 (maximum principle (cf. [25])).
Suppose thatg∈C(Ω¯×R).
(i) Assume thatω∈C2(Ω)∩C1(Ω¯) satisfies (33)Δωx+cxωx≥0inΩ,∂ω∂n≤0on∂Ω. If ω(x0)=maxΩ¯⁡ω, then g(x0,ω(x0))≥0.
(ii) Assume thatω∈C2(Ω)∩C1(Ω¯) satisfies (34)Δωx+cxωx≤0inΩ,∂ω∂n≥0on∂Ω. If ω(x0)=minΩ¯⁡ω, then g(x0,ω(x0))≤0.Lemma 5.
For any positive solution(u,v) of system (8), (35)0<ux<K,0<vx<βρd1d2+r1r2Kfor  anyx∈Ω¯.Proof.
Form Lemma4, u(x)≤K and from the strong maximum principle u(x)<K for all x∈Ω¯. Multiplying the first equation of (8) by β/ρ and adding it to the second equation, we have (36)-βρd1Δu+d2Δv=βρr1u1-uK-r2v≤βρr1K+d1r2βρd2K-r2d2βρd1u+d2v. Then the maximum principle implies that (37)βρd1u+d2v≤βρd2r1r2+d1K. Hence, v(x)<β/ρd1/d2+r1/r2K.In the following, we estimate the positive lower bound of positive solution of (8).Theorem 6.
LetΩ be a bounded smooth domain in Rn. There exist two positive constants C_<C¯ depending possibly on d1, d2, K, β, α, ρ, and Ω, such that such that any positive solution (u(x),v(x)) of system (8) satisfies (38)C_≤ux≤C¯,C_≤vx≤C¯for  anyx∈Ω¯.Proof.
From Lemma5, we obtain (39)ux,vx≤C¯≔max⁡K,βρd1d2+r1r2K, where C¯ depends on d1, d2, K, β, α, and ρ.
From Lemma3, we obtain the fact that there exists a positive constant C2 such that(40)supΩ¯⁡ux≤C2infΩ¯⁡ux,supΩ¯⁡vx≤C2infΩ¯⁡vx.
On the contrary, suppose the result is false. Then there exists a sequence{(un,vn)} of positive solutions to system (8) such that(41)supΩ¯⁡un⟶0orsupΩ¯⁡vn⟶0asn⟶+∞.By the regularity theory for elliptic equations, there exists a subsequence of {(un,vn)}, which will be denoted again by {(un,vn)}, such that {(un,vn)}→(u0,v0) in C2(Ω¯) as n→+∞. Observe that u0≤K and, from (41), either u0≡0 or v0≡0. Therefore, we have the following two cases:
(i)u0≡0, v0≢0; or u0≡0, v0≡0.
(ii)u0≢0, v0≡0.
Since{(un,vn)} is a positive solution of (8), one can obtain the following integral equation by integrating (8) for un and vn over Ω: (42)∫Ωr1un1-unK-ρunαvndx=0,∫Ωβunαvn-r2vndx=0.
(i) In this case,u0≡0; then (43)βunα-r2⟶-r2<0 uniformly as n→∞ and vn>0; then for sufficiently large n, we have (44)∫Ωβunαvn-r2vndx<0, which is a contradiction.
(ii) Ifu0≢0, v0≡0, then this implies that u0 satisfies (8). So u0≡K for large n. Thus (45)βunα-r2⟶βKα-r2>0 for large n since βKα<r2, which derives a contradiction again to the second integral equation of (42). This completes the proof.


## 5. Nonexistence of Nonconstant Positive Steady States

In this section, we can show the nonexistence of nonconstant positive solutions to system (8) when the diffusion coefficients d1 and d2 are large.Theorem 7.
There exists a positive constantd∗ such that elliptic problem (8) has no nonconstant positive solution if min⁡{d1,d2}>d∗.Proof.
Suppose that(u(x),v(x)) is a nonconstant positive solution of system (8). Denote u¯=Ω-1∫Ωu(x)dx≥0, v¯=Ω-1∫Ωu(x)dx≥0. Then (46)∫Ωu-u¯dx=0,∫Ωv-v¯dx=0.
Define(47)hu=u¯-uu¯α-uα for u¯≠u. Indeed, we can prove that h(u)>0 and h′(u)>0. In fact, notice (48)h′u=1-αuα+αu¯uα-1-u¯αu¯α-uα2. Let (49)h1u=1-αuα+αu¯uα-1-u¯α, and we have (50)h1′u=α1-αuα-2u-u¯, which implies that h1(u)>min⁡h1(u)=h1(u¯)=0 for u≠u¯. Therefore, we obtain the fact that h(u)>0.
Furthermore, multiplying the first equation of (8) by β/ρ, adding it to the second equation, and integrating over Ω, we get(51)-∫Ωβρd1Δu+d2Δvdx=∫Ωβρr1u1-uK-r2vdx,and then the Neumann boundary conditions lead to (52)r2∫Ωvdx=∫Ωβρr1u1-uKdx≤Kβr14ρΩ. Thus(53)v¯=1Ω∫Ωvdx≤Kβr14r2ρ.
Multiplying the first equation in (8) by u-u¯, we have(54)d1∫Ω∇u-u¯2dx=∫Ωu-u¯r1u1-uK-ρuαvdx=∫Ωu-u¯r1u-u¯-r1Ku2-u¯2-ρuαv+ρu¯αv¯dx≤3r1∫Ωu-u¯2dx+ρKα∫Ωu-u¯v-v¯dx+ρ∫Ωv¯huu-u¯2dx≤3r1+ρKα2+Kβr14r2h0∫Ωu-u¯2dx+ρKα2∫Ωv-v¯2dx.
Multiplying the second equation in (8) by v-v¯, we have(55)d2∫Ω∇v-v¯2dx=∫Ωv-v¯βuαv-r2vdx=∫Ωv-v¯-r2v-v¯+βuαv-βu¯αv¯dx≤r2+βKα+Kβ2r18r2ρh0∫Ωv-v¯2dx+Kβ2r18r2ρh0∫Ωu-u¯2dx.
From (54) and (55) and the Poincaré inequality, we obtain (56)d1∫Ω∇u-u¯2dx+d2∫Ω∇v-v¯2dx≤A∫Ωu-u¯2dx+B∫Ωv-v¯2dx≤1μ1A∫Ω∇u-u¯2dx+B∫Ω∇v-v¯2dx, where (57)A=3r1+ρKα2+Kβr14r2h0+Kβ2r18r2ρh0,B=r2+βKα+Kβ2r18r2ρh0+ρKα2. Hence, if (58)min⁡d1,d2>d∗≔1μ1max⁡A,B, then (59)∇u-u¯=∇v-v¯=0, and (u,v) must be a constant solution.


## 6. Existence of Nonconstant Positive Steady States

In this subsection, we discuss the existence of nonconstant positive solutions to system (8) when the diffusion coefficients d1 and d2 vary while the parameters r1, K, α, ρ, β, and r2 are fixed by using the Leray-Schauder degree theory. Throughout this section, we assume that the positive constant steady state E∗=(u∗,v∗) exists.For simplicity, denoteu=(u,v) and (60)a1=r1-2r1u∗K-ραu∗α-1v∗,a2=ρu∗α,a3=βv∗u∗α-1,Φu=r1u1-uK-ρuαvβuαv-r2v,ΦuE∗=a1-a2a30,D=diag⁡d1,d2,X=u,v∈C1Ω¯2:∂u∂n=∂v∂n=0on∂Ω,X+=u,v:u,v≥0,u,v∈X,Λ=u,v∈X:C-1<u,v<ConΩ¯. Thus, (8) can be written as(61)-DΔu=Φuin∂u∂n=0onΩ,and, obviously, u is a positive solution of (61) if and only if(62)Fu≔u-I-Δ-1D-1Φu+uinX+,where (I-Δ)-1 is the inverse of I-Δ with the homogeneous Neumann boundary condition. As F(·) is a compact perturbation of the identity operator, the Leray-Schauder degree deg(F(·),Λ,0) is well-defined from Theorem 6. By direct computation, we have(63)FuE∗=I-I-Δ-1D-1ΦuE∗+I. If Fu(E∗) is invertible, the index of F is defined as (64)indexF·,E∗=-1γ, where γ is the number of negative eigenvalues of FuE∗. Note that λ is an eigenvalue of Fu(E∗) on Xj if and only if it is an eigenvalue of the matrix (65)Bj=I-11+μjD-1ΦuE∗+I=11+μjμjI-D-1ΦuE∗. Thus Fu(E∗) is invertible if and only if, for all j≥0, the matrix Bj is nonsingular. Writing(66)Hd1,d2;μ=det⁡μI-D-1ΦuE∗=1d1d2det⁡μD-ΦuE∗,we have that if H(d1,d2;μ)≠0, then H(d1,d2;μ)<0 if and only if the number of negative eigenvalues of Fu(E∗) in Xj is odd. The following lemma gives the explicit formula of calculating the index.Lemma 8.
IfH(d1,d2;μi)≠0 for all i≥0, then(67)index⁡F·,E∗=-1γ,γ=∑i≥0,Hd1,d2;μ<0mμi, where m(μi) is the algebraic multiplicity of μi.To facilitate our computation ofdeg(F(·),E∗), we only need consider the sign of det[μD-Φu(E∗)]. The direct calculation gives(68)det⁡μD-ΦuE∗=d1d2μ2-a1d2μ+a2a3.Obviously, nonnegative roots of (68) exist if and only if a12d2-4d1a2a3>0 and a1>0. Assume that μ+ and μ- are the two roots of (68), we have the following conclusion.Theorem 9.
Assuming thatβKα>r2, (69)d2d1>4a2a3a12,a1>0, and there exist i,j∈N, such that 0≤μj<μ-<μj+1≤μi<μ+<μi+1 and ∑k=j+1im(μk) is odd, then (8) has at least one nonconstant positive solution.Proof.
Fort∈[0,1], we define (70)Atu≜-Δ+I-1u+1-td∗+td1fu,vv+1-td∗+td2gu,v, where d∗ is defined in Theorem 7.
The positive solutions of the problem(71)Atu=uinΩ,∂u∂n=0on∂Ωare contained in Λ. Note that u is a positive solution of system (8) if and only if it is a positive solution of (71) with t=1. u∗ is the unique positive constant solution of (71) for any t∈[0,1]. According to the choice of d∗ in Theorem 7, we have E∗ which is the only fixed point of A0.(72)deg⁡I-A0,Λ,0=index⁡I-A0,Λ,E∗=1.
SinceF=I-H(·,1) and if (8) has no other solutions except the constant one E∗, then we have (73)deg⁡I-A1,Λ,0,0=index⁡F,E∗=-1∑k=j+1imμk=-1. On the other hand, by the homotopy invariance of the topological degree,(74)deg⁡I-A0,Λ,0=deg⁡I-A1,Λ,0,which is a contradiction. Therefore, there exists at least one nonconstant solution of (8).


## 7. Numerical Simulation

### 7.1. Global Stability of EquilibriumE1

Consider system (7) with following parameters: d1=0.8,d2=0.9,r1=0.9,β=0.1,r2=0.2,K=2,ρ=0.1, and α=2/3. According to the discussions in Section 3, the steady state E1 is globally asymptotically stable; see Figure 1.Figure 1
The steady stateE1 is globally asymptotically stable.


### 7.2. Stability of Positive Steady StateE∗

Consider system (7) with following parameters: d1=0.8,d2=0.9,r1=0.9,β=0.3,r2=0.2,K=2,ρ=0.1, and α=2/3. According to the discussions in Section 3, the positive steady state E∗ is locally asymptotically stable; see Figure 2.Figure 2
The positive steady stateE∗ is locally asymptotically stable.


## 7.1. Global Stability of EquilibriumE1

Consider system (7) with following parameters: d1=0.8,d2=0.9,r1=0.9,β=0.1,r2=0.2,K=2,ρ=0.1, and α=2/3. According to the discussions in Section 3, the steady state E1 is globally asymptotically stable; see Figure 1.Figure 1
The steady stateE1 is globally asymptotically stable.


## 7.2. Stability of Positive Steady StateE∗

Consider system (7) with following parameters: d1=0.8,d2=0.9,r1=0.9,β=0.3,r2=0.2,K=2,ρ=0.1, and α=2/3. According to the discussions in Section 3, the positive steady state E∗ is locally asymptotically stable; see Figure 2.Figure 2
The positive steady stateE∗ is locally asymptotically stable.


## 8. Conclusions

In this paper, we have investigated the existence/nonexistence of nonconstant positive steady states for a diffusive predator-prey system with a group defense for prey under Neumann boundary conditions. The existence results provide a theoretical support for pattern formation caused by diffusion. We also study the stability of nonnegative equilibria and obtain the fact thatE1 is globally asymptotically stable when βKα<r2. In fact, the positive steady state does not exist at this time. If βKα>r2 and 1-α-2-α/Ku∗<0, then the positive steady state E∗ is locally asymptotically stable. It is easily obtained that when 1-α-2-α/Ku∗=0, characteristic equation (29) has a pair of purely imaginary roots. Therefore, system (7) occurs with Hopf bifurcation, as shown in Figure 3.Figure 3
Hopf bifurcation occurs with parametersr1=0.9,α=2/3,β=0.1,r2=0.18,K=10;ρ=0.2.



---

*Source: 6519696-2018-01-08.xml*