The Study of a Predator-Prey Model with Fear Effect Based on State-Dependent Harvesting Strategy

In presence of predator population, the prey population may significantly change their behavior. Fear for predator population enhances the survival probability of prey population, and it can greatly reduce the reproduction of prey population. In this study, we propose a predator-prey fishery model introducing the cost of fear into prey reproduction with Holling type-II functional response and prey-dependent harvesting and investigate the global dynamics of the proposed model. For the system without harvest, it is shown that the level of fear may alter the stability of the positive equilibrium, and an expression of fear critical level is characterized. For the harvest system, the existence of the semitrivial order-1 periodic solution and positive order-q (q≥ 1) periodic solution is discussed by the construction of a Poincaré map on the phase set, and the threshold conditions are given, which can not only transform state-dependent harvesting into a cycle one but also provide a possibility to determine the harvest frequency. In addition, to ensure a certain robustness of the adopted harvest policy, the threshold condition for the stability of the order-q periodic solution is given. Meanwhile, to achieve a good economic profit, an optimization problem is formulated and the optimum harvest level is obtained. Mathematical findings have been validated in numerical simulation by MATLAB. Different effects of different harvest levels and different fear levels have been demonstrated by depicting figures in numerical simulation using MATLAB.


Introduction
Prey-predator interaction is a crucial topic in theoretical ecology and evolutionary biology. e history of the study about the prey-predator interactions dates back long. e pioneering work to describe the prey-predator interactions in mathematics belongs to the Lotka-Volterra model [1,2]. Subsequently the model was improved by adding logistic growth term for the prey and variety of population-dependent response functions [3][4][5][6][7][8][9][10][11][12][13][14][15]. A prototype model that captures the prey-predator interaction takes the form dx dt � bx − dx − cx 2 − yp(x, y), dy dt � eϕ(p(x, y))y − my, where x(t) and y(t) represent the densities of prey and predator population, respectively, b, d, (d < b) and c represent the birth rate, natural death rate, and density-dependent decay rate due to the intraspecies competition, respectively, p(x, y) represents the functional response, e is the efficiency of conversion, m is natural mortality of predator, and ϕ is a monotonically increasing function. Due to prey-predator interactions, predators always have an impact, direct, indirect, or both, on prey population. In model (1), the term yp(x, y) models the direct impact of predator on prey by catching and killing behavior. Meanwhile fear of predation risk can be regarded as the indirect impact of predator on prey, and some theoretical ecologists and biologists have realised that a prey-predator model should involve not only direct killing but also the fear [16,17]. e fieldwork of Zanette et al. [18] on song sparrows observed the impact of fear and found a reduction in reproduction by 40% in the number of offspring due to the fear of predation. Based on this phenomenon, Wang et al. [19] incorporated a predator-dependent fear factor into the birth rate of prey in model (1) (i.e., replace b by b(y) � b/(1 + ky)) with linear and Holling type-II functional response to explore the effect of fear on population dynamics. e results show that high level of fear could stabilize the system. Das and Samanta [20] investigated the impact of fear in exponential form on a stochastic prey-predator system when the predator is provided additional food. Sahoo and Samanta [21] investigated a two prey-one predator model by including the cost of fear into prey reproduction and switching mechanism in predation. Das et al. [22] developed and explored a predator-prey model incorporating the cost of perceived fear into the birth and death rates of prey species with Holling type-II functional response. Sarkar and Khajanchi [23] and Kumar and Kumari [24] incorporated a form of fear factor into the birth rate of prey by assuming a nonzero minimum cost of fear. e impact of fear has also been investigated on prey-predator systems with prey refuge [25][26][27], Allee effect [26], hunting cooperation [28], and additional food resource for predator [20,29].
e study of resource management including fisheries, forestry, and wildlife management has great importance. It is necessary to harvest the population but harvesting should be regulated in such a way that ecological sustainability as well as conservation of the species can be implemented in a long run. Besides, it is always hoped that the sustained ability can be achieved at a high level of productivity and good economic profit. In the past decade, scholars considered different kinds of harvest on the dynamics of the predator-prey system such as continuous harvesting [30][31][32][33] and intermittent harvesting [34][35][36][37]. Compared to fixed time harvest strategy, the state-dependent harvesting strategy takes the existing resources of species into full consideration and can maintain the sustainability of species in certain level. Statedependent harvested system can be described by the impulsive semidynamical system [38][39][40][41][42][43][44][45]. Recently, Lai et al. [46] proposed and studied a Lotka-Volterra predator-prey system incorporating both continuous harvesting and fear effect; that is, where k is the level of fear, E is the fishing effort used to harvest, q is the catchability coefficient, and a 1 and a 2 are constants. e harvest in (2) is continuous and in Michaelis-Menten type. However, in reality, the harvest of species should consider the aspect of ecological sustainability as well as conservation. us, in most cases, species are caught intermittently, not continuously.
To the best of our knowledge, to this day, still no scholars investigated the dynamic behavior of the predatorprey system incorporating both fear effect and intermittent harvesting, which motivated us to study a predator-prey model incorporating fear effect based on state-dependent harvesting strategy. e aim of this study is to check the influence of fear level on the stability of the positive steady state of the system without harvest. Meanwhile, for the harvest system, it mainly discusses the existence of the order-q (q > 0) periodic solution, since it provides a possibility to transform the state-dependent harvesting into a cycle one. Meanwhile, in order to make a maximum economic profit in the harvest process, the optimal control problem is discussed. e organization of this study is as follows. In the next section, we introduce the mathematical model for predator-prey system with fear effect based on state-dependent harvesting strategy. In the same section, we present some preliminaries used in the discussion of the system dynamics. Section 3 is dedicated to the existence and stability of semitrivial order-1 and positive order-1 periodic solution. We also study the existence of order-2 and order-3 periodic solution. In Section 4, we demonstrate different effects of different harvest levels and different fear levels by depicting figures in numerical simulation using MATLAB. e paper concludes in Section 5, in which we briefly summarize the biological indications of our analytical findings.

Model Formulation.
In presence of predator population, the prey population may significantly change their behavior. Fear for predator population enhances the survival probability of prey population, and it can greatly reduce the reproduction of prey population [23]. In this study, we consider a predator-prey model introducing the cost of fear into prey reproduction with Holling type-II functional response and a saturation function ϕ in equation (1); that is, Where the variables, model parameters, and their units/ dimensions are given in Table 1. To achieve the commercial purpose of the fishery, it is necessary to harvest the population in such a way that ecological sustainability as well as conservation of the species can be implemented in a long run. e harvest can be continuous or intermittent. In this work, a state-dependent harvest strategy is considered. Let l be the harvest level of prey population; that is, when the density of prey population reaches level l, the harvest is implemented, resulting in a portion of prey and predator being caught. Let E denote the harvest effort, which is dependent on the harvest level l, and let q 1 and q 2 be the catchability coefficients of prey and predator populations.
In addition, to avoid the extinction of predator, it is necessary to release a quantity of predator pups, denoted by τ, which is also dependent on level l. Based on this consideration, the model with state-dependent harvesting takes the following form: 2 Complexity Denote K ≜ (b − d)/c. en K is the carrying capacity of prey population in absence of predator. System (4) is considered in the domain S � (x, y)|0 ≤ x ≤ K, y ≥ 0 for ecological practices. e purpose of this paper is to analyze the dynamics of system (4). Besides, it is always hoped that the harvest can be achieved at a good economic profit, and this requires determining an optimal harvest level l. Next, some preliminaries are listed for the analysis of the harvest model (4).

Preliminaries.
Let us consider a general planar system: where (x, y) ∈ Ω ⊂ R 2 , and χ(x, y) � 0 describes the states at which the harvest is implemented; α and β describe the effects of the harvest strategy. P(x, y) and Q(x, y) are arbitrarily derivative with respect to (x, y) ∈ Ω; χ, α, and β are linearly dependent on x and y; that is, χ x , χ y , α x , α y , β x , and β y are constant. e dynamic system constituted by the solution mapping defined by system (5) is called an impulsive semicontinuous dynamic system, denoted as (Ω, π; I, M IMP ), where π � (π 1 , π 2 ): Ω × R ⟶ Ω, M IMP ≜ (x, y)| χ(x, y) � 0}, and Let z(t) � (x(t), y(t)) be the solution of system (5) with initial value z(0) Definition 1 (priodic solution [47][48][49]). e solution z � z(t) of system (5) is said to be periodic if there exists positive integer m ≥ 1 such that z m � z 0 . Denote k ≜ min m ∈ N, z m � z 0 ; then orbit c(z) is said to be an order-k periodic orbit of system (5).
Definition 2 (orbitally stable [47][48][49]). An orbitc(z) is said to be orbitally stable if, for any ε > 0, there is a neighborhood V of z so that, for all z in V, there is a reparameterization of time (a smooth, monotonic function) t(t) such that |z(t) − z(t(t))| < ε for all t ≥ t 0 .

Remark 1.
If there exists a point L ∈ N ω and q > 0 such that ϕ q N (L) � L and ϕ j N (L) ≠ L (j < q), that is, L is a fixed point of ϕ q N , then system (4) admits an order-qperiodic solution.

Main Works
Define en e M is called the critical value of the conversion; that is, when e ≤ e M , the conversion is not enough to maintain the survival of predator and the predator population will go to extinction. us, in this work, it is reasonable to assume that e > e M .
Besides, denote For system (3), the following result holds.

Theorem 1.
ere are three equilibria for system (3) whene > e M : two boundary saddlesO(0, 0)andE(K, 0)and one positive equilibrium E * (x * , y * ). Moreover, one of the two following cases holds: It is obvious that O(0, 0) and E(K, 0) are saddles. At the equilibrium E * (x * , y * ), the characteristic equation us, the positive equilibrium E * (x * , y * ) is locally asymptotically stable in case of k > k * and unstable in case of k < k * . In this case, there exists a unique stable limit cycle for system (3).
In this case, system (3) is reduced to the following system: en there is x(T) � l and x(T + ) � (1 − q 1 E)l � x 0 by impulse effect. us, the following result holds. , which is orbitally asymptotically stable when R 0 < 1, where Proof. To discuss the stability of (ξ(t), η(t)), let us consider a small disturbance δ 0 . e trajectory starting from . is disturbed trajectory first intersects the harvest set M IMP at point B 2 (l, y 1 ) when t � T + δt, and then it jumps to point Setting δy 1 ≜ y 1 , for 0 < t < T, the variables δx and δy can be expressed by the relation where Φ(t) � (ϕ ij ) 2×2 is the fundamental solution satisfying the variation equation.
According to the first-order Taylor expansion on η(t), there is δy 1 Complexity By impulse effect, there is (10) holds, there is δ 1 < δ 0 . By the arbitrary of δ 0 , it concludes that the order-1 semitrivial periodic solution is orbitally asymptotically stable.

Corollary 1.
e semitrivial order-1 periodic solution(ξ(t), η(t)) is orbitally asymptotically stable if one of the two following cases is satisfied: 3.2. Positive Order-K Periodic Solution for τ > 0. Since the harvest may cause the extinction of predator when τ � 0, in order to keep the predator species from going extinct, it is necessary to reduce the harvest strength and release a certain quantity of predator pups.
For 0 ≤ x ≤ K, define Let N 0 and M 0 denote the intersection point between y � y L (x) and the phase set N PHA and the harvest set M IMP , respectively; G 0 denotes the intersection point between y � τ and the phase set By equation (22), the function y � y(x, y s ) can be expressed as follows: Define l ≜ max l|y(l, y N 0 ) ≤ y L (l) . When l ≤ l, the trajectory of system (4) starting from N 0 will intersect the harvest set M IMP , and denote the intersection point by N 1 ; that is,

Existence of Order-1 Periodic Solution.
By eorem 1, the dynamic behavior of system (3) varies with the model parameter k. us, the discussions will be divided according to parameter k and harvest level l.
is only a function of y. Next, the Poincaré map ϕ N will be characterized and its main property will be analyzed.
there exists a unique y s ′ ∈ (0, y N 0 ) and T s > 0 such that Property 1. For system (4), when l ≤ x * , the Poincaré map ϕ N defined by equation (24) has the following properties: e following result holds.

Theorem 3.
ere exists a unique positive order-1 periodic solution for system (4) Proof. By Remark 1, the existence of order-1 periodic solution is equivalent to the existence of a point L ∈ N PHA such that y L is a fixed point of ϕ N . By Property 1 (i), f N is continuous on [0, +∞). Since f N (0) � ϕ N (0) � τ > 0 and f N (y) � ϕ N (y) − y ⟶ − ∞ as y ⟶ +∞, by the intermediary property of continuous function, there exists at least one y L > 0 such that f N (y L ) � 0; that is, us, the trajectory of system (4) 6 Complexity starting from L((1 − q 1 E)l, y L ) forms an order-1 periodic orbit. Next, the location and uniqueness of the order-1 periodic orbit will be analyzed. By Property 1 (i), ϕ N achieves its . en there is y L ∈ (y L ((1 − q 1 E)l), y min ) and y L is unique, as shown in Figure 1(b)). □ Case I-2: l > x * : in case of x * < l ≤ l, the trajectory of system (4) starting from N 0 will intersect the harvest set M IMP . When l > l, the trajectory starting from point N 0 does not intersect the harvest set M IMP . Denote

Theorem 4.
ere exists a positive order-1 periodic solution for system (4) Proof. When x * < l ≤ l, similar to the proof of eorem 3, system (4) admits an order-1 periodic solution. For l > l, if Combining with ϕ N (y G 0 ) > y G 0 , it can be concluded that system (4) admits an order-1 periodic so- such that ϕ N (y L ) � y L ; that is, system (4) admits an order-1 periodic solution. □ □ Case II: k < k * : in this case, the trajectory of system (4) starting from N 0 will intersect the harvest set M IMP .

Theorem 5.
ere exists a positive order-1 periodic solution for system (4)
us, for any y 0 , a sequence y k k�1,2,... is obtained under ϕ N ; that is, y k � ϕ N (y k− 1 ). If y 0 < y L , then y k is a monotonically increasing sequence with y k < y L , so the limit is y L . Similarly, if y 0 ∈ (y L , y N 0 ], then y k is a monotonically bounded decreasing sequence, and the limit is y L . If y 0 > y N 0 , then y 1 � ϕ N (y 0 ) ∈ (0, y N 0 ); thus y k k�1,2,... is a monotonically bounded sequence with limit y L . To sum up, by the arbitrariness of y 0 , the order-1 periodic solution is globally attractive and so is globally orbitally asymptotically stable. □ 3.2.3. Existence of Order-q (q ≥ 2) Periodic Solution. For l ≤ x * , by eorem 3, if τ ≤ τ f , the order-1 periodic solution is orbitally asymptotically stable and globally attractive, which means that system (4) does not admit order-q (q ≥ 2) periodic solution. For τ > τ f , there exists unique
Proof. "Necessity." Proof by contradiction. Assume that     (4) for l � 60%K and τ 2 � 2 with parameters given in numerical section. e phase portrait diagram shows that the trajectory of system (4) will eventually tend to the positive equilibrium E * . In this case, the positive equilibrium E * is globally asymptotically stable.   Figure 6: e phase portrait diagram of model (4) for τ 2 � 9 and different harvest level l with parameters given in numerical section: (a) l � 60%K; (b) l � 80%K. For l � 60%K, there is τ 2 > τ M 2 , and system (4) admits an order-1 periodic solution; for l � 80%K, there is τ 2 ∈ (τ M 1 , τ M 2 ), and the trajectory of system (4) will tend to the positive equilibrium E * after finite impulses. 10 Complexity   Figure 8: Illustration of the Poincaré maps ϕ N and ϕ 2 N of system (4) for l � 25%K and τ 2 � 8 and different values of q 1 with parameters given in numerical section. For q 1 � 0.8, system (4) admits a unique orbitally asymptotically stable order-1 periodic solution; for q 1 � 0.5, there is ϕ 2 N (4.92) > 4.92 and μ 1 > 1; in this case, system (4) admits a stable order-2 periodic solution; (c) for q 1 � 0.2, there is ϕ 2 N (4.92) < 4.92, and system (4) admits a stable order-2 periodic solution.

) admits a stable order-1 periodic solution or a stable order-2 periodic solution. Moreover, there is
τ 2 =12 v=y y=y N0 Figure 9: Illustration of the Poincaré map ϕ N of system (4) for l � 25%K and different values of τ 2 with parameters given in numerical section. e results indicate that ϕ N has a unique fixed point for any τ 2 ; that is, system (4) admits a unique order-1 periodic solution.
can be concluded that ϕ 3 N (y) � y if and only if y � y L 2 ; that is, the order-3 periodic solution does not exist.
If y N 1 < ϕ 2 N (y N 0 ) < y N 0 , system (4) simultaneously admits an order-1 periodic solution and order-2 periodic solution. Moreover, there is ϕ N (y N 0 ) > y N 2 . If τ ≤ y N 1 , that is, ϕ N (τ) ≤ y N 0 , there exist y R 1 ∈ [0, y N 1 ], y R 2 ∈ (y L 1 , y N 0 ), y R 3 ∈ (y N 0 , y L 2 ), and y R 4 ∈ (y N 2 , +∞) such that   Figure 11: Illustration of functions ϕ 3 N (y) and ϕ N (y) of system (4) for l � 50%K and τ 2 � 2.9 with parameters given in numerical section. e time series evolution of prey population x(t), predator population y(y), and the phase portrait diagram demonstrate the order-3 periodic solution. easily checked that ϕ 3 N is increasing on [0, y R 1 ], [y N 1 , y R 2 ], [y N 0 , y R 3 ], and [y N 2 , y R 4 ], and ϕ 3 N is decreasing on 4 , it can be concluded that ϕ 3 N (y) � y if and only if y � y L 2 ; that is, the order-3 periodic solution does not exist. " If ϕ 2 N (y N 0 ) � y N 1 , then there exists an order-3 periodic solution since ϕ N (y N 1 ) ≠ y N 1 . If ϕ 2 N (y N 0 ) < y N 1 , then there exists at least one y f ∈ (y R 1 , y N 1 ) such that ϕ 3 N (y f ) � y f and ϕ N (y f ) > y f ; that is, system (4) admits an order-3 periodic solution.
Next, the number of order-3 periodic solutions will be discussed: N (y N 0 ). Moreover, there are y N 0 � ϕ N (y N 1 ) and y N 1 � ϕ N (ϕ N (y N 0 )); that is, system (4) admits at least one order-3 periodic solution.

Optimal Harvest Level Determination.
To achieve the commercial purpose of the fishery, it is necessary to harvest the population, and it is always hoped that the sustained ability can be achieved at a good economic profit. For the harvest problem, it is necessary to determine the controlled values E and τ and harvest level l, and this in general involves the optimization theory [48,49]. Let l be the harvest level, which is a decision variable. eorems 3 and 4 show that system (4) admits an order-1 periodic solution when l ≤ l or l > l with τ ≤ y M 1 − (1 − q 2 E)y L (l). Since the harvest effort and yield of released predator are dependent on the harvest level, then it is assumed that E(l) and τ(l) take the following forms: where l 1 ≤ l ≤ l 2 , l 1 and l 2 are minimum and maximum of the harvest level, and E 1 (E 2 ) and τ 1 (τ 2 ) are the harvest effort and yield of released predator at the harvest level l 1 (l 2 ). Let c 1 and c 2 be the unit selling prices of prey and predator, let c 3 be the unit cost of harvest, and let c 4 be the unit cost in breeding predator. en the benefits from harvest can be described as e objective is to maximize the unit benefits; that is,

Numerical Simulations and Optimization
In this section, we compute some numerical simulations regarding the existence and stability of the periodic solution for the predator-prey model (4).  (6), it can be observed that the fear effect factor k only affects the value y * . Figure 2 illustrates the dependence of isoline dx/dt � 0 on x for different fear level k. As illustrated, the positive equilibrium becomes stable from unstable with increasing of k. By equation (6), there is k * � 0.0376. To verify the theoretical results obtained in the above section, the simulations are implemented by considering different combinations of k, q 1 , q 2 , τ 2 , and l. 14 Complexity Case I: k � 0.04, q 1 � 0.8, and q 2 � 0.6. I-(1): τ 2 � 0. By eorem 2, system (4) admits an order-1 semitrivial periodic solution for any l ≤ K � 100, which is expressed by equation (10). Moreover, for l ≤ l � 34.45, by Corollary 1, the semitrivial order-1 periodic solution is orbitally asymptotically stable. Meanwhile, for l � 50%K, by eorem 2, there is R 0 � 0.5673 < 1; that is, the order-1 semitrivial periodic solution is orbitally asymptotically stable, as presented in Figure 3.
I-(2): τ 2 > 0. Firstly, for τ 2 � 2, there is l � 34.45%K. For l � 25%K ≤ x * , eorem 4 and eorem 7 indicate that system (4) admits a globally asymptotically stable positive order-1 periodic solution, as illustrated in Figure 1. Meanwhile, for l > x * , function ϕ N of system (4) for l � 30%K, l � 50%K, and l � 60%K is presented in Figure 4. It can be observed that system (4) admits an order-1 periodic solution for l � 30%K and l � 50%K since the inequality τ < y M 1 − (1 − q 2 E)y L (l) in eorem 4 holds. When l � 60%K, the direction of the inequality has changed, and the trajectory of system (4) will tend to the positive equilibrium E * (25, 4.93) after finite impulses, as shown in Figure 5. But this is not always true; it is dependent on the value of τ 2 . As τ 2 goes up to 9, the inequality τ > y M 2 − (1 − q 2 E)y L (l) in eorem 4 holds; then system (4) admits a positive order-1 periodic solution, as shown in Figure 6(a). However, for l � 80%K, the direction of the inequality is changed again, and the trajectory of system (4) will tend to the positive equilibrium E * (25, 4.93) after finite impulses, as shown in Figure 6(b).
Case II. k � 0.04 and q 2 � 0. II-(1): τ 2 � 0 and q 1 � 0.8. Notice from Figure 1 that a higher catching rate for predators (e.g., q 2 � 0.6) will cause the predator species extinction for l � 50%K. When the catch for the predator is very small or ignored, that is, q 2 � 0, the order-1 semitrivial periodic solution is unstable by Corollary 1 (i.e., R 0 � 1.02 > 1). e function f N (y) � ϕ N (y) − y is presented in Figure 7 and it can be observed that f N (y) > 0, and, for any initial condition, the trajectory of system (4) will eventually tend to the positive equilibrium E * (25, 4.93).
II-(2): τ 2 > 0. For τ 2 � 8 and q 1 � 0.8, by eorem 3, system (4) admits a unique positive order-1 periodic solution for l � 25%K. To show the existence of order-2 periodic solution, the catching rate for prey q 1 is selected as a key parameter to verify how does the dynamic behavior of the system change. Function ϕ N for different catching rate for prey is presented in Figure 8. It can be observed that system (4) admits a unique globally asymptotically stable order-1 periodic solution for a higher catching rate for prey, for example, q 1 � 0.8, as shown in Figure 8(a). As the catching rate for prey goes down, for example   1.2612 > 1 in eorem 8 holds; then system (4) admits a stable order-2 periodic solution, as shown in Figure 8(b). Meanwhile, for q 1 � 0.2, condition (i) ϕ 2 N (y N 0 ) < y N 0 in eorem 8 holds, and then system (4) admits an order-2 periodic solution, as shown in Figure 8(c).
Case III: k � 0.01 and q 1 � q 2 � 0.6. e positive equilibrium E * becomes unstable when k � 0.01, and system (3) admits a limit cycle Γ LC . Since the existence and stability of semitrivial order-1 periodic solution for τ � 0 do not depend on the fear factor k, the results are the same as those in Case I-(1) and are omitted hereby. So it mainly discusses the dynamic behavior of system (4) for τ 2 > 0. It is easily checked that l � 32.7%K.
Next, let us consider l � 50%K > l. It is easily checked that τ M 1 � 1.2495 (i.e., τ 2 � 2.1867); then, by eorem 4, there exists an order-1 periodic solution for system (4) when τ 2 ≤ 2.1867. Here it should be pointed out that the condition given in eorem 4 is only a sufficient one; in fact, as long as τ 2 ≤ 2.55, system (4) admits an order-1 periodic solution, as illustrated in Figure 10.

Optimization.
To achieve a good economic profit, it is necessary to find a level l * at which the benefits from harvest are maximal. Assume that c 1 � 100, c 2 � 5c 1 � 500, and c 3 � 20%c 2 � 100. Denote σ ≜ c 3 /c 1 . In order to sustain the harvest, the releasing yield of predator should not be too large, so, in this part, it is assumed that τ 2 � 1.
For k � 0.04, system (4) admits an order-1 periodic solution when l ≤ 70%K; the dependencies of period T and benefit F benefit on the harvest level l are presented in Figure 16. It can be seen that period T increases as l increases. Meanwhile the benefit function F benefit climbs up and then declines as l increases. For different σ, F benefit achieves its maximum at different l * σ . When the cost of harvest is ignored, that is, σ � 0, there is l * 0 � 40%K. As σ goes up, F benefit goes down. When σ � 8, F benefit achieves its maximum at l * 8 � 50%K.
For the case of k � 0.01, system (4) admits an order-1 periodic solution when l ≤ l 2 , and the dependencies of period T and benefit F benefit on the harvest level l are presented in Figure 17. e benefit funtion F benefit climbs up and then declines as l increases. For different σ, F benefit almost achieves its maximum at l * σ � 50%K.

Conclusion and Discussions
In this paper, we have discussed the dynamics of a harvested prey-predator model, where the prey is provided with fear effect. For the system without harvest (3), there exists a unique positive equilibrium. To verify the stability of the equilibrium, a critical level of fear factor k * is characterized (i.e., equation (6)). When the impact of fear on prey is small, that is, 0 ≤ k < k * , the positive equilibrium is unstable and a limit cycle exists. As the impact of fear grows and exceeds k * , the positive equilibrium becomes stable and the limit cycle disappears (Figure 2). In any case, predators coexist with prey and the system is persistent.
For the system with harvest (4), if we do not consider the release of predator (i.e., τ � 0), system (4) admits a semitrivial order-1 periodic solution for any harvest level (Figure 3). Moreover, the semitrivial order-1 periodic solution is orbitally stability when the harvest level is not higher than the first equilibrium component (i.e., l ≤ x * ). Meanwhile, for the case of l > x * , the semitrivial order-1 periodic solution is orbitally stable when a strong harvest intensity is implemented. is means that the system can be disrupted, and predators will go extinct if the harvest is not properly planned. To maintain the ecological health and avoid the extinction of predator populations, it is necessary to reduce the catch rate of predators (Figure 7) or release a certain quantity of predator pups. In the second case, that is, τ > 0, system (4) admits a positive order-1 periodic solution when l ≤ l (Figures 1, 4, and 9). Moreover, the order-1 periodic solution is orbitally asymptotically stable and globally attractive when l ≤ x * and τ ≤ τ f . Meanwhile, for τ > τ f , system (4) admits an order-2 periodic solution in case of ϕ 2 N (y N 0 ) ≥ y N 0 and μ 1 2 > 1 (Figure 8). In case of l > l, system (4) also admits a positive order-1 periodic solution for τ ≤ τ M 1 (Figures 4 and 10) or τ ≥ τ M 2 ( Figure 6). But, for τ ∈ (τ M 1 , τ M 2 ), when k > k * , the trajectory of system (4) will tend to the positive equilibrium E * (25, 4.93) after finite impulses ( Figure 5). Meanwhile, in case of k < k * , the dynamic behavior of system (4) depends heavily on parameter τ. For different value of τ, system (4) may admit an orderk(k � 1, 2, 3, 4, 5) periodic solution (Figures 11-15).
To achieve a good economic profit, optimization with τ 2 � 1 is carried out and the results show that the benefits from harvest depend on the unit selling prices of prey and predator, as well as the unit cost of harvest. For given c 1 � 100, c 2 � 500, and c 4 � 100, the benefit function first climbs up and then declines as l increases. For k � 0.04, the economic profit F benefit achieves its maximum at l * � 40%K when σ � 2. As σ goes up, F benefit goes down, and F benefit achieves its maximum at l * � 50%K when σ � 8. Meanwhile, in case of k � 0.01, F benefit almost achieves its maximum at l * � 50%K.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that there are no conflicts of interest.