Global Dynamics of the Hastings-Powell System

This paper studies the problem of bounding a domain that contains all compact invariant sets of the Hastings-Powell system. The results were obtained using the first-order extremum conditions and the iterative theorem to a biologically meaningful model. As a result, we calculate the bounds given by a tetrahedron with excisions, described by several inequalities of the state variables and system parameters. Therefore, a region is identified where all the system dynamics are located, that is, its compact invariant sets: equilibrium points, periodic-homoclinic-heteroclinic orbits, and chaotic attractors. It was also possible to formulate a nonexistence condition of the compact invariant sets. Additionally, numerical simulations provide examples of the calculated boundaries for the chaotic attractors or periodic orbits. The results provide insights regarding the global dynamics of the system.


Introduction
Chaos is an interesting phenomenon that has been studied extensively over recent years.Moreover, it continues to be an active research topic, particularly in the area of analysis of nonlinear dynamic systems.Chaotic nature can be found in a huge variety of disciplines; interesting examples can be found in financial systems [1], cancer tumor growth models [2], power systems [3], weather models [4], electromechanical systems [5], and many more.For instance, complex dynamics are detected in several biological models, with many works published in recent years [6][7][8][9][10].In particular, this paper analyzes a chaotic biological model which describes a tritrophic food chain interaction within an ecosystem, commonly referred to as the Hastings-Powell (HP) system [11].Dynamic properties of this model have been studied in some contributions; for instance, see papers [12][13][14].
The main interest in bounding a dynamical system is due to the fact that such bounds can provide relevant information regarding its global dynamics.The bounds are domains given by inequalities, expressed in terms of the parameters and state variables of a system.Within these bounds, it is possible to find every compact invariant set of the system, which characterize its most interesting dynamics which include periodic, homoclinic, or heteroclinic orbits; equilibrium points, and chaotic behavior.Besides the purely theoretical interest of obtaining such bounds, the results may provide valuable insights regarding the evolution and future of the considered system.
The bounds in this work are obtained by the application of first-order extremum conditions, proposed initially by Krishchenko [15] and improved in recent years by Krishchenko and Starkov: they propose the iterative method of localization of compact invariant sets in [16], ellipsoidal estimates in [17], and localization for time-varyng systems in [18], for Hamiltonian systems [19], and they continue to work on the development of the general theory and applications of the method.One of the most interesting improvements of the method is the iterative theorem, which can be used to sharpen the bounds based on previous results of localization [16].This work includes some interesting results obtained by applying the iterative theorem (Theorem 3 described below).The localization method has been successfully applied recently on dynamic system in diverse areas such as permanentmagnet motor system [5], Lanford system [20], Lorenz system [21], the amplitude of a plasma instability [22], and Nonetheless, recently an interest has arisen to understand the global dynamics of a biologically meaningful model.Analyzing the HP system contributes with an alternative way to understanding its behavior.
Similar analytic methods (see [8,[23][24][25] and references therein) calculate the geometric approximation of attractors or special cases of periodic dynamics.We calculate a region that includes all compact invariant sets of the system.Our results can be used to analyze the global dynamics and combined with nonlinear analysis theories may provide important information for engineering applications [26,27].Regarding previous results obtained for the global dynamics analysis of the HP system [28][29][30], we can comment that our approach gives a novel perspective of the boundaries of the system.Besides, various essential properties of dynamics like nonchaoticity and nonexistence of repeating behavior may be analyzed with the help of this approach; see [31,32].
The paper is organized as follows.Some useful assertions and notations are provided in Section 2. Section 3 presents results concerning the bounds of the Hastings-Powell model.This section is divided into several subsections, which presents diverse applications of the localization method.Section 4 includes numerical simulation examples of the bounds with respect to some compact invariant sets.Finally Section 5 contains conclusions.

Preliminaries of Bounding
This section is devoted to describing the general theorems used to obtain the bounds of the system.This method is commonly referred to as the localization of compact invariant sets; however, in this paper we refer to it as the "localization method." The general localization method for a nonlinear system was described by Starkov and Krishchenko in [15,[20][21][22].This section reviews useful results from these works.Consider a nonlinear system where () is a  ∞ -differentiable vector field and  ∈ R  is the state vector.Let ℎ() be a  ∞ -differentiable function such that ℎ is not the first integral of (1).Take ℎ|  the restriction of ℎ on a set  ⊂ R  .(ℎ) denotes the set { ∈ R  |   ℎ() = 0}, where   ℎ() is a Lie derivative with respect to vector field () corresponding to the system (1).Now define Theorem 1.Each compact invariant set Γ of (1) is contained in the localization set [15] If the location of all compact invariant sets is considered inside the domain  ⊂ R  , then the localization set (ℎ) ∩  is valid, with (ℎ) defined in Theorem 1.Let  be a subset in R  .Proposition 2. If  ∩ (ℎ) = 0, then the system (1) has no compact invariant sets located in  [15].
A refinement of the localization set (ℎ) is realized with the help of the iterative theorem, which is stated as follows.

The Hastings-Powell System
The model of interest describes a tritrophic food chain.This system was obtained originally by Hogeweg and Hesper [33].
Later, Hastings and Powell [11] found chaotic dynamics, and is now known as the Hastings-Powell (HP) system in many papers from the field of ecology [12,34,35].The model is given by the following system of equations: with describing a type II functional response of both of the consumer species.Here,  is time,  is the top-predator,  is the prey and  is a predator that consumes the prey .The constant  0 is the intrinsic growth rate. 0 is the carrying capacity of species .Constants  −1 1 and  2 are conversion rates of prey to predator of species  and , respectively. 1 and  2 are constant death rates for species  and  respectively.Constants   and   , for  = 1, 2, parametrize the saturation functional response.  is the prey population level where the predation rate per unit prey is half its maximum value.
Model (5) includes a set of 10 parameters that imply a difficult analysis.As is presented in [11], the parameters can be reduced to 6 with the linear transformation given by 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 8.5 Then, it is possible to have a simpler nondimensional form written as where the new parameters are given by denotes the derivatives of  with respect to .
This model was used to describe the dynamics of classical fishpond management for the tilapia fish culture [12,36].These works consider   ,  ∈ {1, 2}, like the corresponding density of young tilapia (prey), developed tilapia (predator), and variable  3 describes the density of tucunare fish which is considered the top-predator population.Gomes et al. [12] and Varriale and Gomes [36] analyzed the model and found a rich variety of dynamics.Parameters  1 = 5,  2 = 0.1,  2 = 2,  1 = 0.4, and  2 = 0.01, for system (8), are used traditionally; see [11]. Figure 1 shows the named "tea-cup" chaotic attractor obtained for specific parameters of the system (8).From an engineering perspective, the main interest of the HP system is to study secondary effects of the exploitation of biological resources and the harvest in an ecosystem.
In the following section, the bounds for the system (8) are calculated with the application of first-order extremum conditions and the iterative theorem.Numerical simulations are performed to show the bounds of the tea-cup attractor and other interesting compact invariant sets.

Bounding the Hastings-Powell System
This section is devoted to calculating the bounds of system (8).The bounds are calculated with the help of coordinate functions, linear combinations, and rational functions.Therefore, a tetrahedral localization with excisions (cuts over its geometric body) is derived as a result.Additionally, some nonexistence conditions are derived.
4.1.Tetrahedral Localization.Firstly, a compact localization set in R 3 + will be derived.This means that the region of localization delimits all variables of the system.This localization set is given by the following theorem.
Proof.We have Therefore, the set (ℎ 1 ) is defined by taking  1 ≥  2 .Now by using the last expression for (ℎ 1 ) This means that and the localization set (ℎ 1 ) is given by (ℎ This localization means that every compact invariant set of the Hastings-Powell system (8) is bounded in the region given by (ℎ 1 ) ∩ R 3  + .This tetrahedral localization includes the restriction of  1 ≥  2 , which coincides with the biological feasible web-food existence condition ( 1 ≫  2 ) given by [11].The upper limit of the localization set given by ℎ 1,sup varies with resoect to  2 (the constant death rate of the toppredator), where the bound of  1 can be obtained from the localization set (ℎ 1 ), and with  2 = 1, this coincides with that given by [11],  1 ≤ 1 (three lines before Table 1 in [11], page 897), and the (ℎ 2 ) obtained latter in this paper.

Excisions by Means of Additional Linear Functions.
With (ℎ 1 ) we have obtained a compact localization; nevertheless, some additional excisions (cuts over the tetrahedron body) can be obtained with the help of coordinate functions and the application of the iterative theorem.The next theorem gives a localization set from a coordinate function.The first result is given by the next theorem Theorem 5.All compact invariant sets located in R 3  + are contained in the localization set Proof.If the localizing function and Thus, and the localization set is given by (ℎ 2 ) = {0 < ℎ 2 ≤ 1}.
Now, with the application of the iterative theorem, the next theorem can be stated.Theorem 6.All compact invariant sets located in R 3  + are contained in the localization set  3 := (ℎ 2 ) ∩  2,3 , with Proof.If the localizing function Hence solving (ℎ 3 ) respecting  1 and  2 , then So because  1 = ℎ 2 () ≤ 1, we have As a result, one can get The latter formula of  2,3 is an alternative solution of the upper bound for  1 +  2 derived from (ℎ 1 ).In this case, the bound of (ℎ 2 ) was obtained in a similar way by [11], but in this work we obtained a tighter bound by applying the iterative theorem.Here, no restrictions are imposed on the localization set.The variation of the bound ℎ sup 2,3 w.r.t. 1 (the constant death rate of the middle-predator) implies that the minimal value of the bound is ℎ 3,sup = 1 and, as expected, coincides with the bound  1 ≤ 1.Now, an additional bound is given by the following theorem.
This localization gives an additional excision of the localization set.

Rational Functions.
Here, rational functions are applied to sharpen the bounds; this is possible because of the iterative theorem.This improves previous bounds; however, there are restrictions on them.These restrictions are not fulfilled with the parameters given in the literature, but they are included for knowledge in future research.Then with the rational function ℎ 5 =  2  −1 1 it is possible to state the following theorem.

Theorem 8. Suppose that
Then, all compact invariant sets located in R 3 + are contained in the set Proof.The localizing function ℎ 5 =  2  −1 1 is applied respecting R 3  + , then Therefore, on (ℎ 5 ) the formula holds.The lower bound for the coefficient of ℎ 5 in the last formula with respect to the set  3 ∩ (ℎ 5 ) will be calculated with In a similar way, it is possible to state the next.
Theorem 9. Suppose that  3 > 0 and Then, all compact invariant sets located in R 3 + are contained in the set Proof.If the localizing function ℎ 6 =  2  −1 3 is used respecting R 3  + , then It follows from the latter formula that is valid on the set (ℎ 6 ).Indeed, since  4 > 0, then (ℎ 6 ) includes then the localization set  6 is obtained because of Theorem 3 (iterative theorem).Proof completed.

Nonexistence Conditions.
Non-existence conditions in dynamic systems involve the possibility to find some range of parameters where the system lacks compact invariant sets.Therefore, in the case of system (8), if it is possible to find two nonexistence conditions, which implies that if a compact invariant set exists for the system, then, must be located inside the coordinate planes.This conclusion can be obtained due to the main characteristic of the biological system regarding the positiveness and invariance of the first orthant R 3 + .For system (8), two non-existence conditions are calculated and are given by the next theorems.
Theorem 10.If  2 ≤  2  2 , then there are no compact invariant sets located in R 3 + .
Proof.If the function ℎ 7 =  3 is taken, then and can the following set be obtained Therefore, if  2 ≤  2  2 , (ℎ 7 ) = 0.The proof is complete because of the Proposition 2. Now, it is possible to find an additional nonexistence condition.
Theorem 11.If  1 ≤  1  1 , then there are no compact invariant sets located in R 3 + .
Proof.If the function , and the following set can be obtained: Therefore, if The proof is complete based on Proposition 2.

Condition of Survival of Predator Populations.
At this point, it is relevant to mention one noteworthy remark concerning dynamics of system (8).The necessary condition for the population  to survive is Indeed, if  = 2 then we have that The nonexistence conditions coincide with the condition of survival of predators.Here, it is important to mention that they share a similar interpretation but have specific differences between them.A nonexistence condition means that there are no compact invariant sets in R 3  + .The condition of survival for the predator means that if the condition is fulfilled, then the population of predator will be extinct.The combination of both results allows us to assert the attractivity of coordinate planes  2 = 0 and  3 = 0. Therefore, if one of the conditions is fulfilled, the dynamics of the system (8) are trapped inside a specific coordinate plane and a reduced order system is obtained.

Numerical Simulations
This section presents numerical simulations for a few special cases of parameters, given by Kuznetsov and Rinaldi [8]; the localization set (ℎ 1 ) ∩ (ℎ 2 ) ∩  2,3 is valid, and as expected, the localization set is different for each special case.This boundary depends on the parameters  1 and  2 , the constant death rates of the species, meaning that the localization set is described by the parameters of the system.Compared with numerical methods, localization surface is given by a mathematical function of the parameters and state variables of the system.
Inequalities for the existence of other localization sets are not satisfied for these values of parameters.In particular,  3 and  4 are negative.However, if then the localization by the set  5 exists.Besides, as it was indicated in [11], biologically reasonable food chains are realized when  1 ≫  2 , which entails  4 < 0.

Conclusions
This paper describes how to obtain a compact boundary containing all compact invariant sets of a biological system which is given by a three-species food model called the Hastings-Powell system.These bounds are given by a tetrahedron with linear excisions.The linear excisions are obtained by using coordinate and a rational localizing function.The bounds depend on the constant death rates for each species.
Our method requires a candidate with some function, and we have to use heuristic methods for its choice.Here we deal with biological systems that lie within the positive orthant R 3 + ; therefore, our choice of that function is atypical compared with traditional polynomial systems where quadratic-style functions are usually proposed.Here linear, coordinate, and rational functions give better results.Moreover, the iterative method for localization can be applied to sharpen the bounds.Two nonexistence conditions of compact invariant sets are described as well, and two conditions for the survival of predator populations are provided.
Numerical simulations are performed to show the bounds of some compact invariant sets, that is, periodic orbits and chaotic attractors.The bounds provide information regarding the global dynamics of the Hastings-Powell system.The nonexistence conditions allowed us to find an important contribution about the dynamics of the system combined with the conditions of predator survival, about attractivity of coordinate planes  2 = 0 and  3 = 0.Moreover, numerical simulations show a graphical interpretation of the mathematical work.

Figure 4 :
Figure4: Bounds for the system (8) and a case of periodic orbit.

3 Figure 6 :
Figure6: Bounds for the system (8) and a case of periodic orbit.

Figure 7 :
Figure 7: Bounds for the tea-cup chaotic attractor of the system (8).