# Dynamics of Unilateral and Bilateral Control Systems with State Feedback for Renewable Resource Management

Authors: Mingzhan Huang; Shouzong Liu; Xinyu Song; Lansun Chen
Journal: Complexity
 (2020)

---

## Abstract

In this paper, mathematical models for the management of biological resources based on a given predator-prey relationship are proposed, and two types of control strategies, unilateral and bilateral control with impulsive state feedback, are studied. The existence of the order-1 homoclinic orbit, order-1 periodic solution, and bifurcation of homoclinic of the unilateral control system are obtained, and the attraction region of this system is also discussed. Besides, sufficient conditions for the existence and stability of the order-1 and order-2 periodic solutions of the bilateral control system are also gained. A series of numerical simulations including bifurcation diagrams of periodic solution are performed, which not only verify the theoretical results we get but also reveal some peculiar dynamical phenomena, such as the appearance of high-order periodic solutions and existence of parameter intervals with drastic order change of periodic solution. By comparing the two management strategies, our study encourages bilateral control rather than unilateral control for the risk of predator extinction.


---

## Body

## 1. Introduction

Reasonable exploitation of renewable resources has been a topic drawing close attention of researchers, resource owners, and resource management departments [1]. Enormous contributions have been made on sustainable management and exploitation of renewable resources [1–10]. However, the previous studies focused more on the fluctuating environment and different policies. Srinivasu et al. proposed optimization models for resource management and investigated the impact of periodic varying environment and harvest policies on the revenue [1]. Emilie Lindkvist et al. presented a new method to investigate the influence mechanism of several different factors to the management strategies of renewable resources by applying a computational form of adaptive management [11]. Stochastic factors and resource-constrained control tasks were also considered in the management process [12–14].A dynamic model for renewable resource management has played an important role in promoting the sustainable development of a region or a country due to its ability to anticipate future trend. In order to design management strategies of sustainable development, a variety of mathematical models have been constructed to study the impacts of different significant influence factors. Clark studied mathematical models of commercial fishing taking into consideration the economic and environmental factors [2]. Then with some investigations, Leung and Wang [15] presented a commercial fishing model incorporating capital investment. Ragozin and Brown [16] studied a model for a predator-prey system where the prey is supposed to have no commercial value while the predator can benefit humans. In contrast, Mesterton-Gibbons investigated a predator-prey system where the two species are both harvested with an optimal policy [17]. A nonautonomous model with periodic coefficients was proposed by Fan and Wang in [18], and they studied the positive periodic solution and its stability.Feedback control was also used in the general renewable resource management. Sandal and Steinshamn in [19] studied the optimal exploitation of renewable natural resources by a feedback model. They presented some feedback rules in computing optimal yields and depletion rates by applying perturbation theory. They also found the analytical expressions of the optimal harvest under an explicit feedback control law [12]. Fuertes et al. gave an optimal resource allocation policy aiming at maximizing the control performance and also provided experimental comparisons to show the priority of this policy [14].With the rapid development of theories and methods of the semicontinuous dynamic system in recent years, a large number of mathematical models of renewable resource management considering state feedback control are formulated and studied. Xiao and Dai studied a predator-prey system with Allee effect and state feedback and discussed the coexistence mode of the two biological resources [20]. Huang et al. constructed a population system with two types of harvesting, and they exploited the impact of the harvest effort on the development trend of population [21]. Liang et al. formulated a predator-prey system with mutual interference and also investigated the dynamical behaviors of the system [22]. Besides, Chen and Zhao proposed a delayed single-species model with impulsive control, and they mainly discussed the existence and stability of order-1 periodic solutions [23]. Furthermore, Guo et al. in [24] investigated a state feedback for the algae-fish system while Fu and Chen in [25] formulated an ecological system for water hyacinth and studied two types of state-dependent impulsive controls.Overexploitation or excessive use of a renewable resource that has economic value can be a threat to biodiversity. As we all know, all species do not exist in isolation, and they will interact with other species around them. A species may be the prey, or the predator, or a competitor of another species. So when it is overexploited, what will happen to the other one? And what should humans do to protect the biodiversity of the ecological system? Based on these considerations, in this paper by constructing semicontinuous dynamic systems, we present mathematical models for the management of two kinds of biological resources based on a given predator-prey relationship. Since the quantity of the prey population can be easily obtained by some modern devices and measurement techniques, we design two control strategies, unilateral and bilateral control, with the state feedback of the prey population.The paper is organized as follows. In Section2, we formulate two mathematical models for the unilateral and bilateral state feedback control of the two kinds of biological resources. In Section 3, we give qualitative analysis for the unilateral and bilateral control systems, and dynamical behaviors such as the existence and stability of order-1 and order-2 periodic solutions and bifurcation phenomenon corresponding to control parameters are investigated. In Section 4, series of numerical simulations including bifurcation diagrams of periodic solution are performed, which not only verify the theoretical results we get but also reveal some peculiar dynamical phenomena. We finish this paper with a brief discussion in Section 5.


## 2. Model Formulation

### 2.1. Basic Model

The basic model studied in this paper is the following predation system [26]:(1)dxtdt=rxt1−xtK−αxtyte−βxt,dytdt=ytμαxte−βxt−D,where x(t) and y(t) represent the densities of the prey and predator, respectively. r, K, α, β, μ, and D are the positive constants. αx(t)e−βx(t) denotes the functional response when the predation confronts group defense from the prey species.For model (1), there are always an unstable saddle point O(0, 0) and a boundary equilibrium E0(K, 0). To get positive equilibria, we need to solve the following algebra equations:(2)r1−xK−αye−βx=0,μαxe−βx−D=0.By simple calculation, we have the following lemma.Lemma 1.
(i)
Ifαμ < eβD, then system (1) has no positive equilibrium;(ii)
ifαμ = eβD and 1/β<K, then system (1) has a unique positive steady state 1/β,re/α1−1/Kβ;(iii)
ifαμ > eβD, then system (1) has at most two positive equilibria E1(x1, y1) and E2(x2, y2), where xi, i = 1, 2 satisfy μαxie−βxi=D and 0<x1<1/β<x2<K.
For case (iii) in Lemma1, Xiao and Ruan in [26] had discussed two cases: when system (1) has a global asymptotically stable equilibrium E1(x1, y1) and when it has an unstable equilibrium E1(x1, y1) and a saddle E2(x2, y2) (interested readers can refer Theorems 1 and 2 in [26]).Now, we give other conditions under which the existence of the two interior equilibria is guaranteed, andE1(x1, y1) is stable while E2(x2, y2) is still a saddle.Theorem 1.
Ifαμ > eβD and x2<K<1/β+x1, then system (1) has four equilibria: two unstable saddles O(0, 0) and E2(x2, y2), a stable node E0(K, 0), and a stable node or focus E1(x1, y1).Proof.
From Lemma1, we know that if αμ > eβD and x2 < K, the system (1) has two positive steady states E1(x1, y1) and E2(x2, y2) besides two boundary equilibria O(0, 0) and E0(K, 0). To study the stability of these equilibria, we give the Jacobi matrix of the linearized system of (1) in the following:(3)J=r−2rKx−αye−βx1−βx−αxe−βxαμye−βx1−βxαμxe−βx−D.
By simple calculation, we have(4)JO=r00−D,JE0=−r−αKe−βK0αKμe−βK−D,JEi=βrxiKK−xi−1β−αxie−βxiαμyie−βxi1−βxi0.
Obviously,O(0, 0) is an unstable saddle and E0(K, 0) is a stable node if x2 < K. At points Ei(xi, yi), i = 1, 2, the characteristic equation has the form(5)λ2−a1λ+a2=0,where a1=βrxi/KK−xi−1/βanda2=αDyie−βxi1−βxi.
Becausex1<1/β<x2<K<1/β+x1, we can get a2 < 0 when xi = x2 and a2 > 0 when xi = x1; then, E2(x2, y2) is a saddle point. Besides, xi = x1 implies a1 < 0, so E1(x1, y1) is a stable node or focus. The proof is completed.Theorem 2.
Ifrβ−1/K+αμ<2βrD/K, then there are no closed orbits for system (1) provided that the conditions listed in Theorem 1 are satisfied.Proof.
In consideration of the unstable saddleE2(x2, y2), if there is a close orbit (a periodic orbit or a homoclinic loop) of system (1), it must be in the following domain:(6)Ω1=x,y0<x≤x2,0<y<∞.
Let dt = eβxdτ, then system (1) is changed into an equivalent form:(7)dxdτ=rx1−xKeβx−αxy=Px,y,dydτ=yμαx−Deβx=Qx,y.
Applying Dulac functionB(x, y) = x−1, we have(8)∂PB∂x+∂QB∂y=αμ+eβx−βrKx+rβ−1K−Dx.
Note that−βr/Kx−D/x≤−2βrD/K, then we can get(9)∂PB∂x+∂QB∂y≤αμ+eβxrβ−1K−2βrDK.
In addition, we know thateβx ≥ 1, x ∈ (0, x2], and if rβ−1/K+αμ<2βrD/K, there is(10)∂PB∂x+∂QB∂y≤αμ+rβ−1K−2βrDK<0.
By Dulac theory, system (1) does not have a closed orbit in Ω1. This completes the proof.


### 2.2. Model Formulation

Li et al. [27] investigated impulsive control tactics of system (1) at periodic fixed moments, gained the existence and local stability of the prey-extinction periodic solution, and also exploited several kinds of bifurcation phenomena. He [28] studied system (1) incorporating impulsive state feedback controls. However, the study just focused on a special case, that is to say, system (1) has a unique positive equilibrium E1(x1, y1) which is globally stable.In this work, we study system (1) under the conditions provided in Theorems 1 and 2; that is to say, system (1) has two positive equilibria: a hyperbolic saddle E2(x2, y2) and a stable node or focus E1(x1, y1) (see Figure 1).Figure 1
Phase portrait of system (1) if conditions in Theorems 1 and 2 are satisfied.Assume that the prey and predator are two kinds of biological resources which have economic value, while the quantity of the prey population can be easily monitored by some modern devices and measurement techniques. Our research object is providing scientific strategies for the sustainable development and reasonable utilization of these two kinds of resources.From Figure1, we can see that if the prey density is less than xL (xL is the distance between the stable manifold of E2(x2, y2) and y-axis) or greater than x2, there is a high probability that the predator species goes extinct. So in the following, for the biodiversity of this ecological system, we choose a lower control threshold h1 > x1 > xL and an upper control threshold h2 < x2 for the prey density. To avoid overexploitation of the prey population, we first choose h1 as the control threshold. With the state feedback, once the density of the prey drops down to h1, then people release a certain amount of the prey and harvest the predator by a specific ratio at the same time. Then, we formulate a model for the unilateral control with state feedback of the prey population:(11)dxtdt=rxt1−xtK−αxtyte−βxt,dytdt=μαxtyte−βxt−Dyt,xt>h1,xt+=xt+m1,yt+=1−pyt,xt=h1,y>rh11−h1Ke−βh1,where x(0) > xL, y(0) > 0, x1 < h1 < h1 + m1 < x2, and 0 < p < 1.Besides, for any prey sizex = n, x1 < n < x2, the predator species itself has a sustainable range [nmi, nmx]. If the predator resource is overexploited, for example, if its density is less than nmi, then the prey size will exceed x2 and the predator will become extinct. For this reason, we set another upper control threshold h2 besides the lower one. Then, when the density of the predator is less than h2mi and the density of the prey increases to h2, people release a certain amount of the predator and harvest a certain amount of the prey at the same time. For this case, we formulate another management strategy with bilateral control:(12)dxtdt=rxt1−xtK−αxtyte−βxt,dytdt=μαxtyte−βxt−Dyt,h1<xt<h2,xt+=xt+m1,yt+=1−pt,xt=h1,y>rh11−h1Ke−βh1,xt+=xt−m2,yt+=yt+τ,xt=h2,y<rh21−h2Ke−βh2,with initial value h1 < x(0) < h2, y(0) > 0. Here, m1, m2, p, and τ are the positive constants that satisfy h1 + m1 = h2 − m2 ≜ n and x1 < n < x2.The aim of this paper is to investigate the management strategies that can maintain the sustainable development of the two species with rational replenishment and exploitation by studying systems (11) and (12).


## 2.1. Basic Model

The basic model studied in this paper is the following predation system [26]:(1)dxtdt=rxt1−xtK−αxtyte−βxt,dytdt=ytμαxte−βxt−D,where x(t) and y(t) represent the densities of the prey and predator, respectively. r, K, α, β, μ, and D are the positive constants. αx(t)e−βx(t) denotes the functional response when the predation confronts group defense from the prey species.For model (1), there are always an unstable saddle point O(0, 0) and a boundary equilibrium E0(K, 0). To get positive equilibria, we need to solve the following algebra equations:(2)r1−xK−αye−βx=0,μαxe−βx−D=0.By simple calculation, we have the following lemma.Lemma 1.
(i)
Ifαμ < eβD, then system (1) has no positive equilibrium;(ii)
ifαμ = eβD and 1/β<K, then system (1) has a unique positive steady state 1/β,re/α1−1/Kβ;(iii)
ifαμ > eβD, then system (1) has at most two positive equilibria E1(x1, y1) and E2(x2, y2), where xi, i = 1, 2 satisfy μαxie−βxi=D and 0<x1<1/β<x2<K.
For case (iii) in Lemma1, Xiao and Ruan in [26] had discussed two cases: when system (1) has a global asymptotically stable equilibrium E1(x1, y1) and when it has an unstable equilibrium E1(x1, y1) and a saddle E2(x2, y2) (interested readers can refer Theorems 1 and 2 in [26]).Now, we give other conditions under which the existence of the two interior equilibria is guaranteed, andE1(x1, y1) is stable while E2(x2, y2) is still a saddle.Theorem 1.
Ifαμ > eβD and x2<K<1/β+x1, then system (1) has four equilibria: two unstable saddles O(0, 0) and E2(x2, y2), a stable node E0(K, 0), and a stable node or focus E1(x1, y1).Proof.
From Lemma1, we know that if αμ > eβD and x2 < K, the system (1) has two positive steady states E1(x1, y1) and E2(x2, y2) besides two boundary equilibria O(0, 0) and E0(K, 0). To study the stability of these equilibria, we give the Jacobi matrix of the linearized system of (1) in the following:(3)J=r−2rKx−αye−βx1−βx−αxe−βxαμye−βx1−βxαμxe−βx−D.
By simple calculation, we have(4)JO=r00−D,JE0=−r−αKe−βK0αKμe−βK−D,JEi=βrxiKK−xi−1β−αxie−βxiαμyie−βxi1−βxi0.
Obviously,O(0, 0) is an unstable saddle and E0(K, 0) is a stable node if x2 < K. At points Ei(xi, yi), i = 1, 2, the characteristic equation has the form(5)λ2−a1λ+a2=0,where a1=βrxi/KK−xi−1/βanda2=αDyie−βxi1−βxi.
Becausex1<1/β<x2<K<1/β+x1, we can get a2 < 0 when xi = x2 and a2 > 0 when xi = x1; then, E2(x2, y2) is a saddle point. Besides, xi = x1 implies a1 < 0, so E1(x1, y1) is a stable node or focus. The proof is completed.Theorem 2.
Ifrβ−1/K+αμ<2βrD/K, then there are no closed orbits for system (1) provided that the conditions listed in Theorem 1 are satisfied.Proof.
In consideration of the unstable saddleE2(x2, y2), if there is a close orbit (a periodic orbit or a homoclinic loop) of system (1), it must be in the following domain:(6)Ω1=x,y0<x≤x2,0<y<∞.
Let dt = eβxdτ, then system (1) is changed into an equivalent form:(7)dxdτ=rx1−xKeβx−αxy=Px,y,dydτ=yμαx−Deβx=Qx,y.
Applying Dulac functionB(x, y) = x−1, we have(8)∂PB∂x+∂QB∂y=αμ+eβx−βrKx+rβ−1K−Dx.
Note that−βr/Kx−D/x≤−2βrD/K, then we can get(9)∂PB∂x+∂QB∂y≤αμ+eβxrβ−1K−2βrDK.
In addition, we know thateβx ≥ 1, x ∈ (0, x2], and if rβ−1/K+αμ<2βrD/K, there is(10)∂PB∂x+∂QB∂y≤αμ+rβ−1K−2βrDK<0.
By Dulac theory, system (1) does not have a closed orbit in Ω1. This completes the proof.


## 2.2. Model Formulation

Li et al. [27] investigated impulsive control tactics of system (1) at periodic fixed moments, gained the existence and local stability of the prey-extinction periodic solution, and also exploited several kinds of bifurcation phenomena. He [28] studied system (1) incorporating impulsive state feedback controls. However, the study just focused on a special case, that is to say, system (1) has a unique positive equilibrium E1(x1, y1) which is globally stable.In this work, we study system (1) under the conditions provided in Theorems 1 and 2; that is to say, system (1) has two positive equilibria: a hyperbolic saddle E2(x2, y2) and a stable node or focus E1(x1, y1) (see Figure 1).Figure 1
Phase portrait of system (1) if conditions in Theorems 1 and 2 are satisfied.Assume that the prey and predator are two kinds of biological resources which have economic value, while the quantity of the prey population can be easily monitored by some modern devices and measurement techniques. Our research object is providing scientific strategies for the sustainable development and reasonable utilization of these two kinds of resources.From Figure1, we can see that if the prey density is less than xL (xL is the distance between the stable manifold of E2(x2, y2) and y-axis) or greater than x2, there is a high probability that the predator species goes extinct. So in the following, for the biodiversity of this ecological system, we choose a lower control threshold h1 > x1 > xL and an upper control threshold h2 < x2 for the prey density. To avoid overexploitation of the prey population, we first choose h1 as the control threshold. With the state feedback, once the density of the prey drops down to h1, then people release a certain amount of the prey and harvest the predator by a specific ratio at the same time. Then, we formulate a model for the unilateral control with state feedback of the prey population:(11)dxtdt=rxt1−xtK−αxtyte−βxt,dytdt=μαxtyte−βxt−Dyt,xt>h1,xt+=xt+m1,yt+=1−pyt,xt=h1,y>rh11−h1Ke−βh1,where x(0) > xL, y(0) > 0, x1 < h1 < h1 + m1 < x2, and 0 < p < 1.Besides, for any prey sizex = n, x1 < n < x2, the predator species itself has a sustainable range [nmi, nmx]. If the predator resource is overexploited, for example, if its density is less than nmi, then the prey size will exceed x2 and the predator will become extinct. For this reason, we set another upper control threshold h2 besides the lower one. Then, when the density of the predator is less than h2mi and the density of the prey increases to h2, people release a certain amount of the predator and harvest a certain amount of the prey at the same time. For this case, we formulate another management strategy with bilateral control:(12)dxtdt=rxt1−xtK−αxtyte−βxt,dytdt=μαxtyte−βxt−Dyt,h1<xt<h2,xt+=xt+m1,yt+=1−pt,xt=h1,y>rh11−h1Ke−βh1,xt+=xt−m2,yt+=yt+τ,xt=h2,y<rh21−h2Ke−βh2,with initial value h1 < x(0) < h2, y(0) > 0. Here, m1, m2, p, and τ are the positive constants that satisfy h1 + m1 = h2 − m2 ≜ n and x1 < n < x2.The aim of this paper is to investigate the management strategies that can maintain the sustainable development of the two species with rational replenishment and exploitation by studying systems (11) and (12).


## 3. Main Results

### 3.1. Preliminary

Before the study of the above two systems, some notations and definitions of semidynamical systems are introduced. These will be used in the following discussion.Consider the general planar semidynamical system:(13)dxdt=Px,y,dydt=Qx,y,x,y∈ℳx,y,△x=αx,y,△y=βx,y,x,y∈ℳx,y,where (x, y) ∈ R2, P(x, y), Q(x, y), α(x, y), and β(x, y) ∈ C[R2, R]. ℳ⊂R2 is called the impulsive set, while φ: R2 ⟶ R2 is called the impulse mapping, where φZ=Z+=xZ+,yZ+=xZ+αxZ,yZ,yZ+βxZ,yZ for ZxZ,yZ∈ℳ. Z+ is the phase point of point Z, and N=φℳ is the phase set.AssumeX(t) = (x(t), y(t)) be an arbitrary solution of (5). Denote πt(X0, X) = X(t) to be the solution corresponding to initial value condition X(0) = X0. Then, its forward orbit can be described by π(X0, X) = {πt(X0, X) ∣ t ≥ 0}, which is also denoted by π(X) or π+(X) for simplicity. We denote Xk=Xtk+∈πX, where X(tk) ∈ π(X), k = 1, 2, ….Definition 1 (see [29, 30]).
An orbitπ(X) of (13) is said to be a periodic orbit if there is a k ∈ N+ such that Xk = X0. Denote k0 = min{k ∈ N+, Xk = X0}, then the orbit π(X) is called an order-k0 periodic orbit.Definition 2 (see [29, 30]).
Assumeπ(X) and πX˜ are the two orbits of (13), if there is a monotone smooth function ν(t), t > 0 such that πtX−πνtX^ϵ,t0, then we called the two orbits as ϵ-close.Definition 3 (see [29, 30, 31]).
SupposeπX˜ is an orbit of (13), if for any ϵ > 0, there must exist ζ > 0 such that for any point X∈UX˜,ζ, π(X) and πX˜ are ϵ-close, then we call the orbit πX˜ is orbitally stable; if πX˜ is orbitally stable and there exists a ζ0 > 0 such that for all points X∈UX˜,ζ0, we can find a constant δ such that |πt(X) − πt−δ(X)| ⟶ 0, t ⟶ ∞, then we call πX˜ is asymptotically orbitally stable.Lemma 2.
Assume (ξ(t), η(t)) is an order-k periodic solution of (13) with period T. If the Floquet multiplier μ˜ satisfies μ˜<1, where(14)μ˜=∏k=1nΔkexp∫0T∂P∂xξt,ηt+∂Q∂yξt,ηtdt,with(15)Δk=P+∂β/∂y∂φ/∂x−∂β/∂x∂φ/∂y+∂φ/∂x+Q+∂α/∂x∂φ/∂y−∂α/∂y∂φ/∂x+∂φ/∂yP∂φ/∂x+Q∂φ/∂y,and P,Q,∂α/∂x,∂α/∂y,∂β/∂x,∂β/∂y,∂φ/∂x, and ∂φ/∂y are calculated at the point (ξ(tk), η(tk)), P+=Pξtk+,ηtk+,Q+=Qξtk+,ηtk+. Then, (ξ(t), η(t)) is orbitally asymptotically stable.As previously mentioned, if the prey density is less thanxL or greater than x2, there is a high probability that the predator species goes extinct. Therefore, to get a mathematical model which is more nearing to the practice, we assume that the initial pest density is larger than h1 and less than x2.In the following, we will provide a series of qualitative analysis for systems (11) and (12) in the region Ω that is surrounded by x-axis, the threshold line x = h1, and the unstable flow of E2(x2, y2).


### 3.2. Periodic Solution, Homoclinic Bifurcation, and Attraction Region of System (11)

For system (11), we denote the impulsive set by(16)ℳ1=x,yx=h1,y>rh11−h1Ke−βh1.while the corresponding phase set is N1=x,yx=h1+m1,y>0.Theorem 3.
Ifαμ > eβD, x2<K<1/β+x1, rβ−1/K+αμ<2βrD/K, and x1 < h1 < h1 + m1 < x2, then there exists two fixed parameter values 0 < p0 < p∗ < 1 such that(1)
Ifp = p∗, then system (11) admits an order-1 homoclinic cycle;(2)
Ifp0 < p < p∗, then system (11) admits a unique order-1 periodic orbit in Ω;(3)
Ifp∗ < p < 1, then system (11) has no periodic solution in region Ω.Proof.
LetA and D denote the first intersection points of the unstable flow of E2(x2, y2) (in the direction of time increasing) with straight lines x = h1 and x = h1 + m1, respectively. However, point B is the last intersection point of the stable flow of E2(x2, y2) (in the direction of time increasing) and the straight line x = h1 + m1. Besides, lines x = h1 and x = h1 + m1 intersect the isoclinic line x˙ = 0 at points C and G, respectively (see Figure 2). For every point M ∈ Ω, we denote its coordinate by M(xM, yM) and the forward orbit from it by π+(M). Because the phase set N1 is parallel to y-axis, for any point S∈N, we define its coordinate as yS and denote the successor function of point S by fsor(S).
Since 0 <p < 1, there has two fixed parameter values p0 and p∗ such that (1 − p0)yA = yC and (1 − p∗)yA = yB (see the left figure in Figure 2).Figure 2
Existence of order-1 homoclinic cycle and order-1 periodic orbit of system (11).
(a)
(b)(1)
p = p∗. Then, there is fsor(B) = 0, and the orbit BE2^∪E2A^∪AB¯ forms a homoclinic cycle, that is to say, an order-1 circle which has a saddle on it.
(2)
p
0 < p < p∗. To prove that system (11) has an order-1 periodic orbit, we need to find a point S∈N1 such that fsor(S) = 0. Assume that point A jumps to point A+ due to impulsive effect; then, A+ is the successor point of B. Since 1−pyA=yA+>yB=1−p∗yA, we have fsor(B) > 0. Choose point Mh1+m1,yM∈N1 close to point B such that yM > yB; then, the forward orbit π+(M) intersects x = h1 at a point M′ with yA > yM′ and M′ is mapped to point M1 after impulsive effect. Clearly, yM1>yM and fsor(M) > 0. Furthermore, consider the orbit π+(M1), we can easily get fsorM1=yM2−yM1<0. Since fsor is continuous, there must exist a point S∈MM1¯⊂N1 such that fsor(S) = 0 (see the right figure in Figure 2).
In the following, we discuss the uniqueness of the order-1 periodic orbit.
For any pointPxP,yP∈N1,yP<yB, the forward orbit π+(P) tends to the stable node E0(K, 0). Besides, for every point P∈CD¯⊂N1, we have fsor(P) < 0. Therefore, if the order-1 periodic orbit exists, it must begin from BC¯. We arbitrarily choose two points B1,B2∈BC¯,yB<yB1<yB2<yC. Since different trajectories are disjoint, it is easy to obtain fsor(B1) > fsor(B2). That is to say, fsor is a monotonically decreasing function defined on the line segment BC¯, and there is a unique point S∈BC¯⊂N1 such that fsor(S) = 0. From the above, system (11) has a unique order-1 periodic orbit in Ω (see Figure 3).
(3)
p
∗ < p < 1. There is (1 − p)yA < yB, so all the points P on the line segment BD¯⊂N1 satisfy fsor(P) < 0. Besides, the forward orbit π+(P) tends to the stable node E0(K, 0) for any point PxP,yP∈N1,yP<yB. Thus, system (11) has no periodic orbit in Ω. This completes the proof.Figure 3
Uniqueness of order-1 periodic orbit of system (11) when p0 < p < p∗.Remark 1.
Choosep as the bifurcation parameter. According to Theorem 3, we can get that if parameter p equals to the critical value p∗, there is an order-1 homoclinic cycle of system (11); when parameter p gradually decreases from p = p∗to p0 < p < p∗, the homoclinic cycle disappears while an order-1 periodic orbit appears; when parameter p gradually increases from p = p∗ to p∗ < p < 1, the homoclinic cycle disappears but no new order-1 periodic cycle appears.Remark 2.
By Theorem3, system (11) has an order-1 periodic orbit ΓS, S∈BC¯⊂N1. Let point L0+∈U0S,ζ∩N1 and ζ > 0 is a constant that is small enough. Denote the successor point of Ln+ by Ln+1+, n = 0, 1, 2, …. Then, there are two point sequences L2k+ and L2k+1+, k = 0, 1, 2, …, in the phase set N1 which are separated by point S. Obviously, for any n = 0, 1, 2, …, we have yLn+1+=yLn++fsorLn+.
LetdS,Ln+=yS−yLn+ and denote κ=dS,L2+/dS,L0+. According to Definition 3, to prove the orbitally asymptotically stability of the order-1 periodic orbit, we need to show that there exists at least one ζ˜>0 such that κ < 1 for every L0+∈U0ζ˜,S∩N1. If there is always a point L0+∈U0ζ,S∩N1 satisfying κ > 1 for any ζ > 0, then we can get the unstability of the order-1 periodic orbit. If there exists a L0+∈U0ζ,S∩N1 such that κ = 1, then system (11) admits an order-2 periodic orbit passing through L0+. From the analysis in Theorem 3, we can get that function fsor has some good properties, for example, fsor(L) < 0 if L∈N1,yL>yS; fsor(L) > 0 if L∈N1,yL<yS; and fsor(L) = 0 if and only if L = S. However, we cannot ensure that the order-1 periodic orbit is stable.Theorem 4.
For the casep0 < p < p∗ and (1 − p)yG > yB, if αμ > eβD, x2<K<1/β+x1, rβ−1/K+αμ<2βrD/K, and x1 < h1 < h1 + m1 < x2, then system (11) has an attraction region in Ω˜, where Ω˜ is surrounded by the two separatrices of E2(x2, y2) and the threshold line x = h1.Proof.
Let the order-1 periodic orbit ΓS intersect the impulsive set x = h1 at point S′. Denote point B by point L0+. There is fsorL0+>0 and the successor point L1+ is above point S (see Figure 4).
Forπ+L1+, we need to determine the position of its successor point L2+. For this purpose, consider the forward orbit π+L1+ and it meets the impulsive set x = h1 at point L1′. Obviously, L1′ is between two points S′ and G. According to the conditions p0 < p < p∗ and (1 − p)yG > yB, we get 1−pyG<1−pyL1′, so the successor point of L1+ (denoted by L2+) is above the point B and we have L0+<L2+<yS.
Repeating the above steps, two point sequencesL2k+ and L2k+1+, k = 0, 1, 2, …, separated by point S are obtained, where sequence yL2k+k=0,1,2,… is monotonically increasing with yL2k+≤ySwhile yL2k+1+k=0,1,2,… is monotonically decreasing with yS≤yL2k+1+. We further suppose that(17)yL2k⟶yS0≤yS,ask⟶∞,yL2k+1⟶yS1≥yS,ask⟶∞.
And two pointsS0 and point S1 may coincide with each other, that is to say, points S0, S1, and S are the same one.
Let the forward orbitsπ+(S0), π+(S1) meet the impulsive set x = h1 at points S1′ and S1′, respectively. Denote the region circled by the closed curve S0S0′^∪S0′S1′¯∪S1′S1^∪S1S0¯ by △Ω, and we will prove that it is the attraction region of system (11) in Ω˜.
For any pointP0+∈BS0¯, we can pick an integer k > 0 such that yL2k+<yP0<yL2k+1+. Then, similar to the above discussion, we can get a series of successor points P2i+ and P2i+1+, i = 0, 1, 2, …, which satisfy that yP2i+i=0,1,2,… is monotonically increasing with yP2i+≤yS0 and yP2i+1+i=0,1,2,… is monotonically decreasing with yP2i+1+≥yS1. Besides,(18)yP2i+1+⟶yS1,i⟶∞,yP2i+⟶gS0,i⟶∞.
Then, the orbit fromP0 eventually tends to the region △Ω.
Furthermore, for any pointP0+∈S1D¯, the forward orbit π+P0+ will intersect the phase set x = h1 + m1 on BS0¯ after impulsive jumps once or twice and then is eventually attracted to the region △Ω. To sum up, we get that △Ω is the attraction region of system (11) in Ω˜. This completes the proof.Figure 4
Existence of attraction region of system (11).Remark 3.
If three pointsS0, S1, and S coincide with each other, then the attraction region △Ω is exactly the order-1 cycle ΓS we obtained in Theorem 3. In this case, ΓS is globally orbitally asymptotically stable in Ω˜. Of course, these three points S0, S1, and S may not coincide with each other; then, the attraction region △Ω can be circled by an order-2 orbit.


### 3.3. Analysis of System (12)

For system (12), two impulsive sets can be defined as(19)ℳ21=x,yx=h1,y>rh11−h1Ke−βh1,ℳ22=x,yx=h2,y<rh21−h2Ke−βh2.However, the corresponding phase sets are the same which is denoted byN2=x,y∣x=h1+m1=h2−m2≜n,y>0.For the convenience of the following discussion, let the unstable flow ofE2(x2, y2) meet the straight lines x = h1, x = n, and x = h2 for the first time at points A, D, and E, while the stable flow of E2(x2, y2) meet the three lines for the last time at points H, B, and K (in the direction of time increasing). Besides, the three lines intersect the isoclinic line x˙ = 0 at points G, C, and F. For the forward orbit π+(F) and backward orbit π−(F), we denote their intersection points with x = h1 and x = n by N and Q, respectively (see Figure 5).Figure 5
Illustration of the dynamics of system (12).Theorem 5.
IfyQ < (1 − p)yN < yC, αμ > eβD, x2<K<1/β+x1, rβ−1/K+αμ<2βrD/K, and x1 < h1 < n < h2 < x2, then system (12) admits an order-1 periodic orbit in the region Ω.Proof.
Consider the pointL∈N2. If yL < yQ, then the forward orbit π+(L) will firstly meet the second impulsive set ℳ22 and goes through the corresponding impulsive jump, while if yL > yQ, then the forward orbit π+(L) will meet the first impulsive set ℳ21 first and goes through the corresponding impulsive jump.
Similar to the analysis in Theorem3, we can get the existence of the order-1 periodic solution ΓS for system (12), where S∈N2 satisfies yQ < yS < yC. The proof is completed.Remark 4.
Similar to the discussion in Theorem4, we can obtain that if (1 − p)yG > yQ, then system (12) has a local attraction region which is circled by the closed curve AG¯∪GC^∪CQ¯∪QF^∪FE¯∪EDA^.Theorem 6.
Ifr1−n/K/αe−βn<τ<yD−yF, (1 − p)yA < yQ, and x1 < h1 < n < h2 < x2, then system (12) admits an order-2 periodic orbit in Ω provided that conditions in Theorems 1 and 2 are satisfied.Proof.
Sinceτ>r1−n/K/αe−βn, we have yD − yF > τ > yC. Choose a point A1∈N2 next to Q satisfying yA1<yQ, and let the forward orbit π+(A1) meet the impulsive set ℳ22 at point A1′ with yA1′<yF; then, A1′ is mapped to B1∈N2. There is yC<yB1<yD for yD − yF > τ > yC. The forward orbit π+(B1) meets the impulsive set ℳ21 at point C1 with yN<yC1<yA; then, C1 is mapped to D1∈N2. Obviously, we have yD1<yA1. Then, fsor2A1=yD1−yA1<0 (see the left figure in Figure 6).
Choose another pointA2∈N2 satisfying 0<yA2≪1, and let the forward orbit π+(A2) meet the impulsive set ℳ22 at point A2′ with yA2′<yK; then, A2′ is mapped to B2∈N2. There is yC<yB2<yB1 for yD − yF > τ > yC. The forward orbit π+(B2) meets the impulsive set ℳ21 at point C2 with yG<yC2<yC1; then, C2 is mapped to D2∈N2. Obviously, yD2>yA2. Then, fsor2A2=yD2−yA2>0.
Due to the continuity offsor2, there must exist a point S+∈A1A2¯ such that fsor2S=0. Then, system (12) admits an order-2 periodic orbit. The proof is completed.Figure 6
Existence of order-2 periodic orbit of system (12).
(a)
(b)Remark 5.
According to the proof of Theorem5, system (12) may have order-k(k > 2) periodic orbits provided that the parameters satisfy certain conditions.
In the following, we discuss the stability of the order-2 periodic orbitΓ˜S for system (12). Assume Γ˜S passes through four points S+(n, (1 − p)η0), S(h1, η0), W(h2, ζ0), and W+(n, ζ0 + τ) (see the right figure in Figure 6) and its period is T.
Since the analytical expression of this order-2 periodic solution cannot be obtained, we study its stability by using the criterion in Lemma2. According to Lemma 2, we have(20)Px,y=rx1−xK−αxye−βx,Qx,y=μαxye−βx−Dy,α1x,y=m1,β1x,y=−py,φ1x,y=x−h1,α2x,y=−m2,β2x,y=τ,φ2x,y=x−h2.
Then, through simple calculation, we have(21)∂P∂x=r−2rKx−αy1−βxe−βx,∂Q∂y=αμxe−βx−D,∂α1∂x=∂α1∂y=0,∂β1∂x=0,∂β1∂y=−p,∂φ1∂x=1,∂φ1∂y=0,∂α2∂x=∂α2∂y=0,∂β2∂x=∂β2∂y=0,∂φ2∂x=1,∂φ2∂y=0,and then(22)Δ1=P+∂β1/∂y∂φ1/∂x−∂β1/∂x∂φ1/∂y+∂φ1/∂x+Q+∂α1/∂x∂φ1/∂y−∂α1/∂y∂φ1/∂x+∂φ1/∂yP∂φ1/∂x+Q∂φ1/∂y=1−pnr1−n/K−α1−pη0e−βnh1r1−h1/K−αη0e−βh1,Δ2=P+∂β2/∂y∂φ2/∂x−∂β2/∂x∂φ2/∂y+∂φ2/∂x+Q+∂α2/∂x∂φ2/∂y−∂α2/∂y∂φ2/∂x+∂φ2/∂yP∂φ2/∂x+Q∂φ2/∂y=nr1−n/K−αζ0+τe−βnh2r1−h2/K−αζ0e−βh2.
Settingℌt=∂P/∂xξt,ηt+∂Q/∂yξt,ηt, we can get(23)μ˜=Δ1Δ2∫0Tℌtdt=1−pn2r1−n/K−α1−pη0e−βnr1−n/K−αζ0+τe−βnh1h2r1−h1/K−αη0e−βh1r1−h2/K−αζ0e−βh2∫0Tℌtdt.
Then, we have the following result:Theorem 7.
Ifr1−n/K/αe−βn<τ<yD−yF, (1 − p)yA < yQ and x1 < h1 < n < h2 < x2, then the order-2 periodic solution of system (12) is stable provided that μ˜<1.


## 3.1. Preliminary

Before the study of the above two systems, some notations and definitions of semidynamical systems are introduced. These will be used in the following discussion.Consider the general planar semidynamical system:(13)dxdt=Px,y,dydt=Qx,y,x,y∈ℳx,y,△x=αx,y,△y=βx,y,x,y∈ℳx,y,where (x, y) ∈ R2, P(x, y), Q(x, y), α(x, y), and β(x, y) ∈ C[R2, R]. ℳ⊂R2 is called the impulsive set, while φ: R2 ⟶ R2 is called the impulse mapping, where φZ=Z+=xZ+,yZ+=xZ+αxZ,yZ,yZ+βxZ,yZ for ZxZ,yZ∈ℳ. Z+ is the phase point of point Z, and N=φℳ is the phase set.AssumeX(t) = (x(t), y(t)) be an arbitrary solution of (5). Denote πt(X0, X) = X(t) to be the solution corresponding to initial value condition X(0) = X0. Then, its forward orbit can be described by π(X0, X) = {πt(X0, X) ∣ t ≥ 0}, which is also denoted by π(X) or π+(X) for simplicity. We denote Xk=Xtk+∈πX, where X(tk) ∈ π(X), k = 1, 2, ….Definition 1 (see [29, 30]).
An orbitπ(X) of (13) is said to be a periodic orbit if there is a k ∈ N+ such that Xk = X0. Denote k0 = min{k ∈ N+, Xk = X0}, then the orbit π(X) is called an order-k0 periodic orbit.Definition 2 (see [29, 30]).
Assumeπ(X) and πX˜ are the two orbits of (13), if there is a monotone smooth function ν(t), t > 0 such that πtX−πνtX^ϵ,t0, then we called the two orbits as ϵ-close.Definition 3 (see [29, 30, 31]).
SupposeπX˜ is an orbit of (13), if for any ϵ > 0, there must exist ζ > 0 such that for any point X∈UX˜,ζ, π(X) and πX˜ are ϵ-close, then we call the orbit πX˜ is orbitally stable; if πX˜ is orbitally stable and there exists a ζ0 > 0 such that for all points X∈UX˜,ζ0, we can find a constant δ such that |πt(X) − πt−δ(X)| ⟶ 0, t ⟶ ∞, then we call πX˜ is asymptotically orbitally stable.Lemma 2.
Assume (ξ(t), η(t)) is an order-k periodic solution of (13) with period T. If the Floquet multiplier μ˜ satisfies μ˜<1, where(14)μ˜=∏k=1nΔkexp∫0T∂P∂xξt,ηt+∂Q∂yξt,ηtdt,with(15)Δk=P+∂β/∂y∂φ/∂x−∂β/∂x∂φ/∂y+∂φ/∂x+Q+∂α/∂x∂φ/∂y−∂α/∂y∂φ/∂x+∂φ/∂yP∂φ/∂x+Q∂φ/∂y,and P,Q,∂α/∂x,∂α/∂y,∂β/∂x,∂β/∂y,∂φ/∂x, and ∂φ/∂y are calculated at the point (ξ(tk), η(tk)), P+=Pξtk+,ηtk+,Q+=Qξtk+,ηtk+. Then, (ξ(t), η(t)) is orbitally asymptotically stable.As previously mentioned, if the prey density is less thanxL or greater than x2, there is a high probability that the predator species goes extinct. Therefore, to get a mathematical model which is more nearing to the practice, we assume that the initial pest density is larger than h1 and less than x2.In the following, we will provide a series of qualitative analysis for systems (11) and (12) in the region Ω that is surrounded by x-axis, the threshold line x = h1, and the unstable flow of E2(x2, y2).


## 3.2. Periodic Solution, Homoclinic Bifurcation, and Attraction Region of System (11)

For system (11), we denote the impulsive set by(16)ℳ1=x,yx=h1,y>rh11−h1Ke−βh1.while the corresponding phase set is N1=x,yx=h1+m1,y>0.Theorem 3.
Ifαμ > eβD, x2<K<1/β+x1, rβ−1/K+αμ<2βrD/K, and x1 < h1 < h1 + m1 < x2, then there exists two fixed parameter values 0 < p0 < p∗ < 1 such that(1)
Ifp = p∗, then system (11) admits an order-1 homoclinic cycle;(2)
Ifp0 < p < p∗, then system (11) admits a unique order-1 periodic orbit in Ω;(3)
Ifp∗ < p < 1, then system (11) has no periodic solution in region Ω.Proof.
LetA and D denote the first intersection points of the unstable flow of E2(x2, y2) (in the direction of time increasing) with straight lines x = h1 and x = h1 + m1, respectively. However, point B is the last intersection point of the stable flow of E2(x2, y2) (in the direction of time increasing) and the straight line x = h1 + m1. Besides, lines x = h1 and x = h1 + m1 intersect the isoclinic line x˙ = 0 at points C and G, respectively (see Figure 2). For every point M ∈ Ω, we denote its coordinate by M(xM, yM) and the forward orbit from it by π+(M). Because the phase set N1 is parallel to y-axis, for any point S∈N, we define its coordinate as yS and denote the successor function of point S by fsor(S).
Since 0 <p < 1, there has two fixed parameter values p0 and p∗ such that (1 − p0)yA = yC and (1 − p∗)yA = yB (see the left figure in Figure 2).Figure 2
Existence of order-1 homoclinic cycle and order-1 periodic orbit of system (11).
(a)
(b)(1)
p = p∗. Then, there is fsor(B) = 0, and the orbit BE2^∪E2A^∪AB¯ forms a homoclinic cycle, that is to say, an order-1 circle which has a saddle on it.
(2)
p
0 < p < p∗. To prove that system (11) has an order-1 periodic orbit, we need to find a point S∈N1 such that fsor(S) = 0. Assume that point A jumps to point A+ due to impulsive effect; then, A+ is the successor point of B. Since 1−pyA=yA+>yB=1−p∗yA, we have fsor(B) > 0. Choose point Mh1+m1,yM∈N1 close to point B such that yM > yB; then, the forward orbit π+(M) intersects x = h1 at a point M′ with yA > yM′ and M′ is mapped to point M1 after impulsive effect. Clearly, yM1>yM and fsor(M) > 0. Furthermore, consider the orbit π+(M1), we can easily get fsorM1=yM2−yM1<0. Since fsor is continuous, there must exist a point S∈MM1¯⊂N1 such that fsor(S) = 0 (see the right figure in Figure 2).
In the following, we discuss the uniqueness of the order-1 periodic orbit.
For any pointPxP,yP∈N1,yP<yB, the forward orbit π+(P) tends to the stable node E0(K, 0). Besides, for every point P∈CD¯⊂N1, we have fsor(P) < 0. Therefore, if the order-1 periodic orbit exists, it must begin from BC¯. We arbitrarily choose two points B1,B2∈BC¯,yB<yB1<yB2<yC. Since different trajectories are disjoint, it is easy to obtain fsor(B1) > fsor(B2). That is to say, fsor is a monotonically decreasing function defined on the line segment BC¯, and there is a unique point S∈BC¯⊂N1 such that fsor(S) = 0. From the above, system (11) has a unique order-1 periodic orbit in Ω (see Figure 3).
(3)
p
∗ < p < 1. There is (1 − p)yA < yB, so all the points P on the line segment BD¯⊂N1 satisfy fsor(P) < 0. Besides, the forward orbit π+(P) tends to the stable node E0(K, 0) for any point PxP,yP∈N1,yP<yB. Thus, system (11) has no periodic orbit in Ω. This completes the proof.Figure 3
Uniqueness of order-1 periodic orbit of system (11) when p0 < p < p∗.Remark 1.
Choosep as the bifurcation parameter. According to Theorem 3, we can get that if parameter p equals to the critical value p∗, there is an order-1 homoclinic cycle of system (11); when parameter p gradually decreases from p = p∗to p0 < p < p∗, the homoclinic cycle disappears while an order-1 periodic orbit appears; when parameter p gradually increases from p = p∗ to p∗ < p < 1, the homoclinic cycle disappears but no new order-1 periodic cycle appears.Remark 2.
By Theorem3, system (11) has an order-1 periodic orbit ΓS, S∈BC¯⊂N1. Let point L0+∈U0S,ζ∩N1 and ζ > 0 is a constant that is small enough. Denote the successor point of Ln+ by Ln+1+, n = 0, 1, 2, …. Then, there are two point sequences L2k+ and L2k+1+, k = 0, 1, 2, …, in the phase set N1 which are separated by point S. Obviously, for any n = 0, 1, 2, …, we have yLn+1+=yLn++fsorLn+.
LetdS,Ln+=yS−yLn+ and denote κ=dS,L2+/dS,L0+. According to Definition 3, to prove the orbitally asymptotically stability of the order-1 periodic orbit, we need to show that there exists at least one ζ˜>0 such that κ < 1 for every L0+∈U0ζ˜,S∩N1. If there is always a point L0+∈U0ζ,S∩N1 satisfying κ > 1 for any ζ > 0, then we can get the unstability of the order-1 periodic orbit. If there exists a L0+∈U0ζ,S∩N1 such that κ = 1, then system (11) admits an order-2 periodic orbit passing through L0+. From the analysis in Theorem 3, we can get that function fsor has some good properties, for example, fsor(L) < 0 if L∈N1,yL>yS; fsor(L) > 0 if L∈N1,yL<yS; and fsor(L) = 0 if and only if L = S. However, we cannot ensure that the order-1 periodic orbit is stable.Theorem 4.
For the casep0 < p < p∗ and (1 − p)yG > yB, if αμ > eβD, x2<K<1/β+x1, rβ−1/K+αμ<2βrD/K, and x1 < h1 < h1 + m1 < x2, then system (11) has an attraction region in Ω˜, where Ω˜ is surrounded by the two separatrices of E2(x2, y2) and the threshold line x = h1.Proof.
Let the order-1 periodic orbit ΓS intersect the impulsive set x = h1 at point S′. Denote point B by point L0+. There is fsorL0+>0 and the successor point L1+ is above point S (see Figure 4).
Forπ+L1+, we need to determine the position of its successor point L2+. For this purpose, consider the forward orbit π+L1+ and it meets the impulsive set x = h1 at point L1′. Obviously, L1′ is between two points S′ and G. According to the conditions p0 < p < p∗ and (1 − p)yG > yB, we get 1−pyG<1−pyL1′, so the successor point of L1+ (denoted by L2+) is above the point B and we have L0+<L2+<yS.
Repeating the above steps, two point sequencesL2k+ and L2k+1+, k = 0, 1, 2, …, separated by point S are obtained, where sequence yL2k+k=0,1,2,… is monotonically increasing with yL2k+≤ySwhile yL2k+1+k=0,1,2,… is monotonically decreasing with yS≤yL2k+1+. We further suppose that(17)yL2k⟶yS0≤yS,ask⟶∞,yL2k+1⟶yS1≥yS,ask⟶∞.
And two pointsS0 and point S1 may coincide with each other, that is to say, points S0, S1, and S are the same one.
Let the forward orbitsπ+(S0), π+(S1) meet the impulsive set x = h1 at points S1′ and S1′, respectively. Denote the region circled by the closed curve S0S0′^∪S0′S1′¯∪S1′S1^∪S1S0¯ by △Ω, and we will prove that it is the attraction region of system (11) in Ω˜.
For any pointP0+∈BS0¯, we can pick an integer k > 0 such that yL2k+<yP0<yL2k+1+. Then, similar to the above discussion, we can get a series of successor points P2i+ and P2i+1+, i = 0, 1, 2, …, which satisfy that yP2i+i=0,1,2,… is monotonically increasing with yP2i+≤yS0 and yP2i+1+i=0,1,2,… is monotonically decreasing with yP2i+1+≥yS1. Besides,(18)yP2i+1+⟶yS1,i⟶∞,yP2i+⟶gS0,i⟶∞.
Then, the orbit fromP0 eventually tends to the region △Ω.
Furthermore, for any pointP0+∈S1D¯, the forward orbit π+P0+ will intersect the phase set x = h1 + m1 on BS0¯ after impulsive jumps once or twice and then is eventually attracted to the region △Ω. To sum up, we get that △Ω is the attraction region of system (11) in Ω˜. This completes the proof.Figure 4
Existence of attraction region of system (11).Remark 3.
If three pointsS0, S1, and S coincide with each other, then the attraction region △Ω is exactly the order-1 cycle ΓS we obtained in Theorem 3. In this case, ΓS is globally orbitally asymptotically stable in Ω˜. Of course, these three points S0, S1, and S may not coincide with each other; then, the attraction region △Ω can be circled by an order-2 orbit.


## 3.3. Analysis of System (12)

For system (12), two impulsive sets can be defined as(19)ℳ21=x,yx=h1,y>rh11−h1Ke−βh1,ℳ22=x,yx=h2,y<rh21−h2Ke−βh2.However, the corresponding phase sets are the same which is denoted byN2=x,y∣x=h1+m1=h2−m2≜n,y>0.For the convenience of the following discussion, let the unstable flow ofE2(x2, y2) meet the straight lines x = h1, x = n, and x = h2 for the first time at points A, D, and E, while the stable flow of E2(x2, y2) meet the three lines for the last time at points H, B, and K (in the direction of time increasing). Besides, the three lines intersect the isoclinic line x˙ = 0 at points G, C, and F. For the forward orbit π+(F) and backward orbit π−(F), we denote their intersection points with x = h1 and x = n by N and Q, respectively (see Figure 5).Figure 5
Illustration of the dynamics of system (12).Theorem 5.
IfyQ < (1 − p)yN < yC, αμ > eβD, x2<K<1/β+x1, rβ−1/K+αμ<2βrD/K, and x1 < h1 < n < h2 < x2, then system (12) admits an order-1 periodic orbit in the region Ω.Proof.
Consider the pointL∈N2. If yL < yQ, then the forward orbit π+(L) will firstly meet the second impulsive set ℳ22 and goes through the corresponding impulsive jump, while if yL > yQ, then the forward orbit π+(L) will meet the first impulsive set ℳ21 first and goes through the corresponding impulsive jump.
Similar to the analysis in Theorem3, we can get the existence of the order-1 periodic solution ΓS for system (12), where S∈N2 satisfies yQ < yS < yC. The proof is completed.Remark 4.
Similar to the discussion in Theorem4, we can obtain that if (1 − p)yG > yQ, then system (12) has a local attraction region which is circled by the closed curve AG¯∪GC^∪CQ¯∪QF^∪FE¯∪EDA^.Theorem 6.
Ifr1−n/K/αe−βn<τ<yD−yF, (1 − p)yA < yQ, and x1 < h1 < n < h2 < x2, then system (12) admits an order-2 periodic orbit in Ω provided that conditions in Theorems 1 and 2 are satisfied.Proof.
Sinceτ>r1−n/K/αe−βn, we have yD − yF > τ > yC. Choose a point A1∈N2 next to Q satisfying yA1<yQ, and let the forward orbit π+(A1) meet the impulsive set ℳ22 at point A1′ with yA1′<yF; then, A1′ is mapped to B1∈N2. There is yC<yB1<yD for yD − yF > τ > yC. The forward orbit π+(B1) meets the impulsive set ℳ21 at point C1 with yN<yC1<yA; then, C1 is mapped to D1∈N2. Obviously, we have yD1<yA1. Then, fsor2A1=yD1−yA1<0 (see the left figure in Figure 6).
Choose another pointA2∈N2 satisfying 0<yA2≪1, and let the forward orbit π+(A2) meet the impulsive set ℳ22 at point A2′ with yA2′<yK; then, A2′ is mapped to B2∈N2. There is yC<yB2<yB1 for yD − yF > τ > yC. The forward orbit π+(B2) meets the impulsive set ℳ21 at point C2 with yG<yC2<yC1; then, C2 is mapped to D2∈N2. Obviously, yD2>yA2. Then, fsor2A2=yD2−yA2>0.
Due to the continuity offsor2, there must exist a point S+∈A1A2¯ such that fsor2S=0. Then, system (12) admits an order-2 periodic orbit. The proof is completed.Figure 6
Existence of order-2 periodic orbit of system (12).
(a)
(b)Remark 5.
According to the proof of Theorem5, system (12) may have order-k(k > 2) periodic orbits provided that the parameters satisfy certain conditions.
In the following, we discuss the stability of the order-2 periodic orbitΓ˜S for system (12). Assume Γ˜S passes through four points S+(n, (1 − p)η0), S(h1, η0), W(h2, ζ0), and W+(n, ζ0 + τ) (see the right figure in Figure 6) and its period is T.
Since the analytical expression of this order-2 periodic solution cannot be obtained, we study its stability by using the criterion in Lemma2. According to Lemma 2, we have(20)Px,y=rx1−xK−αxye−βx,Qx,y=μαxye−βx−Dy,α1x,y=m1,β1x,y=−py,φ1x,y=x−h1,α2x,y=−m2,β2x,y=τ,φ2x,y=x−h2.
Then, through simple calculation, we have(21)∂P∂x=r−2rKx−αy1−βxe−βx,∂Q∂y=αμxe−βx−D,∂α1∂x=∂α1∂y=0,∂β1∂x=0,∂β1∂y=−p,∂φ1∂x=1,∂φ1∂y=0,∂α2∂x=∂α2∂y=0,∂β2∂x=∂β2∂y=0,∂φ2∂x=1,∂φ2∂y=0,and then(22)Δ1=P+∂β1/∂y∂φ1/∂x−∂β1/∂x∂φ1/∂y+∂φ1/∂x+Q+∂α1/∂x∂φ1/∂y−∂α1/∂y∂φ1/∂x+∂φ1/∂yP∂φ1/∂x+Q∂φ1/∂y=1−pnr1−n/K−α1−pη0e−βnh1r1−h1/K−αη0e−βh1,Δ2=P+∂β2/∂y∂φ2/∂x−∂β2/∂x∂φ2/∂y+∂φ2/∂x+Q+∂α2/∂x∂φ2/∂y−∂α2/∂y∂φ2/∂x+∂φ2/∂yP∂φ2/∂x+Q∂φ2/∂y=nr1−n/K−αζ0+τe−βnh2r1−h2/K−αζ0e−βh2.
Settingℌt=∂P/∂xξt,ηt+∂Q/∂yξt,ηt, we can get(23)μ˜=Δ1Δ2∫0Tℌtdt=1−pn2r1−n/K−α1−pη0e−βnr1−n/K−αζ0+τe−βnh1h2r1−h1/K−αη0e−βh1r1−h2/K−αζ0e−βh2∫0Tℌtdt.
Then, we have the following result:Theorem 7.
Ifr1−n/K/αe−βn<τ<yD−yF, (1 − p)yA < yQ and x1 < h1 < n < h2 < x2, then the order-2 periodic solution of system (12) is stable provided that μ˜<1.


## 4. Numerical Simulations

In this section, a series of numerical simulations are carried out, which will on the one hand confirm the theoretical results and, on the other hand, reveal some peculiar dynamical phenomena.


### 4.1. Existence and Stability of Order-1 and Order-2 Periodic Orbit

To satisfy the conditions listed in Theorems1 and 2, we choose the parameters as follows:(24)r=1.5,K=0.8,α=2,β=1.5,μ=0.615,D=0.3.Then, we can easily getαμ = 1.23 > eβD ≈ 1.223, rβ−1/K+αμ=1.605<2βrD/K≈1.8371, x1 ≈ 0.5995, x2 ≈ 0.7395, and x2<K<1/β+x1. According to Theorems 1 and 2, system (1) admits four equilibria: two boundary equilibria O(0, 0) and E0(K, 0) and two positive steady states E1(x1, y1) and E2(x2, y2), where O(0, 0) and E2(x2, y2) are unstable saddle points while E0(K, 0) and E1(x1, y1) are stable ones.For unilateral control system (11), we choose the lower control threshold h1 = 0.66 > x1 and m1 = 0.05, then there are a fixed parameter values p∗ such that system (11) admits an order-1 periodic orbit for some p < p∗ (see Figure 7), while it has no periodic orbit if p > p∗ (see Figure 8). Obviously, the order-1 periodic orbit in Figure 7 is orbitally stable.Figure 7
The order-1 periodic solution of system (11) with p = 0.3.
(a)
(b)
(c)Figure 8
The predator species tends to be extinct for system (11) with p = 0.8.
(a)
(b)
(c)For bilateral control system (12), we choose the two control thresholds h1 = 0.61 > x1 and h2 = 0.71 < x2, while m1 = 0.05 and m2 = 0.05 such that the shared phase set is x = n = 0.66. By changing the impulse parameter, we show that system (12) can have different types of periodic solution (see Figures 9–12).Figure 9
The order-1 periodic orbit of system (12). Control parameters: p = 0.7 and τ = 0.2.
(a)
(b)
(c)Figure 10
The order-2 periodic orbit of system (12). Control parameters: p = 0.6 and τ = 0.4.
(a)
(b)
(c)Figure 11
The order-3 periodic orbit of system (12). Control parameters: p = 0.6 and τ = 0.5.
(a)
(b)
(c)Figure 12
The order-4 periodic orbit of system (12). Control parameters: p = 0.6 and τ = 0.7.
(a)
(b)
(c)When the order-2 periodic orbit of system (12) exists, we investigate its stability. Choose different initial values; we can see that this periodic solution is orbitally stable (see Figure 13).Figure 13
Orbital stability of the order-2 periodic solution of system (12). Control parameters: p = 0.6 and τ = 0.4.
(a)
(b)
(c)


### 4.2. Bifurcation Diagram of System (12)

According to the discussion in Section4.1, when we fix p = 0.6, system (12) can exhibit different types of periodic solution with the change of parameter τ. We now exploit the bifurcation phenomenon with respect to τ. For this end, we first set p = 0.6, and Figure 14 shows the bifurcation diagram of system (12) with respect to τ (τ ∈ [0, 1]). From Figure 14, we find that system (12) can exhibit several types of periodic solution, such as order-k periodic solution, k = 1, 2, 3, 4, 5. However, the order bifurcation cannot lead to chaos. This is somewhat similar to a discrete dynamical system.Figure 14
Bifurcation diagram of system (12) with respect to τ for τ ∈ [0, 1].Besides, from Figure14, we can see that orders 1, 2, 3, 4, and 5 periodic orbits exist for wide ranges of τ. With the increase of τ, there is an order-1 periodic orbit when τ ∈ [0, 0.2401) and then a positive order-2 periodic solution bifurcates from the order-1 periodic solution at τ ≈ 0.2401 through a fold bifurcation. The order-2 periodic orbit exists when τ ∈ (0.2401, 0.4031); then, a positive order-3 periodic solution bifurcates from the order-2 periodic solution at τ ≈ 0.4031. Analogously, a positive order-4 periodic solution bifurcates from the order-3 periodic solution at τ ≈ 0.6165 and a positive order-5 periodic solution bifurcates from the order-4 periodic solution at τ ≈ 0.8952.In addition, we also setτ = 0.7 and observe the bifurcation phenomenon with respect to p. Figure 15 shows the bifurcation diagram of system (12) with respect to p (p ∈ (0, 1)). Similarly, we find that system (12) also exhibits many types of periodic solution, and the order bifurcation does not lead to chaos. Furthermore, for wide ranges of p, system (12) has periodic solutions of order 4 (for example p ∈ (0.3, 0.7) (see case p = 0.6 in Figure 12) or order 5 (for example, p ∈ (0.18, 0.25] ∪ (0.75, 0.82), see case p = 0.25 in Figure 16(c)).Figure 15
Bifurcation diagram of system (12) with respect to p for p ∈ (0, 1).Figure 16
Change of the order of periodic orbits of system (12) for p ∈ [0, 0.3]: (a) order-8 periodic orbit at p = 0.1; (b) order-6 periodic orbit at p = 0.15; (c) order-5 periodic orbit at p = 0.25; (d) order-4 periodic orbit at p = 0.3.
(a)
(b)
(c)
(d)Besides, from Figure15, we can see that, with the change of τ, there are two small subintervals in which the order of the periodic solution rises and falls dramatically. Figures in Figure 17 are two zoomed parts in the parameter interval p ∈ (0, 0.3] and p ∈ [0.7, 1). From the left one, we can see that the order of the periodic solution is very high when p is near zero and then gradually drops down to 4 (for example, 0.28 < p < 0.3). However, from the right one, we can see that the order of the periodic solution rises from 4 (0.7 < τ < 0.73) to a high level.Figure 17
Bifurcation diagrams of system (12) with respect to p for p ∈ (0, 0.3] and p ∈ [0.7, 1).
(a)
(b)The time series of system (12) are present in Figure 16 for four different values of p ∈ (0, 0.3]. We can see easily that the order of the periodic solution changes greatly. Similarly, the time series of the system (12) for p ∈ [0.9, 0.95] are present in Figure 18 and the order also changes dramatically. Obviously, these profiles of system (12) are associated with Figures 15 and 17.Figure 18
Change of the order of periodic orbits of system (12) for p ∈ [0.9, 0.95]: (a) order-8 periodic orbit at p = 0.9; (b) order-9 periodic orbit at p = 0.92; (c) order-10 periodic orbit at p = 0.93; (d) order-13 periodic orbit at p = 0.95.
(a)
(b)
(c)
(d)


## 4.1. Existence and Stability of Order-1 and Order-2 Periodic Orbit

To satisfy the conditions listed in Theorems1 and 2, we choose the parameters as follows:(24)r=1.5,K=0.8,α=2,β=1.5,μ=0.615,D=0.3.Then, we can easily getαμ = 1.23 > eβD ≈ 1.223, rβ−1/K+αμ=1.605<2βrD/K≈1.8371, x1 ≈ 0.5995, x2 ≈ 0.7395, and x2<K<1/β+x1. According to Theorems 1 and 2, system (1) admits four equilibria: two boundary equilibria O(0, 0) and E0(K, 0) and two positive steady states E1(x1, y1) and E2(x2, y2), where O(0, 0) and E2(x2, y2) are unstable saddle points while E0(K, 0) and E1(x1, y1) are stable ones.For unilateral control system (11), we choose the lower control threshold h1 = 0.66 > x1 and m1 = 0.05, then there are a fixed parameter values p∗ such that system (11) admits an order-1 periodic orbit for some p < p∗ (see Figure 7), while it has no periodic orbit if p > p∗ (see Figure 8). Obviously, the order-1 periodic orbit in Figure 7 is orbitally stable.Figure 7
The order-1 periodic solution of system (11) with p = 0.3.
(a)
(b)
(c)Figure 8
The predator species tends to be extinct for system (11) with p = 0.8.
(a)
(b)
(c)For bilateral control system (12), we choose the two control thresholds h1 = 0.61 > x1 and h2 = 0.71 < x2, while m1 = 0.05 and m2 = 0.05 such that the shared phase set is x = n = 0.66. By changing the impulse parameter, we show that system (12) can have different types of periodic solution (see Figures 9–12).Figure 9
The order-1 periodic orbit of system (12). Control parameters: p = 0.7 and τ = 0.2.
(a)
(b)
(c)Figure 10
The order-2 periodic orbit of system (12). Control parameters: p = 0.6 and τ = 0.4.
(a)
(b)
(c)Figure 11
The order-3 periodic orbit of system (12). Control parameters: p = 0.6 and τ = 0.5.
(a)
(b)
(c)Figure 12
The order-4 periodic orbit of system (12). Control parameters: p = 0.6 and τ = 0.7.
(a)
(b)
(c)When the order-2 periodic orbit of system (12) exists, we investigate its stability. Choose different initial values; we can see that this periodic solution is orbitally stable (see Figure 13).Figure 13
Orbital stability of the order-2 periodic solution of system (12). Control parameters: p = 0.6 and τ = 0.4.
(a)
(b)
(c)


## 4.2. Bifurcation Diagram of System (12)

According to the discussion in Section4.1, when we fix p = 0.6, system (12) can exhibit different types of periodic solution with the change of parameter τ. We now exploit the bifurcation phenomenon with respect to τ. For this end, we first set p = 0.6, and Figure 14 shows the bifurcation diagram of system (12) with respect to τ (τ ∈ [0, 1]). From Figure 14, we find that system (12) can exhibit several types of periodic solution, such as order-k periodic solution, k = 1, 2, 3, 4, 5. However, the order bifurcation cannot lead to chaos. This is somewhat similar to a discrete dynamical system.Figure 14
Bifurcation diagram of system (12) with respect to τ for τ ∈ [0, 1].Besides, from Figure14, we can see that orders 1, 2, 3, 4, and 5 periodic orbits exist for wide ranges of τ. With the increase of τ, there is an order-1 periodic orbit when τ ∈ [0, 0.2401) and then a positive order-2 periodic solution bifurcates from the order-1 periodic solution at τ ≈ 0.2401 through a fold bifurcation. The order-2 periodic orbit exists when τ ∈ (0.2401, 0.4031); then, a positive order-3 periodic solution bifurcates from the order-2 periodic solution at τ ≈ 0.4031. Analogously, a positive order-4 periodic solution bifurcates from the order-3 periodic solution at τ ≈ 0.6165 and a positive order-5 periodic solution bifurcates from the order-4 periodic solution at τ ≈ 0.8952.In addition, we also setτ = 0.7 and observe the bifurcation phenomenon with respect to p. Figure 15 shows the bifurcation diagram of system (12) with respect to p (p ∈ (0, 1)). Similarly, we find that system (12) also exhibits many types of periodic solution, and the order bifurcation does not lead to chaos. Furthermore, for wide ranges of p, system (12) has periodic solutions of order 4 (for example p ∈ (0.3, 0.7) (see case p = 0.6 in Figure 12) or order 5 (for example, p ∈ (0.18, 0.25] ∪ (0.75, 0.82), see case p = 0.25 in Figure 16(c)).Figure 15
Bifurcation diagram of system (12) with respect to p for p ∈ (0, 1).Figure 16
Change of the order of periodic orbits of system (12) for p ∈ [0, 0.3]: (a) order-8 periodic orbit at p = 0.1; (b) order-6 periodic orbit at p = 0.15; (c) order-5 periodic orbit at p = 0.25; (d) order-4 periodic orbit at p = 0.3.
(a)
(b)
(c)
(d)Besides, from Figure15, we can see that, with the change of τ, there are two small subintervals in which the order of the periodic solution rises and falls dramatically. Figures in Figure 17 are two zoomed parts in the parameter interval p ∈ (0, 0.3] and p ∈ [0.7, 1). From the left one, we can see that the order of the periodic solution is very high when p is near zero and then gradually drops down to 4 (for example, 0.28 < p < 0.3). However, from the right one, we can see that the order of the periodic solution rises from 4 (0.7 < τ < 0.73) to a high level.Figure 17
Bifurcation diagrams of system (12) with respect to p for p ∈ (0, 0.3] and p ∈ [0.7, 1).
(a)
(b)The time series of system (12) are present in Figure 16 for four different values of p ∈ (0, 0.3]. We can see easily that the order of the periodic solution changes greatly. Similarly, the time series of the system (12) for p ∈ [0.9, 0.95] are present in Figure 18 and the order also changes dramatically. Obviously, these profiles of system (12) are associated with Figures 15 and 17.Figure 18
Change of the order of periodic orbits of system (12) for p ∈ [0.9, 0.95]: (a) order-8 periodic orbit at p = 0.9; (b) order-9 periodic orbit at p = 0.92; (c) order-10 periodic orbit at p = 0.93; (d) order-13 periodic orbit at p = 0.95.
(a)
(b)
(c)
(d)


## 5. Conclusion

In this paper, mathematical models for the management of two kinds of biological resources based on a given predator-prey relationship are proposed, and two types of control strategies, unilateral and bilateral control with the state feedback of the prey population, are investigated both theoretically and numerically.For the basic predator-prey system (1), we extend and renew the study of [26] and find conditions under which system (1) has two positive equilibria while E1(x1, y1) is a stable saddle and E2(x2, y2) is an unstable saddle (Theorem 1). Besides, we also provide conditions that can guarantee that the system (refer to system (1)) has no closed orbits (refer to Theorem 2).Through analysis of the phase diagram, we find that the basic system has risk of resource depletion without control measures. Then, we introduce two types of control strategies, unilateral and bilateral control, both of which apply impulsive state feedback control technique. By applying the method of successor function, the existence of the order-1 homoclinic orbit and order-1 periodic orbit of the unilateral control system is obtained (refer to Theorem3) and the attraction region of this system is also exploited (refer to Theorem 4). The phenomenon of homoclinic bifurcation is discussed by choosing p as the bifurcation parameter. We first verify the existence of the bifurcation point p∗ and then show that if p = p∗, system (11) admits an order-1 homoclinic cycle; when parameter p gradually decreases from p = p∗ to p < p∗, the homoclinic cycle disappears while an order-1 periodic orbit appears; when parameter p gradually increases from p = p∗ to p∗ < p < 1, the homoclinic cycle disappears but no new order-1 periodic cycle appears. Besides, by using geometry theory and Analogue of Poincaré criterion, we investigate the existence and stability of positive periodic solutions of the bilateral control system (refer to Theorem 5–7). Furthermore, a series of numerical simulations (including bifurcation diagram of periodic solution) are performed which not only confirm the theoretical results but also reveal some peculiar dynamical phenomena, such as the appearance of high-order periodic solutions and existence of parameter intervals with drastic order change of periodic solution.Our study show that both unilateral control and bilateral control are beneficial to the biodiversity of the ecological system we studied in this paper. They can reduce the risk of resource depletion to some extent. From the unilateral control, we can see that the risk of extinction of the predator population still exists. That is, when the density of the predator is low enough while the prey population continues to increase, the predator will go extinct. However, the bilateral control can avoid the extinction of both species. Therefore, when comparing the two management strategies, our study encourages bilateral control rather than unilateral control for the risk of predator extinction.



---

*Source: 9453941-2020-03-10.xml*