Qualitative analysis of Beddington–DeAngelis type impulsive predator–prey models

Nonlinear Analysis: Real World Applications 11 (2010) 1312–1322

Contents lists available at ScienceDirect

Nonlinear Analysis: Real World Applications journal homepage: www.elsevier.com/locate/nonrwa

Qualitative analysis of Beddington–DeAngelis type impulsive predator–prey models H.K. Baek Department of Mathematics, Kyungpook National University, Daegu 702-701, South Korea

article

abstract

info

Article history: Received 29 January 2008 Accepted 20 February 2009

In the paper, the impulsive predator–prey models with Beddington–DeAngelis functional response are studied. Conditions for the existence and stability of a prey-free solution and for the existence of a nontrivial periodic solution have been established. Also, we find a sufficient condition that the model is permanent and show that the model has complex dynamical behaviors via bifurcation diagrams. © 2009 Elsevier Ltd. All rights reserved.

Keywords: Predator–prey model Beddington–DeAngelis functional response Impulsive differential equation Floquet theory Bifurcation theory

1. Introduction Understanding the dynamical relationships between predator and prey is one of central goals in an ecological system. One of the most important components of the predator–prey relationship are so-called functional responses. Functional response refers to the change in the density of prey attached per unit time per predator as the prey density changes. The Beddington–DeAngelis functional response was introduced by Beddington [1] and DeAngelis et al. [2]. It is similar to the wellknown Holling type functional response but contains an extra term describing mutual interference by predators. Actually, there is much significant evidence to suggest that functional responses with predator interference occur quite frequently in laboratory and natural systems [3]. Thus, we can establish a predator–prey model with Beddington–DeAngelis functional response as the following form [1,4,5]:



x0 (t ) = ax(t ) 1 − y0 (t ) = −Dy(t ) +

x(t )

ex(t )y(t )



K

−

by(t ) + x(t ) + c

mx(t )y(t ) by(t ) + x(t ) + c

,

,

(1.1)

(x(0+ ), y(0+ )) = (x0 , y0 ) = x0 , where x(t ), y(t ) represent the population density of the prey and the predator at time t, respectively. Usually, K is called the carrying capacity of the prey. The constant a is called intrinsic growth rate of the prey. The constants m, D are the conversion rate and the death rate of the predator, respectively. The term by measures the mutual interference between predators. As Cushing [6] pointed out that it is necessary and important to consider models with periodic ecological parameters or perturbations which might be quite naturally exposed (for example, those due to seasonal effects of weather, food supply, mating habits, hunting or harvesting seasons and so on). Such perturbations were often treated continually. However, there are still some other perturbations such as fire, flood, etc, that are not suitable to be considered continually. These impulsive

E-mail address: [email protected]. 1468-1218/$ – see front matter © 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.nonrwa.2009.02.021

H.K. Baek / Nonlinear Analysis: Real World Applications 11 (2010) 1312–1322

1313

perturbations bring sudden change to the system. Let’s think of prey as a pest and predator as a natural enemy of prey. There are many ways to beat agricultural pests. For examples, harvesting on prey, spreading pesticides, releasing natural enemies and so on. Such tactics are discontinuous and periodical. With the idea of periodic forcing and impulsive perturbations, in this paper, we consider the following predator–prey model with periodic constant impulsive immigration of the predator and periodic variation in the intrinsic growth rate of the prey.

   ex(t )y(t ) x(t ) 0   − , x ( t ) = ax ( t ) 1 −   K by(t ) + x(t ) + c   0 mx(t )y(t )  y (t ) = −Dy(t ) + , by(t ) +x(t ) + c  +  x(t ) = (1 − p1 )x(t ),   t = nT ,  +  − p2 )y(t ) + q,  y(t )+ = (1 + (x(0 ), y(0 )) = (x0 , y0 ) = x0 ,

   

t 6= nT ,

  

(1.2)

where all parameters are positive constants, T is the period of the impulsive immigration or stock of the predator, 0 ≤ p1 , p2 < 1 present the fraction of prey and predator which die due to harvesting or pesticides etc, q is the size of immigration or stock of the predator. Such a model is an impulsive differential equation whose theory and applications were greatly developed by the efforts of Bainov and Lakshmikantham et al. [7,8] and, moreover, the theory of impulsive differential equations is being recognized to be not only richer than the corresponding theory of differential equations without impulses, but also represents a more natural framework for the mathematical modeling of real world phenomenons. In recent years, models with sudden perturbations have been intensively researched, such as Lotka–Volterra [9], Holling-type [10–13], Ivlevtype [10,14] and Watt-type [15,16]. Most of the models mentioned above have dealt with impulsive harvesting and the immigration of predators at different fixed times. On the contrary, we consider the impulsive harvesting and immigration at the same time in our model which has not been studied well until now. The organization of this paper is as follows. In the next section, we introduce some notations and the properties of the single species equation with impulsive effect which is used in this paper. In Section 3, we show the local stability of preyfree periodic solution and the existence of a positive periodic solution. Moreover, we give a sufficient condition for the permanence of system (1.2) by applying the Floquet theory. In Section 4, we illustrate bifurcation diagrams with respect to various parameters to show dynamical varieties of the model (1.2) by using numerical simulations. Finally, we give a conclusion. 2. Preliminaries Firstly, we give some notations, definitions and Lemmas which will be useful for our main results. Let R+ = [0, ∞) and R2+ = {x = (x(t ), y(t )) ∈ R2 : x(t ), y(t ) ≥ 0}. Denote N the set of all of nonnegative integers and f = (f1 , f2 )T the right hand of (1.2). Let V : R+ × R2+ → R+ , then V is said to be in a class V0 if (1) V is continuous on (nT , (n + 1)T ] × R2+ , and lim (t ,y)→(nT ,x) V (t , y) = V (nT + , x) exists. (2) V is locally Lipschitzian in x.

t >nT

Definition 2.1. Let V ∈ V0 , (t , x) ∈ (nT , (n + 1)T ]× R2+ . The upper right derivatives of V (t , x) with respect to the impulsive differential system (1.2) are defined as D+ V (t , x) = lim sup h→0+

1 h

[V (t + h, x + hf (t , x)) − V (t , x)].

Remark 2.2. (1) The solution of the system (1.2) is a piecewise continuous function x : R+ → R2+ , x(t ) is continuous on (nT , (n + 1)T ], n ∈ N and x(nT + ) = limt →nT + x(t ) exists. (2) The smoothness properties of f guarantees the global existence and uniqueness of solution of the system (1.2). (See [8] for the details.) Since dx = dt = 0 whenever x(t ) = y(t ) = 0, t 6= nT and x(nT + ) = (1 − p1 )x(nT ), y(nT + ) = (1 − p2 )y(nT ) + q(0 ≤ pi < dt 1, i = 1, 2, n ≥ 0). We have the following lemma. dy

Lemma 2.3. Let x(t ) = (x(t ), y(t )) be a solution of (1.2). (1) If x(0+ ) ≥ 0 then x(t ) ≥ 0 for all t ≥ 0. (2) If x(0+ ) > 0 then x(t ) > 0 for all t ≥ 0. We will use the following important comparison theorem on an impulsive differential equation [8].

1314

H.K. Baek / Nonlinear Analysis: Real World Applications 11 (2010) 1312–1322

Lemma 2.4 (Comparison theorem). Suppose V ∈ V0 and



D+ V (t , x) ≤ g (t , V (t , x)), t 6= nτ V (t , x(t + )) ≤ ψn (V (t , x)), t = nτ ,

(2.1)

g : R+ × R+ → R is continuous on (nτ , (n + 1)τ ] × R+ and for u(t ) ∈ R+ , n ∈ N, lim(t ,y)→(nτ + ,u) g (t , y) = g (nτ + , u) exists, ψn : R+ → R+ is non-decreasing. Let r (t ) be the maximal solution of the scalar impulsive differential equation

 0 u (t ) = g (t , u(t )), t 6= nτ , u(t + ) = ψn (u(t )), t = nτ , u(0+ ) = u , 0

(2.2)

existing on [0, ∞). Then V (0+ , x0 ) ≤ u0 implies that V (t , x(t )) ≤ r (t ), t ≥ 0, where x(t ) is any solution of (2.1). Similar results can be obtained when all the conditions of the inequalities in the Lemma 2.4 are reversed. Note that if we have some smoothness conditions of g (t , u(t )) to guarantee the existence and uniqueness of the solutions for (2.2), then r (t ) is exactly the unique solution of (2.2). Now, we give the basic properties of the following impulsive differential equation.

 0 y (t ) = −Dy(t ), t 6= nT y(t + ) = (1 − p2 )y(t ) + q, y(0+ ) = y . 0

t = nT ,

(2.3)

Then we can easily obtain the following results. q exp(−D(t −nT )) q Lemma 2.5. (1) y∗ (t ) = 1−(1−p ) exp(−DT ) , t ∈ (nT , (n + 1)T ], n ∈ N and y∗ (0+ ) = 1−(1−p ) exp(−DT ) is a positive periodic 2 2 solution of (2.3).   (−DT ) (2) y(t ) = (1 − p2 )n+1 y0 − 1−(1q−exp exp(−Dt ) + y∗ (t ) is the solution of (2.3) with y0 ≥ 0, t ∈ (nT , (n + 1)T ] p ) exp(−DT ) 2

and n ∈ N. (3) All solutions y(t ) of (1.2) with y0 ≥ 0 tend to y∗ (t ). i.e., |y(t ) − y∗ (t )| → 0 as t → ∞. It is from Lemma 2.5 that the general solution y(t ) of (2.3) can be synchronized with the positive periodic solution y∗ (t ) of (2.3) for sufficiently large t and we can obtain the complete expression for the prey-free periodic solution of (1.2)



(0, y (t )) = 0, ∗

q exp(−D(t − nT ))



1 − (1 − p2 ) exp(−DT )

for t ∈ (nT , (n + 1)T ].

To study the stability of the prey-free periodic solution as a solution of (1.2) we present the Floquet theory for the linear T -periodic impulsive equation:

(

dx

= A(t )x(t ), t 6= τk , t ∈ R, dt + x(t ) = x(t ) + Bk x(t ), t = τk , k ∈ Z.

(2.4)

Then we introduce the following conditions: (H1 ) A(·) ∈ PC (R, C n×n ) and A(t + T ) = A(t )(t ∈ R), where PC (R, C n×n ) is the set of all piecewise continuous matrix functions which is left continuous at t = τk , and C n×n is the set of all n × n matrices. (H2 ) Bk ∈ C n×n , det(E + Bk ) 6= 0, τk < τk+1 (k ∈ Z), (H3 ) There exists a q ∈ N such that Bk+q = Bk , τk+q = τk + T (k ∈ Z). Let Φ (t ) be a fundamental matrix of (2.4), then there exists a unique nonsingular matrix M ∈ C n×n such that

Φ (t + T ) = Φ (t )M (t ∈ R).

(2.5)

By equality (2.5) there corresponds to the fundamental matrix Φ (t ) and the constant matrix M which we call the monodromy matrix of (2.4) (corresponding to the fundamental matrix of Φ (t )). All monodromy matrices of (2.4) are similar and have the same eigenvalues. The eigenvalues µ1 , . . . , µn of the monodromy matrices are called the Floquet multipliers of (2.4). Lemma 2.6 ([7] Floquet theory). Let conditions (H1 )–(H3 ) hold. Then the linear T -periodic impulsive Eq. (2.4) is (1) stale if and only if all multipliers µj (j = 1, . . . , n) of (2.4) satisfy the inequality |µj | ≤ 1, and moreover, to those µj for which |µj | = 1, there correspond simple elementary divisors; (2) asymptotically stable if and only if all multipliers µj (j = 1, . . . , n) of (2.4) satisfy the inequality |µj | < 1; (3) unstable if |µj | > 1 for some j = 1, . . . , n.

H.K. Baek / Nonlinear Analysis: Real World Applications 11 (2010) 1312–1322

1315

3. Main theorems First, we show that all solutions of (1.2) are uniformly ultimately bounded. Proposition 3.1. There is an M > 0 such that x(t ), y(t ) ≤ M for all t large enough, where (x(t ), y(t )) is a solution of (1.2). Proof. Let x(t ) = (x(t ), y(t )) be a solution of (1.2) and let V (t , x) = D+ V + β V = −

ma eK

m

x(t )2 +

e

m x e

(t ) + y(t ). Then V ∈ V0 , if t 6= nT

(a + β)x(t ) + (β − D)y(t ).

(3.1)

Clearly, the right hand of (3.1), is bounded when 0 < β < D. When t = nT , V (nT + ) = me x(nT + ) + y(nT + ) = (1 − p1 ) me x(nT ) + (1 − p2 )y(nT ) + q ≤ V (nT ) + q. So we can choose 0 < β0 < D and M0 > 0 such that



D+ V ≤ −β0 V + M0 , t 6= nT , V (nτ + ) ≤ V (nτ ) + q.

(3.2)

From Lemma 2.4, we can obtain that V (t ) ≤



V (0+ ) −

M0



β0

exp (−β0 t ) +

q(1 − exp(−(n + 1)β0 T )) 1 − exp(−β0 T )

exp(−β0 (t − nT )) +

M0

β0

for t ∈ (nT , (n + 1)T ]. Therefore, V (t ) is bounded by a constant for sufficiently large t. Hence there is an M > 0 such that x(t ) ≤ M , y(t ) ≤ M for a solution (x(t ), y(t )) with all t large enough.  Now, we present a condition which guarantees locally asymptotical stability of the prey-free periodic solution (0, y∗ (t )). Theorem 3.2. If aT +

e

bq exp(−DT ) + c (1 − (1 − p2 ) exp(−DT ))

ln

bD

bq + c (1 − (1 − p2 ) exp(−DT ))

< ln

1 1 − p1

,

then (0, y∗ (t )) is locally asymptotically stable Proof. The local stability of the periodic solution (0, y∗ (t )) of (1.2) may be determined by considering the behavior of small amplitude perturbations of the solution. Let (x(t ), y(t )) be any solution of (1.2). Define x(t ) = u(t ), y(t ) = y∗ (t ) + v(t ). Then they may be written as



u(t ) v(t )



u(0) = Φ (t ) , v(0)





0 ≤ t ≤ T,

(3.3)

where Φ (t ) satisfies ey∗ (t )

 dΦ dt

a−

 =

0

by∗ (t ) + c my∗ (t )

by∗ (t )

   Φ (t )

(3.4)

−D

+c

and Φ (0) = I, the identity matrix. The linearization of the third and fourth equation of (1.2) becomes



u(nT + ) v(nT + )



 =

1 − p1 0

1 − p2

Note that all eigenvalues of S =



1 − p1 0

p2 ) exp(−dT ) < 1. Since y∗ (t )

T

Z 0

by∗ (t )

+c

dt = −

1 bD

ln

u(nT ) . v(nT )



0



0 1 − p2



(3.5)

R T Φ (T ) are µ1 = (1 − p1 ) exp 0 a −

bq exp(−DT ) + c (1 − (1 − p2 ) exp(−DT )) bq + c (1 − (1 − p2 ) exp(−DT ))

,

we have



µ1 = (1 − p1 ) exp aT +

e bD

ln

bq exp(−DT ) + c (1 − (1 − p2 ) exp(−DT )) bq + c (1 − (1 − p2 ) exp(−DT ))

By Lemma 2.6, (0, y (t )) is locally asymptotically stable if |µ1 | < 1.i.e., ∗

aT +

e bD

ln

bq exp(−DT ) + c (1 − (1 − p2 ) exp(−DT )) bq + c (1 − (1 − p2 ) exp(−DT ))

< ln

1 1 − p1

. 

 .

ey∗ (t ) dt by∗ (t )+c



and µ2 = (1 −

1316

H.K. Baek / Nonlinear Analysis: Real World Applications 11 (2010) 1312–1322

We investigate the permanence of the system (1.2) in the next theorem. We mention the following definition before stating the permanent theorem. Definition 3.3. The system (1.2) is permanent if there exist M ≥ m > 0 such that, for any solution (x(t ), y(t )) of (1.2) with x0 > 0, m ≤ lim inf x(t ) ≤ lim sup x(t ) ≤ M and m ≤ lim inf y(t ) ≤ lim sup y(t ) ≤ M . t →∞

t →∞

t →∞

t →∞

Theorem 3.4. The system (1.2) is permanent if aT +

e bD

ln

bq exp(−DT ) + c (1 − (1 − p2 ) exp(−DT )) bq + c (1 − (1 − p2 ) exp(−DT ))

> ln

1 1 − p1

.

Proof. Let (x(t ), y(t )) be any solution of (1.2)with x0 > 0. From Proposition 3.1, we may assume that x(t ) ≤ M, y(t ) ≤ M, q exp(−DT ) t ≥ 0 and M > ace . Let m2 = 1−(1−p ) exp(−DT ) − 2 , 2 > 0. From Lemma 2.5, clearly we have y(t ) ≥ m2 for all t that are large 2 enough. Now we shall find an m1 > 0 such that x(t ) ≥ m1 for all t that are large enough. We will do this in the following two steps. (Step 1) Since aT >

e bD

ln

bq exp(−DT ) + c (1 − (1 − p2 ) exp(−DT )) bq + c (1 − (1 − p2 ) exp(−DT ))

+ ln

1 1 − p1

we can choose m3 > 0, 1 > 0 to be small enough such that δ = bq exp(−DT )+c (1−(1−p2 ) exp(−DT )) bq+c (1−(1−p2 ) exp(−DT ))

e1 T c

,

mm3 b+m3

 < D and R = (1 − p1 ) exp aT −

a K

Tm3 −



> 1. Suppose that x2 (t ) < m3 for all t. Then we get y0 (t ) ≤ y(t )(−D + δ) from the above assumptions. By Lemma 2.4, we have y(t ) ≤ u(t ) and u(t ) → u∗ (t ), t → ∞, where u(t ) is the solution of  0 u (t ) = (−D + δ)u(t ), t 6= nT , u(t + ) = (1 − p2 )u(t ) + q, t = nT (3.6) u(0+ ) = y , 0 e bD

ln

−

q exp((−D+δ)(t −nT )) and u∗ (t ) = 1−(1−p ) exp((−D+δ)T ) , t ∈ (nT , (n + 1)T ]. Then there exists T1 > 0 such that y(t ) ≤ u(t ) ≤ u∗ (t ) + 1 and 2



ax(t )



K am3



K am3

x0 (t ) = x(t ) a −

≥ x(t ) a − ≥ x(t ) a −

ey(t )

−

− −



by(t ) + x(t ) + c ey(t ) 

by(t ) + c eu∗ (t ) + e1



+ b 1 + c e1  ≥ x(t ) a − − ∗ − , t ≥ T1 . K bu (t ) + c c 

K am3

bu∗ (t )

eu∗ (t )

Let N1 ∈ N and N1 T ≥ T1 . We have, for n ≥ N1

 



x0 (t ) ≥ x(t ) a −

a

m3 −

eu∗ (t )

K bu∗ (t ) + c x(t + ) = (1 − p )x(t ), t = nT . 1

−

e 1  c

,

t 6= nT ,

(3.7)

Integrating (3.7) on (nT , (n + 1)T ](n ≥ N1 ), we obtain

Z

x((n + 1)T ) ≥ x(nT + ) exp

(n+1)T

a− nT

a K

m3 −

eu∗ (t ) bu∗ (t ) + c

−

e 1 c

dt



= x(nT )R.

Then x((N1 + k)T ) ≥ x(N1 T )Rk → ∞ as k → ∞ which is a contradiction. Hence there exists a t1 > 0 such that x(t1 ) ≥ m3 . (Step 2) If x(t ) ≥ m3 for all t ≥ t1 , then we are done. If not, we may let t ∗ = inft >t1 {x(t ) < m3 }. Then x(t ) ≥ m3 for t ∈ [t1 , t ∗ ] and, by the continuity of x(t ), we have x(t ∗ ) = m3 . In this step, we have only to consider two possible cases. Case 1) t ∗ = n1 T for some n1 ∈ N. Then (1 − p1 )m3 ≤ x(t ∗+ ) = (1 − p1 )x(t ∗ ) < m3 . Select n2 , n3 ∈ N such that

(n2 − 1)T >

1 ln( M + ) q

−D+δ

and (1 − p1 )n2 Rn3 exp(n2 σ T ) > (1 − p1 )n2 Rn3 exp((n2 + 1)σ T ) > 1, where σ = a −

a K

m3 − ce M < 0.

H.K. Baek / Nonlinear Analysis: Real World Applications 11 (2010) 1312–1322

1317

Let T 0 = n2 T + n3 T . In this case we will show that there exists t2 ∈ (t ∗ , t ∗ + T 0 ] such that x(t2 ) ≥ m3 . Otherwise, by (3.6) with u(t ∗+ ) = y(t ∗+ ), we have u(t ) = (1 − p2 )

n 1 +1



u(t

∗+

q exp((−D + δ)T )

)−

1 − (1 − p2 ) exp((−D + δ)T )



exp((−D + δ)(t − t ∗ )) + u∗ (t )

for (n − 1)T < t ≤ nT and n1 + 1 ≤ n ≤ n1 + 1 + n2 + n3 . So we get |u(t )− u∗ (t )| ≤ (M + q) exp((−D +δ)(t − t ∗ )) < 1 and y(t ) ≤ u(t ) ≤ u∗ (t ) + 1 for t ∗ + n2 T ≤ t ≤ t ∗ + T 0 . Then we know that the Eq. (3.7) also holds for t ∈ [t ∗ + n2 T , t ∗ + T 0 ]. As in step 1, we have x(t ∗ + T 0 ) ≥ x(t ∗ + n2 T )Rn3 . Since y(t ) ≤ M, we have

(

x0 (t ) ≥ x(t ) a −

a

m3 −

K x(t ) = (1 − p1 )x(t ), +

e c

M = σ x(t ),




t 6= nT ,

(3.8)

t = nT ,

for t ∈ [t ∗ , t ∗ + n2 T ]. Integrating (3.8) on [t ∗ , t ∗ + n2 T ] we have x((t ∗ + n2 T )) ≥ m3 exp(σ n2 T )

≥ m3 (1 − p1 )n2 exp(σ n2 T ) > m3 . Thus x(t ∗ + T 0 ) ≥ m3 (1 − p1 )n2 exp(σ n2 T )Rn3 > m3 which is a contradiction. Now, let t¯ = inft >t ∗ {x(t ) ≥ m3 }. Then x(t ) ≤ m3 for t ∗ ≤ t < t¯ and x(t¯) = m3 . For t ∈ [t ∗ , t¯), suppose t ∈ (t ∗ + (k − 1)T , t ∗ + kT ], k ∈ N and k ≤ n2 + n3 . So, we have, for t ∈ [t ∗ , t¯), from (3.8) we obtain x(t ) ≥ x(t ∗+ )(1 − p1 )k−1 exp((k − 1)σ T ) exp(σ (t − (t ∗ + (k − 1)T ))) ≥ m3 (1 − p1 )k exp(kσ T ) ≥ m3 (1 − p1 )n2 +n3 exp(σ (n2 + n3 )T ) ≡ m01 . Case (2) t ∗ 6= nT , n ∈ N. Then x(t ) ≥ m3 for t ∈ [t1 , t ∗ ) and x(t ∗ ) = m3 . Suppose that t ∗ ∈ (n01 T , (n01 + 1)T ) for some 0 n1 ∈ N. There are two possible cases. Case(2(a)) x(t ) < m3 for all t ∈ (t ∗ , (n01 + 1)T ]. In this case we will show that there exists t2 ∈ [(n01 + 1)T , (n01 + 1)T + T 0 ] such that x(t2 ) ≥ m3 . Suppose not. i.e., x(t ) < m3 , for all t ∈ [(n01 + 1)T , (n01 + 1 + n2 + n3 )T ]. Then x(t ) < m3 for all t ∈ (t ∗ , (n01 + 1 + n2 + n3 )T ]. By (3.6) with u((n01 + 1)T + ) = y((n01 + 1)T + ), we have



u(t ) = u((n01 + 1)T + ) −

q exp(−D + δ)



1 − (1 − p2 ) exp(−D + δ)

exp((−D + δ)(t − (n01 + 1)T )) + u∗ (t )

for t ∈ (nT , (n + 1)T ], n01 + 1 ≤ n ≤ n01 + n2 + n3 . A similar argument as in step 1, we have x((n01 + 1 + n2 + n3 )T ) ≥ x((n01 + 1 + n2 )T )Rn3 . It follows from (3.8) that x((n01 + 1 + n2 )T ) ≥ m3 (1 − p1 )n2 +1 exp(σ (n2 + 1)T ). Thus x((n01 + 1 + n2 + n3 )T ) ≥ m3 (1 − p1 )n2 +1 exp(σ (n2 + 1)T )Rn3 > m3 which is a contradiction. Now, let t¯ = inft >t ∗ {x(t ) ≥ m3 }. Then x(t ) ≤ m3 for t ∗ ≤ t < t¯ and x(t¯) = m3 . For t ∈ [t ∗ , t¯), suppose t ∈ (n01 T + (k0 + 1)T , n01 T + k0 T ], k0 ∈ N, k0 ≤ 1 + n2 + n2 , we have x(t ) ≥ m3 (1 − p)1+n2 +n3 exp(σ (1 + n2 + n3 )T ) ≡ m1 . Since m1 < m01 , so x(t ) ≥ m1 for t ∈ (t ∗ , t¯). Case (2(b)) There is a t 0 ∈ (t ∗ , (n01 + 1)T ] such that x2 (t 0 ) ≥ m3 . Let tˆ = inft >t ∗ {x(t ) ≥ m3 }. Then x(t ) ≤ m3 for t ∈ [t ∗ , tˆ) and x(tˆ) = m3 . Also, (3.8) holds for t ∈ [t ∗ , tˆ). Integrating the equation on [t ∗ , t )(t ∗ ≤ t ≤ tˆ), we can get that x(t ) ≥ x(t ∗ ) exp(σ (t − t ∗ )) ≥ m3 exp(σ T ) ≥ m1 . Thus in both case the similar argument can be continued since x(t ) ≥ m1 for some t > t1 . This completes the proof. 

R T

Remark 3.5. Let F (T ) = (1 − p1 ) exp ec (y∗ (0)−y∗ (T )) (by∗ (T )+c )(by∗ (0)+c )



0

a −

ey∗ (t ) dt by∗ (t )+c



R T − 1. Then F 0 (T ) = (1 − p1 ) exp 0 a −

ey∗ (t ) dt by∗ (t )+c



a +

q(1−exp(−DT ))

. Since y∗ (0) − y∗ (T ) = 1−(1−p ) exp(−DT ) > 0, F (T ) is strictly increasing. So, there exists a unique solution 2 T ∗ of the equation F (T ) = 0 because F (0) = −p1 < 0 and limT →∞ F (T ) = ∞. From Theorems 3.2 and 3.4, we know that T ∗ plays a role in a critical value that discriminates between stability and permanence. Now, we deal with problem of the bifurcation of a nontrivial periodic solution of the system (1.2), near (0, y∗ (t )). The following theorem establishes the existence of a positive periodic solution of the system (1.2) near (0, y∗ (t )). Theorem 3.6. When T = T ∗ , the prey-free solution (0, y∗ (t )) bifurcates to another periodic solution which is supercritical if Ke ≤ 4abc.

1318

H.K. Baek / Nonlinear Analysis: Real World Applications 11 (2010) 1312–1322

Proof. To prove this theorem, we will use the results of [17,18,10]. It is convenient for the computation to exchange the subscripts of x and y and change the period T to τ . Thus the system (1.2) becomes as follows:

 my(t )x(t )  , x0 (t ) = −Dx(t ) +    bx ( t ) + y(t ) + c    y ( t ) cy ( t )x(t ) y0 (t ) = ay(t ) 1 − − ,  K bx(t ) + y(t ) + c    +   x(t +) = (1 − p2 )x(t ) + q, t = nτ . y(t ) = (1 − p1 )y(t ),

   

t 6= nτ , (3.9)

  

Let Φ be the flow associated with (3.9). We have X (t ) = Φ (t , x0 ), 0 < t ≤ τ , where x0 = (x(0), y(0)). The flow Φ applies up to time τ . So, X (τ ) = Φ (τ , x0 ). We will use all notations in [18]. Note that  cyx myx y − F1 (x, y) = −Dx + , F2 (x, y) = ay 1 − , bx + y + c K bx + y + c

Θ1 (x, y) = (1 − p2 )x + q, Θ2 (x, y) = (1 − p1 )y, ζ (t ) = (y∗ (t ), 0), Z t  Z t  ∂ Φ1 (t , x0 ) ∂ F1 (ζ (r )) ∂ F2 (ζ (r )) ∂ Φ2 (t , x0 ) = exp dr , = exp dr , ∂x ∂x ∂y ∂y 0 0 Z t  Z ν  Z t ∂ Φ1 (t , x0 ) ∂ F1 (ζ (r )) ∂ F1 (ζ (ν)) ∂ F2 (ζ (r )) = exp dr exp dr dν. ∂y ∂x ∂y ∂y 0 ν 0 Then

∂ Θ1 ∂ Θ2 = = 0, ∂y ∂x

∂ Θ1 = 1 − p2 , ∂x

∂ Θ2 = 1 − p1 , ∂y

∂ 2 Θi = 0, ∂ x∂ y

i = 1, 2

and

∂ 2 Θ2 = 0. ∂ y2

Now, we can compute d00 = 1 −



∂ Θ2 ∂ Φ2 · ∂y ∂y

 (τ0 ,x0 )

τ0

Z = 1 − (1 − p1 ) exp

a− 0

ey∗ (r ) by∗ (r ) + c



dr ,

where τ0 is the root of d00 = 0. Actually, we easily see that τ0 = T ∗ .



0

a0 = 1 −

∂ Θ1 ∂ Φ1 · ∂x ∂x

 (τ0 ,x0 )

= 1 − (1 − p2 ) exp(−Dτ0 ) > 0,

Clearly.|1 − a00 | < 1.

   ∂ Θ1 ∂ Φ2 ∂ Φ1 ∂ Θ1 ∂ Φ1 · + · = −(1 − p2 ) ∂x ∂y ∂y ∂ y (τ0 ,x0 ) ∂ y (τ0 ,x0 ) Z τ0  Z ν  Z τ0 ∂ F1 (ζ (r )) ey∗ (ν) ∂ F2 (ζ (r )) = −(1 − p2 ) exp dr exp dr dν < 0. ∂x by∗ (ν) + c ∂y 0 ν 0 Z  Z  Z t t ν ∂ 2 Φ2 (t , x0 ) ∂ F2 (ζ (r )) ∂ 2 F2 (ζ (ν)) ∂ F2 (ζ (r )) = exp dr exp dr dν, ∂ y∂ x ∂ y ∂ y ∂ x ∂y 0 ν 0 Z Z   Z ν τ0 τ0 ∂ F2 (ζ (r )) ∂ F2 (ζ (r )) ∂ 2 Φ2 (τ0 , x0 ) eby∗ (ν) exp =− exp dr dr dν < 0. ∂ y∂ x ∂y (by∗ (ν) + c )3 ∂y 0 0 ν Z ν Z t  2  Z t ∂ 2 Φ2 (t , x0 ) ∂ F2 (ζ (r )) ∂ F2 (ζ (ν)) ∂ F2 (ζ (r )) = exp dr exp dr dν 2 ∂ y2 ∂y ∂y 0  ∂ y2  0 Z t  νZ t ∂ F2 (ζ (r )) ∂ F2 (ζ (ν)) + exp dr ∂ y Z0 ν Zν ν  ∂ y∂ x Z θ   ∂ F1 (ζ (r )) ∂ F1 (ζ (θ )) ∂ F2 (ζ (r )) × dr exp dr dθ dν. exp ∂x ∂y ∂y 0 0 θ Z   Z  Z τ0 τ0 ν 2a 2ey∗ (ν) ∂ 2 Φ2 (τ0 , x0 ) ∂ F2 (ζ (r )) ∂ F2 (ζ (r )) = exp dr − + exp dr dν ∂ y2 (by∗ (ν) + c )2 ∂y 0 0 K Z τ0  νZ τ0 ∂ y ∗ ∂ F2 (ζ (r )) eby (ν) − exp dr ∗ 3 ∂y Z0 ν Z νν  (by∗ (ν) + c ) Z θ   ey (θ ) ∂ F2 (ζ (r )) ∂ F1 (ζ (r )) × exp dr exp dr dθ dν < 0 . ∂x by∗ (θ ) + c ∂y 0 0 θ Z t  ∂ 2 Φ2 (t , x0 ) ∂ F2 (ζ (t )) ∂ F2 (xp (r ), 0) = exp dr . ∂ y∂τ ∂y ∂y 0 

b00 = −

H.K. Baek / Nonlinear Analysis: Real World Applications 11 (2010) 1312–1322

1319

τ0 ∂ F2 (xp (r ), 0) ∂ 2 Φ2 (τ0 , x0 ) ∂ F2 (ζ (τ0 )) = exp dr ∂ y∂τ ∂y ∂y 0 Z   τ0 ey∗ (r ) ey∗ (τ0 )  a− ∗ exp = a− ∗ dr . by (τ0 ) + c by (r ) + c 0 ∂ Φ1 (τ0 , x0 ) ) q exp (− D τ 0 = y˙∗ (τ0 ) = − < 0. ∂τ D(1 − (1 − p2 ) exp(−Dτ0 ))   ∂ 2 Θ2 ∂ Φ1 (τ0 , x0 ) ∂ Φ1 (τ0 , x0 ) 1 ∂ Θ1 ∂ Φ1 (τ0 , X0 ) ∂ Φ2 (τ0 , x0 ) + · 0 · · B = − ∂ x∂ y ∂τ ∂x a0 ∂x ∂τ ∂y  2  ∂ Θ2 ∂ Φ2 (τ0 , x0 ) 1 ∂ Θ1 ∂ Φ1 (τ0 , x0 ) ∂ 2 Φ2 (τ0 , x0 ) − · 0 · · + ∂y ∂ x∂ y a0 ∂x ∂τ ∂τ ∂ y   Z τ0   2 ey∗ (r ) ∂ Φ1 (τ0 , x0 ) ey∗ (τ0 )  ∂ Φ2 (τ0 , x0 ) 1 a− ∗ · 0 · (1 − p2 ) · + a− ∗ exp dr . = −(1 − p1 ) ∂ x∂ y a0 ∂τ by (τ0 ) + c by (r ) + c 0  0  2 ∂ Θ2 b ∂ Φ1 (τ0 , x0 ) ∂ Φ1 (τ0 , x0 ) ∂ Φ2 (τ0 , x0 ) C = −2 − 00 · + ∂ x∂ y a0 ∂x ∂y ∂y  2 ∂ 2 Θ2 ∂ Φ2 (τ0 , x0 ) ∂ Θ2 b00 ∂ 2 Φ2 (τ0 , x0 ) ∂ Θ2 ∂ 2 Φ2 (τ0 , x0 ) − +2 · − · · ∂ y2 ∂y ∂ y a00 ∂ y∂ x ∂y ∂ y2 0 2 2 ∂ Φ2 b ∂ Φ2 − (1 − p1 ) 2 . = 2(1 − p1 ) 00 a0 ∂ x ∂ y ∂y

Z



ey∗ (t )

ecq exp(−Dt )

First we determine the sign of B. For this, let φ(t ) = a− by∗ (t )+c . Then we obtain that φ 0 (t ) = D(by∗ (t )+c )2 (1−(1−p ) exp(−Dτ )) > 0 2 0 Rτ bq exp(−Dτ0 )+c (1−(1−p2 ) exp(−Dτ0 )) e 1 and so φ(t ) is strictly increasing. Since 0 0 φ(t )dt = aτ0 + bD ln = ln > 0, we have bq+c (1−(1−p ) exp(−Dτ )) 1−p 2

0

1

(ν) ∗ φ(τ0 ) > 0. This implies that B < 0. To show C > 0 we need to prove that − Ka + (byey ∗ (ν)+c )2 ≤ 0. Let t = y (ν) > 0. Then t a et does not depend on the variables a, b, c , e and K . Moreover, it is easily shown that − K + (bt +c )2 ≤ 0 if Ke ≤ 4abc. So C > 0. Hence BC < 0. Thus it is from Theorem 2 of [18] that the theorem holds.  ∗

4. Numerical simulations In Section 3, we have shown that there exists a locally asymptotically stable prey-free periodic solution of the model (1.2) under the condition of Theorem 3.2 and established the conditions in Theorem 3.4 that the model (1.2) is permanent. Now, consider the following choice of parametric values in order to substantiate the our theoretical results: a = 7, p2 = 0.01

b = 0.001, and

q = 5.

c = 1,

D = 0.2,

e = 1.1,

K = 1.65,

m = 1.045,

p1 = 0.1, (4.1)

Throughout this section, we take an initial value (x0 , y0 ) = (1, 1). It follows from Theorems 3.2 and 3.4 that the prey-free periodic solution is locally asymptotically stable if T < T ∗ ≈ 3.9142 (Fig. 1) and the prey and the predator can coexist if T > T ∗ (Fig. 2(a)). Moreover, we display typical bifurcation diagrams for the model (1.2) with respect to parameter T in Fig. 3. It is clear from Fig. 3 that the model undergoes a wide variety of dynamical behaviors including cycles (Fig. 2(a)), periodic doubling bifurcations (Fig. 4), chaotic bands (Fig. 2(b)), periodic windows, crises (the phenomenon of crisis in which chaotic attractors can suddenly appear or disappear, or change size discontinuously as a parameter smoothly varies), period-halving bifurcations (Fig. 5), etc. In fact, there is a crisis at T ≈ 8.645, which leads to a chaos band and then, the chaos disappears, a 4T -period solution of the model (1.2) appears (Fig. 6) Now, to exhibit the bifurcation diagrams of the model (1.2) with respect to a, we assume that all parameters except a are the same as in (4.1) and T = 6. Fig. 7 shows the bifurcation diagrams with respect to a in the range 4 < a < 20. From Theorems 3.2 and 3.4, we know that if a < a∗ ≈ 4.5777, then the prey-free periodic solution is locally stable, while the model (1.2) is permanent if a > a∗ . Further, Fig. 7 indicates that the growth rate a of the prey also have profound effects on the dynamics of the model (1.2). In Figs. 8 and 9, we illustrate the bifurcation diagrams of the model (1.2) with T = 6 and (4.1) with respect to parameters b and c, respectively. These bifurcation diagrams also show the rich dynamics including periodic windows, chaotic bands, crises, periodic doubling bifurcation, periodic halving bifurcation and so on. 5. Conclusions In this paper, we have investigated effects of impulsive perturbations on a predator–prey model with Beddington–DeAngelis functional response. We have proven that there exists a stable prey-free periodic solution when the impulsive

1320

H.K. Baek / Nonlinear Analysis: Real World Applications 11 (2010) 1312–1322

Fig. 1. Time series of a solution for T = 3.9. (a) x(t ). (b) y(t ).

Fig. 2. Phase portraits of the model (1.2). (a) A T -period solution for T = 3.92. (b) A chaotic solution for T = 6.

Fig. 3. Bifurcation diagrams of the model (1.2) with respect to T . (a) x(t ) is plotted. (b) y(t ) is plotted.

Fig. 4. Periodic doubling cascade. (a) Phase portrait of a 2T -period solution for T = 5.2. (b) Phase portrait of a 4T -period solution for T = 5.3.

Fig. 5. Periodic halving cascade. (a) Phase portrait of a 2T -period solution for T = 12. (b) Phase portrait of a T -period solution for T = 12.5.

Fig. 6. Crisis phenomenon. (a) A chaotic solution for T = 8.64. (b) A 4T -period solution for T = 8.65.

H.K. Baek / Nonlinear Analysis: Real World Applications 11 (2010) 1312–1322

1321

Fig. 7. Bifurcation diagrams of the model (1.2) with respect to a. (a) x(t ) is plotted. (b) y(t ) is plotted.

Fig. 8. Bifurcation diagrams of the model (1.2) with respect to b. (a) x(t ) is plotted. (b) y(t ) is plotted.

Fig. 9. Bifurcation diagrams of the model (1.2) with respect to c. (a) x(t ) is plotted. (b) y(t ) is plotted.

period T is less than the critical value T ∗ by using the Floquet theory of impulsive differential equation and small amplitude perturbation skills. In addition, it has been shown that the model (1.2) is permanent when the period T is greater than T ∗ via the comparison theorem. We have illustrated some bifurcation diagrams of the model (1.2) with respect to several parameters. Under the impulsive controlled model, our results are useful tools to control the prey population compared to the classical one. Actually, the model (1.1) is one of models with a classical biological control technique. It is easy to see that the preyfree equilibrium (0, y+ )(y+ > 0) of the model (1.1) does not exist and the equilibrium (0, 0) of the model (1.1) is unstable, which indicates that one cannot stamp out the prey steadily. However, from the impulsive controlled model (1.2), the prey can be eradicated if the impulsive period T is less than a critical value T ∗ . Thus, if one regards the prey as a pest, then the impulsive control strategy is more suitable than the classical one. Acknowledgement We would thank the referee for carefully reading the manuscript and suggesting improvements. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11]

J.R. Beddington, Mutual interference between parasites or predator and its effect on searching efficiency, J. Animal Ecol. 44 (1975) 331–340. D.L. DeAngelis, R.A. Goldstein, R.V. O’Neill, A model for trophic interaction, Ecology 56 (1975) 881–892. G.T. Skalski, J.F. Gilliam, Functional responses with predator interference: Viable alternatives to the Holling type II mode, Ecology 82 (2001) 3083–3092. M. Fan, Y. Kuang, Dynamics of a nonautonomous predator–prey system with the Beddington–DeAngelis functional response, J. Math. Anal. Appl. 295 (2004) 15–39. T.-W. Hwang, Uniqueness of limit cycles of the predator–prey system with Beddington–DeAngelis functional response, J. Math. Anal. Appl. 290 (2004) 113–122. J.M. Cushing, Periodic time-dependent predator–prey systems, SIAM J. Appl. Math. 32 (1977) 82–95. D.D. Bainov, P.S. Simeonov, Impulsive Differential Equations: Periodic Solutions and Application, in: Pitman Monographs and Surveys in Pure and Applied Mathematics, vol. 66, Longman Science & Technical, Harlo, UK, 1993. V. Lakshmikantham, D. Bainov, P. Simeonov, Theory of Impulsive Differential Equations, World Scientific Publisher, Singapore, 1989. B. Liu, Y. Zhang, L. Chen, Dynamic complexities in a Lotka–Volterra predator–prey model concerning impulsive control strategy, Internat. J. Bifur. Chaos 15 (2) (2005) 517–531. B. Liu, Y. Zhang, L. Chen, The dynamical behaviors of a Lotka–Volterra predator–prey model concerning integrated pest management, Nonlinearity Anal. 6 (2005) 227–243. O.D. Makinde, Solving ratio-dependent predator–prey system with constant effort harvesting using Adomian decomposition method, Appl. Math. Comput. 186 (2007) 17–22.

1322

H.K. Baek / Nonlinear Analysis: Real World Applications 11 (2010) 1312–1322

[12] D. Xiao, L.S. Jennings, Bifurcations of a ratio-dependent predator–prey system with constant rate harvesting, SIAM J. Appl. Math. 65 (3) (2005) 737–753. [13] S. Zhang, L. Chen, A study of predator–prey models with the Beddington–DeAngelis functional response and impulsive effect, Chaos Solitons Fractals 27 (2006) 237–248. [14] H. Wang, W. Wang, The dynamical complexity of a Ivlev-type prey-predator system with impulsive effect, Chaos Solitons Fractals (2007) doi:10.1016/j.chaos.2007.02.008. [15] H. Wang, W. Wang, The dynamical complexity of a Ivlev-type prey-predator system with impulsive effect, Chaos Solitons Fractals (2007) doi:10.1016/j.chaos.2007.02.008. [16] W. Wang, H. Wang, Z. Li, The dynamic complexity of a three-species Beddington-type food chain with impulsive control strategy, Chaos Solitons Fractals 32 (2007) 1772–1785. [17] S. Gakkhar, K. Negi, Pulse vaccination in SIRS epidemic model with non-monotonic incidence rate, Chaos Solitons Fractals 35 (2008) 626–638, 2007. [18] A. Lakmeche, O. Arino, Bifurcation of non trivial periodic solutions of impulsive differential equations arising chemotherapeutic treatment, Dyn. Contin. Discrete Impuls. Syst. 7 (2000) 265–287.