Global Dynamics of the Hastings-Powell System

Abstract

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.

1. 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–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–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 permanent-magnet 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–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–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.

2. 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–22]. This section reviews useful results from these works.

Consider a nonlinear system 𝑥̇=𝑓(𝑥), (1) where 𝑓(𝑥) is a 𝐶^∞-differentiable vector field and 𝑥∈𝐑^𝑛 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 𝐵⊂𝐑^𝑛. 𝑆(ℎ) denotes the set {𝑥∈𝐑^𝑛∣𝐿_𝑓ℎ(𝑥)=0}, where 𝐿_𝑓ℎ(𝑥) is a Lie derivative with respect to vector field 𝑓(𝑥) corresponding to the system (1). Now define ℎ_inf∶=inf{ℎ(𝑥)∣𝑥∈𝑆(ℎ)}; ℎ_sup∶=sup{ℎ(𝑥)∣𝑥∈𝑆(ℎ)}.

Theorem 1. Each compact invariant set Γ of (1) is contained in the localization set [15] 𝐾(ℎ)={ℎ_inf≤ℎ(𝑥)≤ℎ_sup}.

If the location of all compact invariant sets is considered inside the domain 𝑈⊂𝐑^𝑛, then the localization set 𝐾(ℎ)∩𝑈 is valid, with 𝐾(ℎ) defined in Theorem 1. Let 𝑄 be a subset in 𝐑^𝑛.

Proposition 2. If 𝑄∩𝑆(ℎ)=∅, 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.

Theorem 3. Let ℎ_𝑚(𝑥),𝑚=0,1,2,…, be a sequence of functions from 𝐶^∞(𝐑^𝑛) [16]. Sets 𝐾₀=𝐾ℎ₀,𝐾_𝑚=𝐾_𝑚−1∩𝐾_{𝑚−1,𝑚},𝑚>0, (2) with 𝐾_{𝑚−1,𝑚}=𝑥∶ℎ_𝑚,inf≤ℎ_𝑚(𝑥)≤ℎ_𝑚,sup,ℎ_𝑚,sup=sup_{𝑆ℎ𝑚∩𝐾𝑚−1}ℎ_𝑚(𝑥),ℎ_𝑚,inf=inf_{𝑆ℎ𝑚∩𝐾𝑚−1}ℎ_𝑚(𝑥), (3) contain any compact invariant set of the system (1) and 𝐾₀⊇𝐾₁⊇⋯⊇𝐾_𝑚⊇⋯. (4)

3. 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: 𝑑𝑋𝑑𝑇=𝑅₀𝑋1−𝐾₀⁻¹𝑋−𝐶₁𝐹₁(𝑋)𝑌,𝑑𝑌𝑑𝑇=𝐹₁(𝑋)𝑌−𝐹₂(𝑌)𝑍−𝐷₁𝑌,𝑑𝑍𝑑𝑇=𝐶₂𝐹₂(𝑌)𝑍−𝐷₂𝑍, (5) with 𝐹_𝑖(𝑈)=𝐴_𝑖𝑈𝐵_𝑖+𝑈,𝑖=1,2, (6) 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 𝑅₀ is the intrinsic growth rate. 𝐾₀ is the carrying capacity of species 𝑋. Constants 𝐶₁⁻¹ and 𝐶₂ are conversion rates of prey to predator of species 𝑌 and 𝑍, respectively. 𝐷₁ and 𝐷₂ 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 𝑥₁=𝑋𝐾₀;𝑥₂=𝐶₁𝑌𝐾₀;𝑥₃=𝐶₁𝑍𝐶₂𝐾₀,𝑡=𝑅₀𝑇. (7) Then, it is possible to have a simpler nondimensional form written as 𝑥̇₁=𝑥₁1−𝑥₁−𝑎₁𝑥₁𝑥₂1+𝑏₁𝑥₁;𝑥̇₂=𝑎₁𝑥₁𝑥₂1+𝑏₁𝑥₁−𝑎₂𝑥₂𝑥₃1+𝑏₂𝑥₂−𝑑₁𝑥₂;𝑥̇₃=𝑎₂𝑥₂𝑥₃1+𝑏₂𝑥₂−𝑑₂𝑥₃, (8) where the new parameters are given by 𝑎₁=(𝐾₀𝐴₁)/(𝑅₀𝐵₁), 𝑏₁=𝐾₀/𝐵₁, 𝑎₂=(𝐶₂𝐴₂𝐾₀)/(𝐶₁𝑅₀𝐵₂), 𝑏₂=𝐾₀/(𝐶₁𝐵₂), 𝑑₁=𝐷₁/𝑅₀, 𝑑₂=𝐷₂/𝑅₀; finally 𝑥̇_𝑖,𝑖 with 𝑖∈{1,2,3} 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 𝑥₃ 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 𝑎₁=5, 𝑎₂=0.1, 𝑏₂=2, 𝑑₁=0.4, and 𝑑₂=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.

Figure 1 

Chaotic attractor for the Hastings-Powell system (8) (Figure 2 from [11], page 898).

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.

4. 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.

Since variables 𝑥_𝑖,𝑖∈{1,2,3}, in (8), have a biological interpretation we examine compact invariant sets that lie within the positive orthant 𝐑³₊∶={𝑥_𝑖>0;𝑖∈{1,2,3}}. In addition, all parameters of the model are supposed to be positive and 𝐑³₊ is an invariant set. In what follows, 𝑓 denotes the vector field of the system under consideration. Also, for the sake of simplicity of notations, 𝐾(ℎ)∶=𝐾(ℎ)∩𝐑³₊ and 𝑆(ℎ)∶=𝑆(ℎ)∩𝐑³₊.

4.1. Tetrahedral Localization

Firstly, a compact localization set in 𝐑³₊ 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.

Theorem 4. Let ℎ₁=𝑥₁+𝑥₂+𝑥₃; if 𝑑₁≥𝑑₂ then all compact invariant sets are located in the tetrahedron 𝐾ℎ₁=𝑥₁+𝑥₂+𝑥₃≤𝑑₂+1²4𝑑₂. (9)

Proof. We have 𝐿_𝑓ℎ₁=𝑥₁(1−𝑥₁)−𝑑₁𝑥₂−𝑑₂𝑥₃. Therefore, the set 𝑆(ℎ₁) is defined by 𝑆ℎ₁∶=𝑥₃=𝑥₁1−𝑥₁𝑑₂−𝑑₁𝑑₂𝑥₂, (10) taking 𝑑₁≥𝑑₂. Now by using the last expression for 𝑆(ℎ₁)ℎ_1𝑆(ℎ1)1≤−𝑑₂𝑥₁−𝑑₂+12²+𝑑₂+1²4𝑑₂. (11) This means that ℎ_1,sup≤𝑑₂+1²4𝑑₂, (12) and the localization set 𝐾(ℎ₁) is given by 𝐾(ℎ₁)={ℎ₁≤ℎ_1,sup}.

This localization means that every compact invariant set of the Hastings-Powell system (8) is bounded in the region given by 𝐾(ℎ₁)∩𝐑³₊. This tetrahedral localization includes the restriction of 𝑑₁≥𝑑₂, which coincides with the biological feasible web-food existence condition (𝑑₁≫𝑑₂) given by [11]. The upper limit of the localization set given by ℎ_1,sup varies with resoect to 𝑑₂ (the constant death rate of the top-predator), where the bound of 𝑥₁ can be obtained from the localization set 𝐾(ℎ₁), and with 𝑑₂=1, this coincides with that given by [11], 𝑥₁≤1 (three lines before Table 1 in [11], page 897), and the 𝐾(ℎ₂) obtained latter in this paper.

4.2. Excisions by Means of Additional Linear Functions

With 𝐾(ℎ₁) 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 𝐑³₊ are contained in the localization set 𝐾ℎ₂=0<𝑥₁≤1. (13)

Proof. If the localizing function ℎ₂=𝑥₁ is taken, then 𝐿_𝑓ℎ₂=𝑥₁(1−𝑥₁)−𝑎₁𝑥₁𝑥₂(1+𝑏₁𝑥₁)⁻¹, and 𝑆ℎ₂∶=1−𝑥₁−𝑎₁𝑥₂1+𝑏₁𝑥₁=0. (14) Thus, ℎ_2𝑆(ℎ2)𝑎=1−¹𝑥₂1+𝑏₁𝑥₁≤1, (15) and the localization set is given by 𝐾(ℎ₂)={0<ℎ₂≤1}.

Now, with the application of the iterative theorem, the next theorem can be stated.

Theorem 6. All compact invariant sets located in 𝐑³₊ are contained in the localization set 𝐾₃∶=𝐾(ℎ₂)∩𝐾_2,3, with 𝐾_2,3∶=𝑥₁+𝑥₂≤1+𝑑₁⁻¹. (16)

Proof. If the localizing function ℎ₃=𝑥₁+𝑥₂ is used, then 𝐿_𝑓ℎ₃=𝑥₁(1−𝑥₁)−𝑎₂𝑥₂𝑥₃(1+𝑏₂𝑥₂)⁻¹−𝑑₁𝑥₂ and 𝑆ℎ₃=𝑥₁−𝑥²₁=𝑥₂𝑎₂𝑥₃1+𝑏₂𝑥₂⁻¹+𝑑₁. (17) Hence solving 𝑆(ℎ₃) respecting 𝑥₁ and 𝑥₂, then 𝑆ℎ₃∩𝐾ℎ₂⊂𝑥₂𝑎₂𝑥₃1+𝑏₂𝑥₂⁻¹+𝑑₁⊂≤10<𝑥₂≤𝑑₁⁻¹. (18) So because 𝑥₁=ℎ₂(𝑥)≤1, we have ℎ_3,sup=sup_{𝑆ℎ𝑥∈3∩𝐾2}ℎ₃(𝑥)=sup_{𝑆ℎ3∩𝐾2}𝑥₁+𝑥₂≤1+𝑑₁⁻¹. (19) As a result, one can get 𝐾_2,3∶=𝑥∣ℎ₃(𝑥)≤ℎ_3,sup. (20)

The latter formula of 𝐾_2,3 is an alternative solution of the upper bound for 𝑥₁+𝑥₂ derived from 𝐾(ℎ₁). In this case, the bound of 𝐾(ℎ₂) 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 ℎ_sup2,3 w.r.t. 𝑑₁ (the constant death rate of the middle-predator) implies that the minimal value of the bound is ℎ_3,sup=1and, as expected, coincides with the bound 𝑥₁≤1. Now, an additional bound is given by the following theorem.

Theorem 7. If 𝑑₁≥𝑑₂ then all compact invariant sets are located in 𝐾ℎ₄=𝑥₁+𝑥₂−𝑥₃≤𝑑₁+1²4𝑑₁. (21)

Proof. Proposing the localizing function ℎ₄=𝑥₁+𝑥₂−𝑥₃, then 𝐿_𝑓ℎ₄=𝑥₁−𝑥²₁−𝑑₁𝑥₂+𝑑₂𝑥₃−2𝑎₂𝑥₂𝑥₃𝑏₂𝑥₂+1⁻¹. (22) Therefore by taking ℎ₄|_𝑆(ℎ4) by completing squares in a similar way to Theorem 4, taking 𝑑₁≥𝑑₂, then ℎ_4,sup≤𝑑₁+1²4𝑑₁, (23) and the localization set 𝐾(ℎ₄) is given by 𝐾(ℎ₄)={ℎ₄≤ℎ_4,sup}.

This localization gives an additional excision of the localization set.

4.3. 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 ℎ₅=𝑥₂𝑥₁⁻¹ it is possible to state the following theorem.

Theorem 8. Suppose that 𝜂₃𝑑∶=²¹𝑏₁+1−𝑎₁𝑑₁−𝑎₁𝑑₁1+𝑏₁>0. (24) Then, all compact invariant sets located in 𝐑³₊ are contained in the set 𝐾₅=𝑥₂≤𝑎₁𝑑₁+1𝑥₁𝑑₁𝜂₃∩𝐾₃. (25)

Proof. The localizing function ℎ₅=𝑥₂𝑥₁⁻¹ is applied respecting 𝐑³₊, then 𝐿_𝑓ℎ₅=𝑎₁𝑥₂1+𝑏₁𝑥₁−𝑎₂𝑥₂𝑥₃𝑥₁1+𝑏₂𝑥₂−𝑑₁𝑥₂𝑥₁+𝑎₁𝑥²₂𝑥₁1+𝑏₁𝑥₁−𝑥₂1−𝑥₁𝑥₁. (26) Therefore, on 𝑆(ℎ₅) the formula ℎ₅𝑑₁+1−𝑥₁−𝑎₁𝑥₂1+𝑏₁𝑥₁+𝑎₂𝑥₃1+𝑏₂𝑥₂=𝑎₁𝑥₂1+𝑏₁𝑥₁ (27) holds. The lower bound for the coefficient of ℎ₅ in the last formula with respect to the set 𝐾₃∩𝑆(ℎ₅) will be calculated with inf𝑑₁+1−𝑥₁−𝑎₁𝑥₂1+𝑏₁𝑥₁+𝑎₂𝑥₃1+𝑏₂𝑥₂≥𝑑₁+1+inf−𝑥₁−𝑎₁𝑥₂1+𝑏₁𝑥₁≥𝜂₃. (28)

In a similar way, it is possible to state the next.

Theorem 9. Suppose that 𝜂₃>0 and 𝜂₄∶=𝑑₂−𝑑₁−𝑎₂𝑏₂+𝑎₂𝑑₁𝑏₂𝑑₁+𝑏₂𝑑₁+𝑏₂>0. (29) Then, all compact invariant sets located in 𝐑³₊ are contained in the set 𝐾₆=𝑥₂≤𝑎₂𝑑₁+1𝑥₃𝑑₁+𝑏₂𝑑₁+𝑏₂𝜂₄∩𝐾₅. (30)

Proof. If the localizing function ℎ₆=𝑥₂𝑥₃⁻¹ is used respecting 𝐑³₊, then 𝐿_𝑓ℎ₆=𝑎₁𝑥₁𝑥₂𝑥₃1+𝑏₁𝑥₁−𝑎₂𝑥₂1+𝑏₂𝑥₂−𝑑₁𝑥₂𝑥₃−𝑎₂𝑥²₂𝑥₃1+𝑏₂𝑥₂+𝑑₂𝑥₂𝑥₃. (31) It follows from the latter formula that ℎ₆𝑑₂−𝑑₁+𝑎₁𝑥₁1+𝑏₁𝑥₁−𝑎₂𝑥₂1+𝑏₂𝑥₂=𝑎₂𝑥₂1+𝑏₂𝑥₂ (32) is valid on the set 𝑆(ℎ₆). Indeed, since 𝜂₄>0, then 𝑆(ℎ₆) includes ℎ₆𝑑₂−𝑑₁+𝑎₁𝑥₁1+𝑏₁𝑥₁−𝑎₂𝑥₂1+𝑏₂𝑥₂=𝑎₂𝑥₂1+𝑏₂𝑥₂; (33) then the localization set 𝐾₆ is obtained because of Theorem 3 (iterative theorem). Proof completed.

4.4. 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 𝐑³₊. For system (8), two non-existence conditions are calculated and are given by the next theorems.

Theorem 10. If 𝑎₂≤𝑏₂𝑑₂, then there are no compact invariant sets located in 𝐑³₊.

Proof. If the function ℎ₇=𝑥₃ is taken, then 𝐿_𝑓ℎ₇=𝑎₂𝑥₂𝑥₃1+𝑏₂𝑥₂−𝑑₂𝑥₃, (34) and can the following set be obtained 𝑆ℎ₇=𝑎₂−𝑏₂𝑑₂𝑥₂=𝑑₂. (35) Therefore, if 𝑎₂≤𝑏₂𝑑₂, 𝑆(ℎ₇)=∅. The proof is complete because of the Proposition 2.

Now, it is possible to find an additional nonexistence condition.

Theorem 11. If 𝑎₁≤𝑏₁𝑑₁, then there are no compact invariant sets located in 𝐑³₊.

Proof. If the function ℎ₈=𝑥₂+𝑥₃ is taken, then 𝐿_𝑓ℎ₈=𝑎₁𝑥₁𝑥₂(1+𝑏₁𝑥₁)⁻¹−𝑑₁𝑥₂−𝑑₂𝑥₃, and the following set can be obtained: 𝑆ℎ₈=𝑎₁−𝑑₁𝑏₁𝑥₁𝑥₂=𝑑₁𝑥₂+𝑑₂𝑥₃1+𝑏₁𝑥₁. (36) Therefore, if 𝑎₁≤𝑏₁𝑑₁, 𝑆(ℎ₈)=∅. The proof is complete based on Proposition 2.

4.5. 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 𝑎_𝑗−1≥𝑏_𝑗−1𝑑_𝑗−1, 𝑗=2,3. Indeed, if 𝑗=2 then we have that 𝑥̇₂≤𝑥₂−𝑑₁+𝑎₁𝑥₁⁻¹+𝑏₁≤𝑥₂−𝑑₁+𝑎₁𝑏₁, (37) which entails lim_𝑡→∞𝑥₂(𝑡)=0 in case 𝑎₁<𝑏₁𝑑₁. The case 𝑗=3 can be considered similarly and lim_𝑡→∞𝑥₃(𝑡)=0 in case 𝑎₂<𝑏₂𝑑₂.

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 𝐑³₊. 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 𝑥₂=0 and 𝑥₃=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.

5. Numerical Simulations

This section presents numerical simulations for a few special cases of parameters, given by Kuznetsov and Rinaldi [8]; the localization set 𝐾(ℎ₁)∩𝐾(ℎ₂)∩𝐾_2,3 is valid, and as expected, the localization set is different for each special case. This boundary depends on the parameters 𝑑₁ and 𝑑₂, 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, 𝜂₃ and 𝜂₄ are negative. However, if 𝑏₁>𝑎₁𝑑₁+1𝑑²₁−1, (38) then the localization by the set 𝐾₅ exists. Besides, as it was indicated in [11], biologically reasonable food chains are realized when 𝑑₁≫𝑑₂, which entails 𝜂₄<0.

5.1. Periodic Orbits

5.1.1. Periodic Orbit [8]

With the parameters 𝑎₁=1.81; 𝑎₂=0.181; 𝑏₁=4.5; 𝑏₂=0.45; 𝑑₁=0.16; 𝑑₂=0.08; [8] it is a periodic orbit behavior is obtained. The bounds are given by the set 𝐾(ℎ₁)∩𝐾(ℎ₂)∩𝐾(ℎ₄); in this case 𝐾_2,3do not cut the tetrahedron 𝐾(ℎ₁). The periodic orbit and the bounds for these parameters are shown in Figure 2, where 𝐾(ℎ₁)={𝑥₁+𝑥₂+𝑥₃≤3.65}, and 𝐾(ℎ₂)={0<𝑥₁≤1}, and 𝐾(ℎ₄)={𝑥₁+𝑥₂−𝑥₃≤2.102}.

Figure 2 

Bounds for the system (8) and a periodic orbit.

5.1.2. Two Periodic Orbits [8]

With the set of parameters 𝑎₁=8; 𝑎₂=2.88; 𝑏₁=6.66; 𝑏₂=2.4; 𝑑₁=0.87; 𝑑₂=0.25; [8] is possible to obtain different periodic orbits for the system (8).Figure 3 shows two periodic orbits obtained for this set of parameters and the boundary surfaces. In this case 𝐾_2,3do not cut the tetrahedron 𝐾(ℎ₁). Here 𝐾(ℎ₁)={𝑥₁+𝑥₂+𝑥₃≤1.56};𝐾(ℎ₂)={0<𝑥₁≤1};𝐾(ℎ₄)={𝑥₁+𝑥₂−𝑥₃≤1.005}.

Figure 3 

Bounds for of the system (8) and two periodic orbits.

5.1.3. Periodic Orbit [8]

For parameters 𝑎₁=0.22; 𝑎₂=0.011; 𝑏₁=5.0; 𝑏₂=0.25; 𝑑₁=0.02; 𝑑₂=0.02; [8], an interesting behavior is given by a periodic orbit of the system (8). Figure 4 shows a periodic orbit obtained for this set of parameters and the bounds given by 𝐾(ℎ₁)∩𝐾(ℎ₂). In this case 𝐾_2,3and 𝐾(ℎ₄) do not cut the tetrahedron 𝐾(ℎ₁).Here 𝐾(ℎ₁)={𝑥₁+𝑥₂+𝑥₃≤13.05};𝐾(ℎ₂)={0<𝑥₁≤1}.

Figure 4 

Bounds for the system (8) and a case of periodic orbit.

5.1.4. Periodic Orbit [8]

Taking the set of parameters 𝑎₁=1.75, 𝑎₂=0.35, 𝑏₁=2.62, 𝑏₂=1.75, and 𝑑₁=0.1, 𝑑₂=0.04 [8], it is possible to obtain another periodic orbit. Similar to the previous localizations, 𝐾_2,3do not cut the tetrahedron 𝐾(ℎ₁). Figure 5 shows the boundary given by 𝐾(ℎ₁)∩𝐾(ℎ₂)∩𝐾(ℎ₄), and the periodic orbit for these parameters. Here 𝐾(ℎ₁)={𝑥₁+𝑥₂+𝑥₃≤6.76}; 𝐾(ℎ₂)={0<𝑥₁≤1}; and 𝐾(ℎ₄)={𝑥₁+𝑥₂−𝑥₃≤3.025}.

Figure 5 

Bounds for the system (8) and a case of periodic orbit.

5.1.5. Periodic Orbit [8]

Consider 𝑎₁=2.5, 𝑎₂=0.25, 𝑏₁=3, 𝑏₂=2, 𝑑₁=0.2, and 𝑑₂=0.025 [8] and again 𝐾_2,3do not cut the tetrahedron 𝐾(ℎ₁). Figure 6 shows the boundary set 𝐾(ℎ₁)∩𝐾(ℎ₂)∩𝐾(ℎ₄) and this periodic orbit, localization sets are given by 𝐾(ℎ₁)={𝑥₁+𝑥₂+𝑥₃≤10.506};𝐾(ℎ₂)={0<𝑥₁≤1}, and 𝐾(ℎ₄)={𝑥₁+𝑥₂−𝑥₃≤1.8}.

Figure 6 

Bounds for the system (8) and a case of periodic orbit.

5.2. Chaotic Attractor

Parameters 𝑎₁=5, 𝑎₂=0.1, 𝑏₂=2, 𝑑₁=0.4, 𝑑₂=0.01, for system (8), are used traditionally, see [7, 8, 11, 12, 36]. In this case, the bifurcation parameter 𝑏₁ is varied within the segment [2.0,6.2]. A chaotic attractor can be observed at 𝑏₁=2.75; see in [12]. Figure 7 provides a localization of this attractor by the set 𝐾(ℎ₁)∩𝐾(ℎ₂)∩𝐾_2,3∩𝐾(ℎ₄), with 𝐾(ℎ₁)={𝑥₁+𝑥₂+𝑥₃≤25.503}; 𝐾(ℎ₂)={0<𝑥₁≤1}; 𝐾_2,3={𝑥₁+𝑥₂≤3.5};𝐾(ℎ₄)={𝑥₁+𝑥₂−𝑥₃≤1.225}.

Figure 7 

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

The upper side of the plot in Figure 7 is a piece of the face of the tetrahedral localization 𝐾(ℎ₁); the latter has three excisions corresponding to the set 𝐾(ℎ₂)∩𝐾_2,3∩𝐾(ℎ₄).

6. 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 𝐑³₊; 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 𝑥₂=0 and 𝑥₃=0. Moreover, numerical simulations show a graphical interpretation of the mathematical work.

Conflict of Interests

The author declare that there is no conflict of interests.

Acknowledgments

This work was supported by DGEST Project TIJ-ING-2012-111 and the PROMEP Research Group ITTIJ-CA-6. The author would like to acknowledge Dr. Leonardo Trujillo and Dr. Konstantin E. Starkov for their contributions in various aspects of the preparation in this work. The author is pleased to acknowledge the helpful and insightful comments of the reviewers. A part of this work was presented in the IFAC CHAOS 2009 congress [37].