Bifurcation analysis and chaos control in discrete-time system of three competing species

In this paper, we investigate the complex dynamics of three-dimensional Ricker-type discrete-time competition model. We discuss the existence and uniqueness, and find parametric conditions for local asymptotic stability of positive fixed point of this model. It is also proved that the system undergoes Neimark–Sacker (NS) and period-doubling bifurcation (PDB) at certain parametric values for positive fixed point with the help of an explicit criterion for NS and PDB. The system shows chaotic dynamics at increasing values of bifurcation parameter. To control the chaos, we apply the hybrid control methodology. Finally, numerical simulations are provided to illustrate the theoretical discussions. These results of numerical simulations show chaotic long-term behavior over a broad range of parameters. The computation of the maximum Lyapunov exponents confirms the presence of chaotic behavior in the model.


Introduction
Difference and differential equations have been used to study a wide range of population models. Discrete-time systems may have favorable properties for simulations and mathematical analysis [1,2]. On the other hand, continuous systems are more realistic for ecological applications in which there is a strong overlap between generations (such as in microbial ecosystems). When the population remains constant over a generation or remains small over a number of generations, it would seem that the dynamics of the population is best described by discrete-time models [3,4]. Difference equations are used in modeling the interactions of species with nonoverlapping generations which are more suitable to study the behavior of population models [5][6][7][8]. Recently, discrete-time models are frequently implemented to study the complex and chaotic behavior of dynamical systems. For more detail, we refer the interested readers to [9][10][11][12][13]. Due to extremely complex dynamics, it may be very difficult to predict the comprehensive qualitative behavior even for single species [14,15]. In this article, we discuss the dynamics of a discrete-time system of three competing species, which is analogous to the May-Leonard differential equation model and is given as follows [16]: x n+1 = x n e r (1−x n −ay n −bz n ) , y n+1 = y n e r (1−y n −az n −bx n ) , z n+1 = z n e r (1−z n −ax n −by n ) . (1) In this symmetric model, three species have the same intrinsic growth rate r > 0. Moreover, we assume that the competition coefficients a and b for the discrete system (1) satisfy the condition that 0 < a < 1 < b. This assumption guarantees the competition among the three species. The general form of symmetric model (1) was proposed by Hofbauer et al. [17]. The aforementioned System (1) models the scramble competition among the competing species. The exponential terms model the effects of population density on vital rates. For example, the exponential term of the first equation of (1), that is, e r (1−x n −ay n −bz n ) , denotes the per unit time production of species x n in the presence of other species. Competition comes from the fact that the increase of other species will decrease one species production. When the intrinsic growth rate r is high, then System (1) shows more complex behavior due to the complicated nature of the Ricker equations. To help us understand the biological meaning of the competition parameters a and b of System (1), we consider B = b − 1 to be the dominant factor and A = 1 − a to be the sub-dominant factor with respect to all three species. For example, in the first equation of (1), B = b − 1 > 0 is the dominant factor of z n with respect to x n and A = 1 − a > 0 is the sub-dominant factor of y n with respect to x n . Therefore, when the dominant factor is less than the sub-dominant factor, none of the three species are dominant and the three species approach the fixed point. On the other hand, if the dominant factor is greater than the sub-dominant factor, the three species take turns becoming the dominant species. Our purpose is to study the local asymptotic stability, and NS and PDB of unique positive fixed point of System (1). To control the chaos due to emergence of PDB and NS bifurcations, hybrid control method is implemented to System (1).

Existence of positive fixed point
In this section, we study the existence and uniqueness of interior fixed point of System (1). First, we consider the possible fixed points of System (1), which can easily obtained by solving the equations of given system.

Linearized stability of unique positive fixed point
Let E * = (x * , y * , z * ) be a positive fixed point of System (1). Then the Jacobian matrix for (1) evaluated at this point is given by The characteristic polynomial of J (x * , y * , z * ) is given by: where To study the linear stability analysis of unique positive fixed point of System (1), we have the following Theorem 3.1 [18].

Theorem 3.1 [Jury condition]
Consider the third-degree polynomial equation: where C 1 , C 2 , and C 3 are real numbers. Then, necessary and sufficient conditions that all the roots of Eq. (6) lie in an open disk|λ| < 1 are:

Theorem 3.2 The unique positive fixed point E
is locally asymptotically stable for 0 < a < 1 < b if the following conditions are satisfied: where Proof The Jacobian matrix for System (1) evaluated at 1 1+a+b , 1 1+a+b , 1 1+a+b is given by: The characteristic polynomial for aforementioned matrix is computed as follows: where Now, applying Theorem 3.1, the positive fixed point 1 1+a+b , 1 1+a+b , 1 1+a+b is locally asymptotically stable for 0 < a < 1 < b, if the following conditions are satisfied:

Bifurcation analysis
In this section, we investigate the parametric conditions for the existence of NS and PDB at the positive fixed point (x * , y * , z * ) of System (1). Recently, many authors have discussed similar type of bifurcation analysis for discrete-time dynamical systems. For more detail on bifurcation analysis, we refer the interested readers to [19][20][21]. In dynamical systems, when a particular parameter passes through its critical value, various types of bifurcations emerge from its fixed point. Many dynamical properties of a system can be discussed owing to emergence of NS and PDB. Bifurcation usually occurs when the stability of a fixed point changes, i.e., qualitative properties of a dynamical system change.
We discuss NS and PDB for the positive fixed point 1 1+a+b , 1 1+a+b , 1 1+a+b of System (1) using an explicit criterion and taking r as a bifurcation parameter. Due to emergence of Neimark-Sacker bifurcation, closed invariant circles are produced. Equivalently, one can find some isolated orbits of periodic behavior along with trajectories that cover the invariant circle densely [22]. The bifurcation can be super-critical or sub-critical, resulting in a stable or unstable closed invariant curve, respectively. A period-doubling bifurcation in a discrete dynamical system is a bifurcation in which a slight change in a parameter value leads to the system switching to a new behavior with twice the period of the original system.

Neimark-Sacker bifurcation
To study the NSB in System (1), we need the following explicit criterion of Hopf bifurcation [23,24].

Lemma 4.1
Consider an n-dimensional discrete dynamical system X k+1 = f r (X k ), where r ∈ IR is a bifurcation parameter. Let X * be a fixed point of f r and the characteristic polynomial for Jacobian matrix where c i = c i (r, u), i = 1, 2, 3, . . . , n and u is control parameter or another parameter to be determined. Let Moreover, the following conditions hold: , when n is even or odd, respectively.
H2 Transversality condition: u). Then, Neimark-Sacker bifurcation occurs at r o . The following result shows that System (1) undergoes NSB if we take r as bifurcation parameter.
Proof According to Lemma (4.1), for three-dimensional system (n = 3), we have the characteristic polynomial (9) of System (1) evaluated at its unique positive equilibrium point (x * , y * , z * ). Then we obtain the following equalities and inequalities:

Period-doubling bifurcation
An explicit critical criterion for the existence of PDB is proposed for higher dimensional discrete-time systems [25].
Consider an n-dimensional discrete dynamical system X k+1 = f r (X k ), where X k+1 , X k ∈ IR n are the state vectors, k is the iterative index, and r ∈ IR is a bifurcation parameter. Let X * be a fixed point of f r and the characteristic polynomial for Jacobian matrix J (X * ) = (c i j ) n×n of n-dimensional map f r (X k ) is given by: where c i = c i (r ), i = 1, 2, . . . n. Consider a series of determinants: ± 0 (r ) = 1, ± 1 (r ), . . . , ± n (r ), which are defined as ± i (r ) = det (M 1 ± M 2 ), i = 1, 2, 3, . . . , n, where M 1 and M 2 are given in (12). Lemma 4.2 Assume that f r has a fixed point X * . Then a PDB takes place at r = r o iff the following conditions (H1) and (H2) are satisfied: (H1) Eigenvalue assignment: the following inequalities or equalities hold: . . , 1 (or 2); when n is odd or even; respectively; (H2)Transversality condition: where c i stands for first derivative of c(r ) with respect to r at r = r o .

Theorem 4.2
The unique positive fixed point (x * , y * , z * ) of System (1) undergoes PDB for 0 < a < 1 < b if the following conditions hold: where C 1 , C 2 , and C 3 are given in (10).
Proof According to Lemma (4.2), for three-dimensional system (n = 3), we have the characteristic polynomial (9) of System (1) evaluated at its unique positive equilibrium point (x * , y * , z * ). Then we obtain the following equalities and inequalities:

Chaos control
In dynamical systems, it is desired that the system be optimized with respect to some performance criterion and chaos be avoided. In recent times, controlling chaos in discrete-time systems is a topic of great interest for many researchers, and practical methods can be used in many fields such as communications, physics laboratories, cardiology, and turbulence [26]. In discrete-time models, chaos control can be obtained using various methods, some of them are the state feedback method, pole-placement technique, and hybrid control method [27][28][29][30]. Among all these chaos control methods, hybrid control technique is very simple to implement. Moreover, this method can be implemented for controlling chaos under the influence of both NS bifurcation and PDB.
To control the chaos in System (1), we apply the hybrid control feedback methodology [31][32][33][34]. Assume that System (1) undergoes NS bifurcation or PDB at unstable fixed point (x * , y * , z * ). Then the corresponding controlled system can be written as: where 0 < ρ < 1 denotes the control parameter, and in (18), the controlled strategy is a combination of both parameter perturbation and feedback control. Chaos for controlled system (18) can be delayed, advanced, or even completely eliminated by suitable choice of controlled parameter ρ. The Jacobian matrix of controlled system (18) evaluated at E * = ( 1 1+a+b , 1 1+a+b , 1 1+a+b ) is given by: The unique positive fixed point (x * , y * , z * ) of the controlled system (18) is locally asymptotically stable if roots of the characteristic polynomial of (19) lie in an open unit disk. The roots of (20) are λ 1 = 0.7099999986 and λ 2,3 = 0.9949999995 ± 0.09526279388i with λ 2,3 = 1.

Numerical simulations and discussions
According to Theorem 4.1, the conditions of NSB are satisfied near the fixed point (0.35, 0.35, 0.35) at the critical value of bifurcation parameter r = 0.29. Bifurcation diagrams and corresponding maximum Lyapunov exponent are shown in Fig. 1. Figure 1a-c shows that all three species undergo NSB. MLE confirms the existence of the chaotic sets (Fig. 1d). The phase portraits corresponding to Fig. 1 for different values of r are shown in Fig. 2, also the plots of x n , y n , and z n are shown in Fig. 3. Figure 3 shows that the unique positive fixed point (0.35, 0.35, 0.35) is unstable at the critical value of bifurcation parameter r = 0.29. The roots of (22) are λ 1 = −1 and λ 2,3 = 0.6666666675 ± 0.03849001803i with λ 2,3 = 1. Moreover, we have According to Theorem 4.2, the conditions of PDB are satisfied near the unique positive fixed point (0.444444, 0.444444, 0.444444) at the critical value of bifurcation parameter r = 2. Bifurcation diagrams and corresponding maximum Lyapunov exponent are shown in Fig. 4. Figure 4a-c shows that all three species undergo PDB. MLE confirms the existence of the chaotic sets (Fig. 4d). The phase portraits corresponding to Fig. 4 of System (1) for different values of r are shown in Fig. 5; also the plots of x n , y n and z n are shown in Fig. 6. Figure 6 shows that the fixed point  (24). b Plot of y n for system (24). c Plot of z n for system (24). d Phase portrait for system (24) period of the original system. Often, the fixed point loses its stability and a stable cycle of period 2 appears. Continued parameter changes may result in a cascade of PDB and the onset of chaos.  (1) undergoes NSB. Figure 2c shows that closed invariant circles appear at r = 0.35 enclosing this unstable fixed point.

Conclusion
We study the dynamics of three-dimensional discrete-time model of three competing species. We find the results about the existence, uniqueness, and conditions for local asymptotic stability of positive fixed point using the Jury condition. To confirm the complex and chaotic behavior in System (1), existence of Neimark-Sacker and period-doubling bifurcations for the positive fixed point is proved mathematically, as well as numerical simulations are provided. The computation of the maximum Lyapunov exponents confirms the presence of chaotic behavior in the model. Chaotic behavior is successfully controlled using the hybrid control strategy. From our numerical investigation, it is clear that the hybrid control method based on feedback control and parameter perturbation works effectively. To study PDB and NS bifurcation, intrinsic growth rate r is taken as bifurcation parameter. From our investigation, it is clear that the intrinsic growth rate r for species has strong stability effects or vice versa. Moreover, it is observed that System (1) exhibits chaotic trajectories, and intrinsic growth rate parameter r leading to extinction of three species, or coexistence of one or three species. This gives a high sensitivity of biodiversity with respect to parameter variations in the chaotic regions. Numerical simulations show that the feedback control law can suppress chaos to unstable equilibrium points successfully.