Dynamics and Optimal Harvesting Control for a Stochastic One-Predator-Two-Prey Time Delay System with Jumps

We consider a stochastic one-predator-two-prey harvesting model with time delays and Lévy jumps in this paper. Using the comparison theorem of stochastic differential equations and asymptotic approaches, sufficient conditions for persistence in mean and extinction of three species are derived. By analyzing the asymptotic invariant distribution, we study the variation of the persistent level of a population. Then we obtain the conditions of global attractivity and stability in distribution. Furthermore, making use of Hessian matrix method and optimal harvesting theory of differential equations, the explicit forms of optimal harvesting effort and maximum expectation of sustainable yield are obtained. Some numerical simulations are given to illustrate the theoretical results.


Introduction
Many researchers are widely focused on the complex dynamics of biological systems such as delay population systems [1][2][3], stochastic population systems [4][5][6][7][8][9][10][11][12], and impulsive population systems [13][14][15].Recently, many scholars have investigated two-species models and studied the extinction and persistence [16][17][18].However, in the real world, it is a common phenomenon that one predator catches two or more kinds of preys [19].Consequently, models with three or more species which can explain the dynamical behaviors of the population accurately are investigated [20,21].On the other hand, it is necessary and important to consider time delay caused by the competition and predation of species.Delayed differential equations can exhibit much more complex dynamics than differential equations without delay, and stable equilibrium can become unstable with the effects of a time delay.Therefore, many researchers have studied the Lotka-Volterra time delay models with two competitive preys and one predator [22,23].Notice that the composite population systems with stochastic effects and time delays present some complex dynamics; thus this causes widespread researchers concern [24][25][26][27][28][29][30].
In order to rationalize the model, some coefficients should be modified into the existing models.The onepredator-two-prey model with time delays is described by with initial value   () =   () ,   () > 0,   () ∈  [− 0 , 0] , (2) where   () denotes the size of the -th species of the prey at time  (=1, 2),  3 () denotes the size of the predator at Ñ (, dV) =  (, V) −  (dV) Ñ(, dV) is a compensated Poisson process and  is a Poisson counting measure with characteristic measure  on a measurable subset Y of (0, +∞) with (Y ) < ∞.The distribution of Lévy jump   () can be completely parameterized by (  ,   , ) and satisfies the property of infinite divisibility.It is characterized by its characteristic function Φ   () ; we can get a detailed explanation from the following Lévy-Khintchine formula [37].There are many other papers about stochastic models with Lévy jump; the readers could refer to [38,39] and references therein.Considering the inevitable situations in the real world, we assume that the intrinsic growth rates  1 and  2 and the death rate  3 of the model are perturbed by the Lévy jump to signify the sudden climate change, so we introduce the Lévy jump into the underlying stochastic model (1).
Taking the economic factors into account, reasonable natural resources management can increase sustainable production and profits.Therefore, harvesting models have been already used to exploited the optimal harvesting policies of renewable resources [40][41][42][43].We only consider to harvest preys  1 and  2 and ℎ 1 > 0 and ℎ 2 > 0 are the harvesting effort rates of  1 and  2 , respectively [42,43].From system (1), we can ( This paper is organized as follows.In Section 1, we formulate the model.We show the solution of model ( 5) is global and positive.We give the main theorems for persistence in mean and extinction under the model ( 5) and its proof in Section 3. In Section 4, we prove the global attractivity and stability in distribution.Our main aim of this paper is to investigate the optimal strategy of the proposed model; we give the conclusions in Section 5. Finally, we carry out numerical simulations and some figures to support the main conclusions in Section 6.

Global Positive Solution
In order to explore the dynamical behaviors of ecological population, we first study the positivity of the solutions of system (5).
with initial value   () = ln  () .It is not difficult to know that the coefficients of model ( 5) satisfy the local Lipschitz condition; therefore, model ( 5) has a unique local solution () on [0,   ], where   is the explosion time.Applying Itô's formula, we can get the following solution: It is the unique positive local solution to model (5).Now let us prove   = +∞.Thus, we introduce the following auxiliary model: with initial value Obviously Before starting proving, we state several hypotheses.We assume that which implies that the persistent ability of species 1 is better than that of species 2.
Assumption 4.  > 0,   > 0,  = 1, 2, 3 express that all the population could coexist under the condition that the model frees from stochastic noises. 33 > 0, Γ  > 0,  = 1, 2, which imply that the two prey populations could coexist when there is no environmental noises and the predators are absent, where   is the complement minor of   in the determinant , ,  = 1, 2, 3.
Therefore, if If Lemma 9.For model (20), we have the following.
Dividing both sides of ( 38), (39), and ( 40) by , we can obtain that Note that lim Firstly, we prove (a).Since  1 < 0,  2 < 0, using Lemma 7, yields Substituting ( 45) and ( 46) into (40) gives where  is small enough satisfying  3 +  < 0. Applying (i) in Lemma 7 gives Secondly, we prove (b) and (c).Using Lemma 7, since  1 ≥ 0 and  2 < 0, it is easy to obtain that lim and As a consequence, we can study and discuss the following model: We know that It is not difficult to obtain the system The proofs of (d) and (e) are similar to these of (b) and (c).
Then, we give the proofs of (f) and (g).Using Lemma Theorem 11.For system (5), define From Assumptions 4, 5, and 6, we can obtain the following results.Proof.Now we will give the proof of (I), and the proof of (II) is parallel to (I).We can get and, as a consequence, we can also compute that Firstly, we prove (i).Since  11 ,  12 , and  13 are positive, we can get Note that  1 < 0,  2 < 0. By Lemma 9, we get lim →+∞   () = 0, ..,  = 1, 2, 3.
Taking advantage of the above identity, we will get the following two-species models: It is easy to get the result with the similar proof of that in [40][41][42] lim Thirdly, we prove (iii).Denote ,  as the solution of the following equations: Consequently, According to Lemma 7, for arbitrarily given  > 0, there exists a  1 > 0, for all  >  1 , Multiplying ( 75), (76), and (77) by (−1), , and , respectively, and adding them, one can observe that for sufficiently large  such that  >  1 , Substituting ( 60) and ( 81) into (89) yields For sufficiently large , by  1 / Ã1 >  1 >  2 > 1, and choosing  > 0 to be sufficiently small, then we have Similarly, denote ,  as the solution of the following equations: Then we have (93) In the same way, we can choose a  2 > 0, for arbitrarily given  > 0, such that It follows that, for any sufficiently small , there exists  3 and  4 such that multiplying (75), (76), and (77) by , (−1) and , respectively, and adding them, we can obtain that Applying ( 60) and ( 81) into (95) yields that Note that  2 / Ã2 >  3 / Ã3 > 1.According to the arbitrariness of  and Lemma 7, we have Substituting ( 81), (91), and (97) into (77), For sufficiently large , then according to  3 < Ã3 , Lemma 7, and the arbitrariness of , it is not difficult to show lim →+∞  3 () = 0, .. (99) Consequently, model ( 5) reduces to the following model: which has already been investigated in [5].Then similarly to the proof of Theorem 5.1 in [5], the following identities can be derived: Fourthly, we prove (iv).Since  3 − Ã3 > 0, for arbitrariness of , the application of Lemma 7 to (98) yields Substituting (103) into (75), we can obtain that, for sufficiently large , According to the arbitrariness of  and Lemma 7, we have Analogously, we can prove that ⟨ 2 ()⟩ * ≥ ( 2 − Ã2 )/ is established.Substituting ( 103) and ( 105) into (77) yields that when  is large enough, According to the arbitrariness of  and Lemma 7, we have Subsequently, we have lim This completes the proof.

Stability in Distribution
For the convenience, we define the following notations: Lemma 13.For any  > 1, there exists a constant  = () which makes the solution () of model ( 5) with any given initial value satisfy the property that Proof.The proof is rather standard and hence is omitted.
Firstly we prove (i).It is easy to see that Λ 1 ∈ , so  is not empty.By (iv) of Theorem 11, for any  ∈ , we have Applying to Lemma 14, model ( 18) has a unique invariant measure (⋅).According to Corollary 3.4.3 in [46], we can get that (⋅) is strong mixing.At the same time, it is ergodic by Theorem 3.2.6 in [46].Hence, it then follows from (3.3.2) in [46] that Let () represent the stationary probability density of model (18); we obtain Since the invariant measure of model ( 18) is unique, then, according to the one-to-one correspondence between () and its corresponding invariant measure (⋅), we obtain That is to say, Assume that Λ 1 = ( 1 ,  2 ) T is the unique stagnation point of the following equation:  It holds that Λ 1 = [( −1 ) T + ] −1 .We can take use of the following Hessian matrix [42,43]: is negatively defined, so Λ 1 is the unique extreme point of ().In other words, if Λ 1 ∈ , i.e.,   ≥ 0 and Θ  > 0,  = 1, 2, then the optimal harvesting effort is  * = Λ 1 and  * is the maximum value of ESY.We are now in the position to prove (ii).Suppose that the optimal harvesting effort  consider the optimal harvesting of preys [50].To begin with, we establish the modified model.In Theorem 11, we obtain the sufficient criteria for extinction and persistence in mean of each species.Lévy jump is important to study.When  1 < 1, all species go to extinction.When  1 > 1 >  2 , then the predator  3 goes to extinction and  2 is persistent in mean,  1 goes to extinction.When  2 > 1 >  3 , then  3 goes to extinction and  1 ,  2 are persistent in mean.When  3 > 1, then three species are persistent in mean.The main purpose of this paper is to study the optimal harvesting.After discussing the stability of distribution, we study the optimal harvesting and obtain the maximum yield of two preys.
From the numerical simulations, we list the following biological meanings.
(1) We find the noise can cause the variation of species.When the noise is large, it in reality can suppress the increase of population, then it dies out.
(2) Time delay and Lévy jump have important effects on the persistence in mean and the harvesting yield.
In traditional papers, scholars consider two species or three species without optimal harvesting.We consider a three-species model with Lévy jump and optimal harvesting.Time delay is ineluctable in the ecological environment and is necessary to consider delays.Recently, stochastic models with the telephone noise have been studied by many authors  The optimal harvesting effort and the maximum of ESY [45].Red line is with ℎ 1 = 0.9778, ℎ 2 = 0.9243, blue line is with ℎ 1 = 0.5, ℎ 2 = 0.4, and green line is with ℎ 1 = 0.2, ℎ 2 = 0.1.[51,52].In the future research, we hope to add more realistic conditions and study more interesting topics, for example, pulse process, Markov Chain, telephone noise, and partial differential system [53][54][55].

Figure 4 :
Figure 4: The persistence in mean of two preys and of three species are given in Theorem 11.(b) stands for the phase portrait of (a).The density of white noises is taken:  1 = 0.8,  2 = 0.65,  3 = 0.42.