Dynamic Complexity in a Prey-Predator Model with State-Dependent Impulsive Control Strategy

In this paper, an ecological model described by a couple of state-dependent impulsive equations is studied analytically and numerically. The theoretical analysis suggests that there exists a semitrivial periodic solution under some conditions and it is globally orbitally asymptotically stable. Furthermore, using the successor function, we study the existence, uniqueness, and stability of order-1 periodic solution, and the boundedness of solution is also presented. The relationship between order- k successor function and order- k periodic solution is discussed as well, thereby giving the existence condition of an order-3 periodic solution. In addition, a series of numerical simulations are carried out, which not only support the theoretical results but also show the complex dynamics in the model further, for example, the coexistence of multiple periodic solutions, chaos, and period-doubling bifurcation.


Introduction
Since mathematical models represented by differential equations were introduced into biological and ecological systems, it has been proved that they are very useful tools to deal with practical problems such as emerging disease [1] and population dynamics between plankton and nutrient [2][3][4], which also promote the development of the theory of differential equation further. Especially, in the process, there is a significant development in impulsive differential equations [5][6][7][8][9] because many reality systems exhibit the abrupt jumps phenomena in population sizes during the evolution processes. Additionally, a theoretical method is urgently needed to provide a guide for the management of some practical problems such as pest outbreaks, phytoplankton blooms, vaccination, and so on. Accordingly, impulsive differential equation theory may prove to be one of the most potential theories to strategize about the approach of management for these problems, as impulsive differential equation can model abrupt jump behaviors induced by management, for example, some management occurs at fixed time, which can be modeled using the models with impulses at fixed times [10][11][12][13].
State-dependent impulsive differential equation, as a branch of impulsive differential equation, has the advantage to model the control behaviors depending on the population state. In recent years, the studies on state-dependent impulsive differential equation have been paid increasing attention [13][14][15][16]. Many studies devote to the properties of periodic solutions [17][18][19][20], including the existence, stability and periodicity, etc. It has been confirmed that the Poincaré map is a very useful tool to prove the existence of periodic solution in state-dependent impulsive differential equations [21][22][23]. In addition, Simeonov and Bainov [24] gave analogue of Poincaré criterion in ref. [25], which contributes to the stability of the periodic solution. Based on the Poincaré map, it can be proved theoretically that the positive periodic solution can bifurcate from the trivial periodic solution through a transcritical bifurcation [25]. Particularly, Chen [26] introduced the successor function a in semicontinuous dynamical system described by state-dependent impulsive differential equation to prove the existence of periodic solution. Using the successor function, Dai et al. [27] studied homoclinic bifurcation in a semicontinuous dynamical system [28]. e successor function has been extensively employed to investigate the existence of periodic solution and even stability and bifurcation in state-dependent impulsive differential equation [28][29][30][31][32][33][34][35][36][37][38]. ese studies enrich the theory of state-dependent impulsive differential equation further.
In this paper, we will present a prey-predator two-species model with state-dependent impulsive control strategy to study the dynamics induced by impulsive control. Firstly, a model describing population dynamics of prey and predator is given below.
where x and y denote the densities of prey and predator, respectively; r is the intrinsic growth rate of the prey; K is the carrying capacity of the prey; b is the maximum growth rate; a represents the half-saturation concentration; α is the assimilation efficiency of predator; and m is the specific mortality rate of predator. Here, we employ a Holling II functional response to present the interaction between prey and predator because some results suggest that the population dynamic following Holling II may exist in reality [2,3,39,40]. For some species, when their densities are beyond a certain value, some seriously negative effects will emerge, for instance, algal bloom. For this reason, their densities should be kept below the critical value by harvest or/and released natural enemy. Obviously, the occurrence of control will rely on the density of the species. Here, under the assumption that the interaction between prey x and predator y follows model (1), we investigate the dynamics of impulsive control and control strategies for prey x using state-dependent impulsive differential equation, and the model can be described as dx dt where h > 0 denotes the critical value below which the density of prey x should be kept. 0 < p < 1 is the harvest rate of prey x; τ ≥ 0 is the releasing amounts of predator y; and 0 < q < 1 represents the harvest rate of predator y or the level of negative effect caused by the harvest of prey x on predator y. e rest of the paper is organized as follows. In Section 2, some notations, basic definitions, and lemmas are given. In Section 3, a semitrivial solution and its stability are studied; we also study the existence, uniqueness, and stability of an order-1 periodic solution; moreover, the boundedness of solution and some propositions of periodic solution are discussed as well. Whereafter, a series of numerical results are carried out to study the dynamics of model (2) further in Section 4. Finally, the paper ends with the conclusion in Section 5.

Preliminaries
In model (1), there exist three equilibria, E 0 � (0, 0), and y * � r(K − x * )(a + x * )/(bK). When the condition max(0, αb(K − a)/(K + a)) < m < Kαb/(a + K) holds, E * is a stable positive focus or node. We assume that the following condition holds throughout the paper: e vertical isocline and the horizontal isocline can be defined as follows: Let z(t) � (x(t), y(t)) be any solution of model (2); then, we define the positive trajectory initializing at point In this paper, some notions are always mentioned, including phased point, impulsive point, and successor point. From works in [20,23,26], we state the following definition. Definition 1. Take the point A 0 ∈ N; set π(t; A 0 , t 0 ) as the trajectory initializing at point A 0 in model (2). If there exists t k > t 0 such that the trajectory π(t; A 0 , t 0 ) intersects the impulsive set M at A k � (h, y k ) when t � t k (k  | and A + k is the order-k successor point of point A 0 . e map F k is called the order-k successor function, that is, Definition 3. Take a point S in the phase set N; then, π(t; S, t 0 ) is a solution of model (2). If there exists a positive integer k ≥ 1 such that k is the smallest integer for F k (S) � 0, then the solution π(t; S, t 0 ) is called order-k periodic solution.
Lemma 1 (see [26]). e successor function F k (S) is continuous, where the point S ∈ N.
Obviously, we can obtain the following lemma by Lemma 1.
For the point O, we always assume the point

Main Results
When the values of parameters in model (1) are given, the dynamics of model (2) only depend on the parameters h, p, q, and τ. In this paper, we are interested in the effect of impulsive control on dynamics of model (2). Hence, we will investigate the dynamical behaviors of model (2) under the condition that the parameters of model (1) are fixed.
For points O and H, it is obvious that In this section, we just show the proof with respect to case (K − a)/2 < (1 − p)h < (1 − p/2)h in the figure of illustration. e proof corresponding to other cases is similar, which has been omitted.
Proof. For τ � 0, solutions of model (2) will enter the space Π eventually. Hence, we study the stability of the semitrivial periodic solution in this space.
. en, we assume that the trajectory π(t; N ε , t 0 ) intersects with impulsive set M at point P � (h, y P ), and y P < y N ε � ε under condition (H). Set δ � ε; then, for any point G � Figure 1(a)), which suggests that the distance holds for all t > t 0 , where Γ � (ξ(t), η(t)). Hence, according to the definition of orbital stability in ref. [40], the semitrivial periodic solution Γ is orbitally stable. For any G 0 ∈ N, take a sequence of number en, the sequence of number y k is a strictly decreasing positive sequence because h < x * and q ∈ (0, 1). Hence, there exists a y * ≥ 0, such that en, for any point Q ∈Q * Q # (i.e., y Q * < y Q < y Q # ), its order-1 successor point is below the point Q * (see Figure 1(b)), that is, y Q + ≤ y * .
Because lim k⟶+∞ y k � y * , there exists a positive integer Θ 1 such that y Q * < y k < y Q # for k > Θ 1 . According to above analysis, there exists a positive integer Θ 2 > Θ 1 such that y k < y * when k > Θ 2 , which contradicts with y k ≥ y * . Hence, lim k⟶+∞ y k � y * � 0, which means that the semitrivial periodic solution Γ is orbitally attractive. So, the semitrivial periodic solution Γ is orbitally asymptotically stable. Since the attraction domain is Ω, the semitrivial periodic solution Γ is globally orbitally asymptotically stable.
is completes the proof. Suppose that M is a positive constant, and then we have the following theorem for any M ∈ (0, +∞).
Proof. For any τ > 0, Figure 2(a) displays the relation between point O and its order- According to Lemmas 1 and 2, there exists an order-1 periodic solution in model (2). If in phase set N. en, let the point R + 1 be the order-1 successor point of point R and point R 1 be the impulsive point of point R + 1 , and there exists an order-1 periodic solution in model (2) using Lemmas 1 and 2.
Suppose there are two order-1 periodic solutions in model (2), which are Γ 1 and Γ 2 , respectively. Let point A + � ((1 − p)h, y A + ) and A � (h, y A ) be the phase point and impulsive point of Γ 1 , respectively, and point B + � ((1 − p)h, y B + ) and B � (h, y B ) be the phase point and impulsive point of Γ 2 , respectively. Without loss of generality, set y A + > y B + . en, the relation among A + , B + , and O is (i) Figure 2(b) shows case (i); it is obvious that Γ 1 and Γ 2 cannot coexist by disjoint property of trajectory of model (1). Hence, the order-1 periodic solution is unique for this case.
For case (ii), because A + ∈B + O, we only need to focus on the region Π. We rewrite the model (1) in the region Π as the following initial-value problem: Obviously, the function Ψ(x, y) is continuously differentiable and monotone decreasing in y on the space Π. en, the solution of model (6) is equivalent to the following integral equation: By equation (7), set  (2).
On the other hand, because both Γ 1 and Γ 2 are order-1 periodic solutions, we have (2), which is a contradiction with inequality (9). Hence, the order-1 periodic solution is unique for this case.
For case (iii), the proof is similar with case (ii), so it is omitted.
is completes the proof. From eorem 2, we can find that the existence and stability of order-1 periodic solution depend on parameter τ.
Proof. When τ ∈ (0, τ c ), the order-1 periodic solution of model (2) must belong to space Π under condition (H). Γ + denotes an order-1 periodic solution of model (2), whose phased point and impulsive point are point C � (x C , y C ) and point D � (x D , y D ), respectively. In addition, the solutions of model (2) will enter the space Π eventually when τ ∈ (0, τ c ). Hence, we only need to study the stability of order-1 periodic solution on the space Π.
For case (1) (see Figure 3(b)), we have lim k⟶+∞ y k � y * ≥ y C and y k ≥ y * . Let Q * � ((1 − p)h, y * ); there exists a point Q # � ((1 − p)h, y # ) when y * > y C , where y # > y * , such that point Q * is the order-1 successor point of point Q # . For any point Q ∈Q * Q # (i.e., y Q * < y Q < y Q # ), its order-1 successor point is below the point Q * . Let point Q + be the order-1 successor point of Q, then y Q + < y * .
Because lim k⟶+∞ y k � y * and y k ≥ y * , there exists a positive integer Θ 1 such that y Q * < y k < y Q # for k > Θ 1 . According to above analysis, there exists a positive integer Θ 2 > Θ 1 such that y k < y * when k > Θ 2 , which contradicts with y k ≥ y * . So, lim k⟶+∞ y k � y C . Likewise, lim k⟶+∞ y k � y C also holds for case (2), which says the order-1 periodic solution Γ + is orbitally attractive.
Hence, the order-1 periodic solution Γ + is orbitally asymptotically stable. Since the attraction domain is Ω, the order-1 periodic solution Γ + is globally orbitally asymptotically stable.
is completes the proof.

Theorem 4. Under condition (H), if h < x * and τ ≥ y O , then there exists a globally orbitally asymptotically stable and unique order-1 periodic solution or a unique order-1 periodic solution and an order-2 solution coexists in model
Proof. By Proposition 2, it is obvious that there is no orderk (k ≥ 3) periodic solution in model (2) τ, q, p, and h. is completes the proof. From eorem 5, it is not difficult to find that there is a rough estimate of the bound, and we have the following remark.

Remark 1. Under condition (H), if
h < x * , then 0 < lim t⟶+∞ sup x(t) ≤ h and 0 ≤ lim t⟶+∞ sup y(t) < (1 − q)y H + τ. (2)  Proof. Suppose positive integer k f is not a factor of the positive integer k, then let k r be the remainder of k/k f and k q be the quotient of k/k f . Because there exists an order-k f periodic solution in model (2), without loss of generality, assume F k f (S) � 0. Obviously, we have F k q k f (S) � 0, but F k r (S) ≠ 0. Hence, the equation F k q k f +k r (S) ≠ 0, which is a contradiction with F k (S) � 0 because k � k q k f + k r . So, the positive integer k f is a factor of the positive integer k. When the positive integer k f is a factor of the positive integer k, F k (S) � 0 is obvious if there exists an order-k f periodic solution in model (2).

Proposition 3. Under condition (H), if there exists a point S such that F k (S) � 0 in model
is completes the proof.

Proposition 4. Under condition (H), if there exists a point
, and q ∈ (0, 1), then model (2) has an order-3 periodic solution when F 1 (O) > 0, and the order-3 periodic solution coexists with an order-1 periodic solution. is completes the proof. Although it is not proved that point O 0 must exist in Proposition 4, it provides a method to search an order-3 periodic solution for numerical analysis. Furthermore, the existence of order-3 periodic solution means that orderk (k > 3) periodic solution and chaos exist in model (2) by Sarkovskii's theorem [42] and Li and Yorke's theorem [43].

Numerical Results
In the previous section, we have analyzed the properties of periodic solution in model (2), including the existence, uniqueness, stability, and boundedness. In this section, we will further show the complex dynamics induced by impulsive control using the numerical simulations. Since the focus of this paper is the role of impulsive control in the population dynamics, the parameters of model (1) are fixed: r � 0.9, K � 100, b � 1.2, a � 7, α � 0.7, and m � 0.7, and we can obtain the equilibrium E * � (35, 20.475) with index +1.
Obviously, condition (H) always holds under this parameter set.
By eorem 1, there only exists a semitrivial periodic solution in model (2) when τ � 0 and h < x * . Given h � 3, p � 0.6, q � 0.5, and τ � 0, we get a semitrivial periodic solution as shown in Figure 4(a), and its period is about 1.309. In addition, we also obtain the theoretical value of the period by T � (1/r)ln(K − (1 − p)h/ (1 − p)(K − h)), which is the same as the numerical one. When h � 3 and p � 0.6, the point O is fixed at (1.2, 6.0762). e numerical simulation shows that all the successor points of any trajectory except the semitrivial solution are below initial point in model (2) and these trajectories converge towards the semitrivial periodic solution (see Figure 4(b)).
When τ > 0, however, the semitrivial periodic solution disappears. en, an order-1 periodic solution emerges and it should be unique by eorem 2. In phase set N, because the order-1 successor point of any point above point O is below point O + 1 , we compute the order-1 successor function for different values of τ (see Figure 5(a)), where h � 32, p � 0.1, and q � 0.1. In the rest of this section, we will employ this control parameter set unless otherwise specified.
From Figure 5(a), it is not difficult to find that the equation F 1 (y) � 0 has a unique root, which means that there exists a unique order-1 periodic solution. Hence, we carry out some numerical simulations with respect to these roots of F 1 (y) � 0, and the results verify that these solutions are order-1 periodic solutions as shown in Figure 5(b). But we find that some of these order-1 periodic solutions are unstable, which indicates that there may exist order-k periodic solution for some values of τ.
In order to illustrate the existence of order-k periodic solution, we further compute order-1 successor function, order-2 successor function, and order-3 successor function, and the results are very interesting and are shown in Figure 6(a). From Figure 6(a), it is easy to find that the equation F 1 (y) � 0 has a unique root, the equation F 2 (y) � 0 has three roots, and the equation F 3 (y) � 0 has seven roots. Furthermore, we can find that the root of F 1 (y) � 0 is the common root of three equations: F 1 (y) � 0, F 2 (y) � 0, and F 3 (y) � 0, which means that there exists an order-1 periodic solution corresponding to this root. Except this root, the equation F 2 (y) � 0 still has two roots which indicate an order-2 periodic solution exists in virtue of the uniqueness of order-1 periodic solution. For order-3 successor function F 3 (y), the equation F 3 (y) � 0 has six roots except for the root of F 1 (y) � 0.
Except for the root of F 1 (y) � 0, it is not difficult to find that F 1 (y) and F 2 (y) are positive or negative definite in a small neighborhood of each root of F 3 (y) � 0, which means that there may exist two order-3 periodic solutions in model (2). us, we perform some numerical simulations to show the solutions corresponding to these roots. e results show that a stable order-3 periodic solution coexists with an unstable order-3 periodic solution, an unstable order-2 periodic solution, and an unstable order-1 periodic solution (see Figure 6(b)). Additionally, the numerical solutions for Complexity 7 prey y are given to further exhibit these periodic solutions, as shown in Figure 6(c). e above results indicate that the order-1 periodic solutions are not always stable. Using order-1 and order-2 successor functions, we compute y − value with respect to the phased point of order-1 periodic solution to simulate its variation as parameter τ increases (see Figure 7(a)). From Figure 7(a), the order-1 periodic solution is stable before parameter τ reaches line l 1 , and the order-1 periodic solution is still stable near l 1 when parameter τ is beyond l 1 . is result is in accord with eorem 3. Keeping parameter τ increasing, the order-1 periodic solution loses its stability until it is beyond l 3 . It is well worth noting that model (2) only has order-1 and order-2 periodic solution when parameter τ is between l 2 and l 3 , which agrees with eorem 4. Furthermore, Figure 7(a) demonstrates that the order-1 periodic solution is unique, which supports the result in eorem 2. Additionally, Figure 7(a) also shows that y � y(τ) is a monotonously increasing function with respect to τ, which is consistent with Proposition 1.
In addition, our numerical results also indicate that the phased point of order-1 periodic solution rises along phase set as parameter τ increases. In order to present the stability of order-1 periodic solution, we compute two groups of order-k successor points whose initial point is on the sides of the phased point of order-1 periodic solution, respectively. It is revealed that these points monotonously converge towards the phased point of order-1 periodic solution (see Figure 7(b)). Hence, this order-1 periodic solution is stable due to the nonexistence of other periodic solutions.
In order to display complex dynamics, we chose parameter τ as the controlled parameter and ran numerous simulations using a wide range of parameter τ to show the bifurcation diagram (see Figure 8(a)). Since the stability of order-1 periodic solution has been discussed above (Figure 7(a)), we will present the dynamics when the order-1 periodic solution is unstable. From Figure 8    chaos, which suggests that all states described by these timeseries are in disorder.

Conclusion
In this paper, we proposed an ecological model with statedependent impulsive control strategy to investigate the effect of impulsive control on prey-predator dynamics. eoretical results suggest that model (2) has a globally orbitally asymptotically stable semitrivial periodic solution when parameter τ � 0 and h < x * , where condition (H) holds to guarantee existence of positive equilibrium E * . When τ > 0, the semitrivial periodic solution disappears, but the theoretical analysis shows that a unique order-1 periodic solution emerges for any q ∈ (0, 1) and p ∈ (0, 1). Furthermore, we find that the order-1 periodic solution is globally orbitally asymptotically stable when the value of τ is below τ c . Additionally, when the value of τ is beyond y o , theoretical analysis reveals that a globally orbitally asymptotically stable order-1 periodic solution exists or coexistence of a unique order-1 periodic solution and order-2 periodic solution occurs.
Additionally, the existence of order-k is theoretically analyzed based on equation F k (y) � 0, which demonstrates that an order-k f periodic solution exists when F k (y) � 0, and k f is a factor of the positive integer k, including 1 and k.
According to this result, we discuss the existence of order-3 periodic solution. Although we cannot get the exact range of parameter τ with respect to the existence of order-3 periodic solution, our results imply that the order-3 periodic solution may exist for τ ∈ (τ c , y O ), which even means that chaos may occur.
In order to study dynamics induced by impulsive control further, a series of numerical simulations are carried out, which agreed with the theoretical results well. Numerical simulations show that order-1, order-2, and order-3 periodic solutions can coexist. Specially, we find the coexistence of double order-3 periodic solutions. Moreover, numerical bifurcation analysis shows that an unstable order-1 periodic solution and a stable order-2 periodic solution can coexist for τ > y o , but there exists only a stable order-1 periodic solution when τ is beyond some critical value. In addition, numerical simulations show that chaos occurs. Obviously, dynamics induced by impulsive control is much richer. e theoretical and numerical results predict that the population x can be controlled using the state-dependent impulsive control strategies. Because h < x * , the method by impulsive harvesting population x was a unique choice which can keep the density of population x below the critical value. Our results shows that the control period is much lower when the released amount τ is much smaller, which means that the control frequency is much higher. But from the viewpoint of management for some ecological problems (e.g., pest outbreak, phytoplankton blooms, etc.), it is always hoped that the control period is much longer. us, much more predators are needed to be released. However, our results also indicate that chaos may occur for some values of τ, which is unfavourable for prediction of control period. Hence, the released amount of predator should be chosen carefully. Obviously, our results are much more interesting, and we expect that these results are helpful to the management of ecological problems.

Data Availability
e data used to support the findings of this study are included within the article.

Conflicts of Interest
e author declares that there are no conflicts of interest.