Dynamics of a Guanaco-sheep Competitive System with Unilateral and Bilateral Control

In this paper, based on a guanaco-sheep competitive system, we develop and analyze mathematical models with unilateral and bilateral control for the management of overgrazing. We ﬁrst analyze the dynamics of the uncontrolled system. Then, we investigate the system with impulsive control by diﬀerential equation geometry theory. And we mainly prove the existence and stability of order-1 periodic solution for unilateral control system and order-2 periodic solution for bilateral control system. Comparing the unilateral and bilateral control strategy, we encourage bilateral control rather than unilateral control for the management of sheep species.

only 26% of its original scope. And the results of [4] suggested that the range distribution of the guanacos has been reduced by 60% in Argentina, 75% in Chile.
This phenomenon is mainly attributed to habitat change and loss, competition with livestock (sheep, cattle and horses) for pasture, and the lack of a proper management for livestock [5].
Literature study shows that there is a competitive relationship between the guanaco and sheep species. Puig et al. [6] indicated that competition between guanaco and domestic animals is due to the fact that their diets overlap up to 83%. Baldi et al. [7] and Gordon [8] suggested that guanaco and sheep compete for food and water resources. Ranchers believe that guanaco depletes food supplies, hence reducing productivity of sheep species [9]. So it is very meaningful to investigate the interaction of these two species. In this regard, Laguna [10] developed a three-dimensional mathematical model with two competing herbivores (sheep and guanaco) and one predator, and analyzed the steady state of the occupation of patches for the three species. Daza et al. [11] proposed a delayed differential model capturing the interaction between sheep and guanaco, and investigated coexistence and extinction of these two species.
However, overgrazing or domestic animal overstock has caused vegetation degradation, and the reduction of wildlife is the most serious consequences of vegetation degradation. Pedrana et al. [1] indicated that sheep species seems to have a higher impact on vegetation. Marino et al. [12] showed that grassland degradation and the decline of guanaco species will continue due to overgrazing. Bharucha et al. [13] pointed out that overgrazing modified the plant cover and impaired its productive capacity. So what management programs should be developed by human beings to change this situation and protect guanaco species and grassland? With human intervention, the population size of sheep and guanaco species will change in a short time, and this phenomenon could be well described by impulsive differential equation [14] which has been extensively applied in applications [15][16][17][18][19][20][21][22][23][24][25][26][27]. For example, integrated pest management models were studied in [15][16][17] where the existence of order-1 periodic solution is proved, and the authors of [18] and [19] also considered the optimal control problem. For resource management, Guo et al. [22] analyzed a algae-fish with unilateral control. And bilateral control strategy is applied in [23,25,27] where the existence and stability of order-2 periodic solution are mainly proved. And Huang et al. [27] also analyzed the attraction region of the system, and performed bifurcation diagrams of periodic solution. In addition, bilateral control strategy has also been applied to biodiversity conservation and commodity price adjustment [24,26].
Motivated by the above research, in this paper, by constructing semi-continuous dynamical systems, we develop mathematical models for the management of sheep and guanaco species. On the one hand, the number of sheep species will be limited to avoid overgrazing; on the other hand, the economic interests of farmers should also be protected while keeping sheep from overgrazing. Based on these two aspects, we consider unilateral and bilateral control strategies to keep the sheep species in a suitable range.
The rest of the paper is organised as follows. In section 2, we develop guanacosheep competitive systems with unilateral and bilateral control strategy. In section 3, the existence of the order-1 and order-2 periodic orbit is proved by successor function. And stability of the order-1 and order-2 periodic orbit is verified by analogue of the Poincaré criterion. In section 4, the numerical simulations are carried out with a special case to validate the theoretical results. Finally, a brief conclusion is presented in section 5.

Model formulation
For the guanaco and sheep species, their competitive ability is related to their body size. Tilman et al. [28] considered that the competitiveness of these two competing herbivores is not in the same hierarchi. And Daza et al. [11] pointed out that there exists hierarchical competition between guanaco and sheep, the guanaco is assumed to be the superior species and sheep is the inferior one. So we assume that the effect of inferior species on superior species satisfies Holling II functional response, and the influence of higher competitors on lower competitors satisfies a linear function. Based on the classical two-species competitive system [29], we get a modified competitive system as follows: where x(t) and y(t) stand for the density of guanaco and sheep species, respectively.
The parameters a, b, c, d, e and f have the same biological meaning as in [29].
As the application of semi-continuous dynamic system in overgrazing manage- and When two kinds of control are considered simultaneously, we have the following competitive system with bilateral control: All the parameters in systems (2), (3) and (4) are positive, and 0 < α 1 , α 2 < 1.

Main results
In this section, we mainly discuss the dynamical properties of system (1) and the existence and stability of periodic solution of systems (2), (3) and (4).

Dynamics of free system (1)
For system (1), it is easy to determine that trivial equilibrium O(0, 0) is always unstable. And the equilibrium E 1 ( a b , 0) is locally stable if ae > bd, otherwise it is unstable. The equilibrium E 2 (0, d f ) is locally stable if a(f + d) > cd, otherwise it is unstable. Positive equilibrium (x * , y * ) satisfies the following equation: After simple calculation, we have x * = d − f y * e , and y * satisfies the following second-order algebraic equation: Obviously, the curve l 2 := (x, y)|f (x, y) = 0 has a horizontal asymptote y = −1 and a vertical asymptote y = a−c b , respectively. Combining with differential equation geometry theory, we analyze the existence and stability of the positive equilibrium of system (1) from the following two cases.
lies on the right of y-axis.
The location of the two isoclines l 1 and l 2 is shown in Figure 1.
(b) Fig. 1 The location of the two isoclines l 1 and l 2 a > c.
exists unique positive equilibrium E * 1 (x * 1 , y * 1 ) (see Figure 1(a)). The Jacobian So E * 1 is locally stable. By constructing Dulac function D(x, y) = 1 xy , we find that system (1) has no limit cycle. So we conclude that E * 1 is globally stable.
lies on the left of y-axis.

Existence and stability of the order-1 periodic solution
In this work, we mainly prove the existence and stability of periodic solution when the unique positive equilibrium E * 1 (x * 1 , y * 1 ) is a saddle point. It then follows that the proof the existence of order-1 periodic solution of system (2). (2) As is shown in Figure 3

Existence of the order-1 periodic solution of system
The location of the two isoclines l 1 and l 2 for a < c.
respectively. We denote the coordinates of points A 2 and B 2 as (x A 2 , (1 + α 1 )h 1 ) and (x B 2 , h 1 ). We now give the existence of the order-1 periodic solution of system (2). (2) exists an order-1 periodic solution.
Proof Denote the intersection point of the curve l 2 and the line M 1 by B 1 . According to the phase diagram of system (2), there must exist a point such that the trajectory from A 1 intersects with M 1 at point B 1 , then jumps to only need to find another point A 2 ∈ N 1 such that G(A 2 ) > 0. Consider the intersection point A 2 of the dashed line l 3 and N 1 , the trajectory from A 2 intersects with N 1 at point B 2 , and then jumps to A + 2 due to the impulsive effect. Based as is shown in Figure 4. According to the continuity of successor functions, there must exist a point A ∈ N 1 between A 1 and A 2 satisfying that G(A) = 0, i.e., the trajectory AB ∪ BA forms an order-1 periodic solution of system (2). Fig. 3 The existence of order-1 periodic solution of system (2).

Stability of order-1 periodic solution of system
Theorem 2 The order-1 periodic solution of system (2) is orbitally asymptotically By simple calculation, we have Then we have 1, then the order-1 periodic solution is orbitally asymptotically stable.

3.2.3
Existence and stability of the order-1 periodic solution of system (3) Next, we mainly investigate the existence of the order-1 periodic solution of system As is shown in Figure 4, line l 4 intersects with the phase set N 2 = {(x, y) ∈ R 2 + |y = (1 − α 2 )h 2 } and impulsive set M 2 = {(x, y) ∈ R 2 + |y = h 2 } at points A 4 and B 4 , respectively. We denote the coordinates of points A 4 and B 4 as (x A 4 , (1 − α 2 )h 2 ) and (x B 4 , h 2 ). In the following, we prove the existence of order-1 periodic solution of system (3).
small, the trajectory from A 3 intersects with M 2 at point B 3 , and then jumps to A + 3 due to the impulsive effect. Since x A 3 is sufficiently small, then the the successor function of point A 3 satisfies that G( we only need to find another point A 4 ∈ N 2 such that G(A 4 ) < 0. Consider the intersection point A 4 of the dashed line l 4 and N 2 , the trajectory from A 4 intersects with M 2 at point B 4 , and then jumps to A + 4 due to the impulsive effect. Figure   4. By the continuity of successor function, there must exist a point A between A 3 and A 4 satisfying that G(A) = 0, i.e., there exist an order-1 periodic solution of system (3). Fig. 4 The existence of order-1 periodic solution of system (3).

Remark 1
The proof of the stability of order-1 periodic solution of system (3) is similar to the proof of theorem 2, so we omit it here.

Existence and stability of the order-2 periodic solution
In this subsection, we mainly prove the existence and stability of order-2 periodic solution of system (4). (4) Next, we will prove the existence of order-2 periodic solution of system (4). For convenience, we suppose that the isocline l 1 intersects with impulsive set M 2 and phase set N 1 at points G 1 (x G 1 , h 2 ) and G 2 (x G 2 , (1 + α 1 )h 1 ), respectively. Isocline l 2 and the line l 3 intersect with impulsive set M 1 at point D 3 (x D 3 , h 1 ) and point H(x H , h 1 ), respectively. l 5 and l 6 intersect with N 1 and N 2 at point

Existence of the order-2 periodic solution of system
, respectively. Besides, there must exist a point S(x S , (1 − α 2 )h 2 ) ∈ N 2 , such that the trajectory from point S intersects with impulsive set M 1 at point D 3 .
an order-2 periodic solution.
Proof Choose a point A 1 ∈ N 1 , where A 1 is very closed to y-axis, i.e., x A 1 is positive and small enough. The trajectory starting from point A 1 will intersect with the impulsive set M 2 at point B 1 (x B 1 , h 2 ). Based on the second impulsive function, In addition, x A 1 is positive and small enough, Figure 5(a)).
Choose another point A 2 ∈ N 1 , and A 2 is sufficiently closed to point F 2 and it is on the left side of point F 2 . Similar to the discussion above, the trajectory By the continuity of successor function, there is a point A between A 1 and A 2 satisfying G(A) = 0, then AB ∪ BC ∪ CD ∪ DA forms an order-2 periodic solution, as is shown in Figure 5(b).
It then follows that the stability of order-2 periodic solution of system (4).

Examples and numerical simulations
In this section, some numerical simulations will be presented to verify the the- and In order to get the order-1 periodic solution of system (6), we select control parameters h 1 = 0.15, α 1 = 0.25 and τ 1 = 0.02. We find that the trajectory starting from the initial value (x 0 , y 0 ) = (0.4, 0.5) will tend to be an order-1 periodic solution (see Figure 6). And Figure 6(d) is a partial enlargement of Figure   6(c). We keep other parameters unchanged but τ 1 = 0.2, the trajectory from the same initial value also tends to be an order-1 periodic solution (see Figure 7).
Comparing Figures 6(d) and 7(c), we find that system (6)  In order to get the order-1 periodic solution of system (7), we select control parameters: h 2 = 0.75, α 2 = 0.4 and τ 2 = 0.12. And we find that the trajectory starting from initial value (x 0 , y 0 ) = (0.2, 0.5) tends to be an order-1 periodic solution (see Figure 8). Keep other parameters unchanged, choosing the control parameter τ 2 = 0.05, and we find the trajectory from same initial value also tends to be an order-1 periodic solution (see Figure 9). Comparing Figures 9(c) and 10(c), we find that system (7) has different types of order-1 solutions with different values of τ 2 . Example 2 In order to illustrate the existence of order-2 periodic solution of system (4) in section 3.3, by using the same model parameters, we get the following specific In order to get the order-2 periodic solution of system (8) 3) tends to be an order-2 periodic solution (see Figure 10). Keep other parameters unchanged but τ 1 = 0.075, , and we find that the trajectory from the same initial value tends to be an order-2 periodic solution (see Figure 11).
From Figure 12, we can see that system (4) can exhibit order 1, 2, 3, and 4, periodic orbits with τ 2 ∈ (0, 0.9). More precisely, there exists an order-1 periodic orbit for  with the change of parameter τ 2 . And if we select τ 1 , α 1 or α 21 as bifurcation parameter, we have similar results, so omitted here.

Conclusions
In this paper, we develop and analyze a mathematical model for the management of overgrazing based on a guanaco-sheep competitive system with unilateral and bilateral control. And the case that the unique positive equilibrium E * 1 of the basic ODE model is unstable is mainly discussed. More precisely, the existence of order-1 and order-2 periodic solution is verified by the method of successor function.
And the stability of order-1 and order-2 periodic solution is proved by applying the geometry theory and analogue of Poincaré criterion.
Numerical results show that if the control parameter α 1 , τ 1 are given suitably, then system (2) exists an order-1 periodic solution. It means that guanaco and sheep species can be maintained periodically in certain level. And if we change the value of τ 1 and keep other parameters unchanged, then system (2) also has an order-1 periodic solution. Thus, system (2) exhibits different types of order-1 periodic solution with different values of τ 1 . There exists similar results for system (3). In addition, if the control parameters α 1 , τ 1 , α 2 , τ 2 are chosen suitably, then system (4) has an order-2 periodic solution. By drawing bifurcation diagram of periodic solution with τ 2 ∈ (0, 0.9), we find that system (4) exists high order periodic solution, such as order-3 and order-4 periodic solutions.
Our study shows that if unilateral control measure is adopted, then system (2) or (3) has order-1 periodic solution, it means that population size of these two species changes periodically. For system (2), this unilateral control strategy can only protect farmers' interests, but it may not alleviate desertification caused by overgrazing. For system (3), this unilateral control strategy can alleviate desertification caused by overgrazing, but it does not guarantee farmers' benefits.
Combining these two aspects, then bilateral control strategy is adopted, and the results show that the guanaco species can be protected and sheep species can be limited at a certain level. It means that the farmer's profits could be guaranteed and the desertification caused by overgrazing may be mitigated. Comparing the unilateral control and bilateral control strategies, our study encourages bilateral control rather than unilateral control for the management of sheep species.  (x, y)). And if (x, y) / ∈ M (x, y), the system develops under the regulation of f (x, y) = ( dx dt = P (x, y), dy dt = Q(x, y)), this part is similar to a continuous system.
2. For a certain point P 2 ∈ N 2 , the trajectory from P 2 intersects with M 2 at point Q 2 , then the impulsive mapping φ 2 maps it to the phase point Q + 2 . Q + 2 is called the order-1 successor point of P 2 , and G(P 1 ) = y Q + 1 − y P 1 is called the successor function of point P 2 . If Q + 2 coincides with P 2 , then the trajectory P 2 Q 2 ∪ Q 2 P 2 forms an order-1 periodic solution, as is shown in Figure 14(a).
3. For a certain point P 1 ∈ N 2 , the trajectory from P 2 intersects with M 1 at point Q 1 , then the impulsive mapping φ 1 maps it to the phase point Q + 1 ∈ N 1 . And, the trajectory from Q + 1 intersects with M 2 at point S 1 , then the impulsive mapping φ 2 maps it to the phase point S + 1 ∈ N 2 . S + 1 is called order-2 successor point of P 1 . If S + 1 coincides with P 1 , then the trajectory P 1 Q 1 ∪ Q 1 Q + 1 ∪ Q + 1 S 1 ∪ S 1 P 1 forms an order-2 periodic solution, as is shown in Figure 14(b). Besides, point S + 1 is called the order-2 successor point of P 1 , and G(P 1 ) = y S + 1 − y P 1 is called the order-2 successor function of point P 2 .