Abstract
This study examines the complexity of a discrete-time predator-prey system with ratio-dependent functional response. We establish algebraically the conditions for existence of fixed points and their stability. We show that under some parametric conditions the system passes through a bifurcation (flip or Neimark-Sacker). Numerical simulations are presented not only to justify theoretical results but also to exhibit new complex behaviors which include phase portraits, orbits of periods 9, 19, and 26, invariant closed circle, and attracting chaotic sets. Moreover, we measure numerically the Lyapunov exponents and fractal dimension to confirm the chaotic dynamics of the system. Finally, a state feedback control method is applied to control chaos which exists in the system.
1. Introduction
The interaction between predator and prey is one of the most studied topics in ecology and mathematical biology. The Lotka-Volterra model [1, 2] has received more attention from mathematician and ecologist. After them, many authors qualitatively analysed other types of predator-prey models in which response function of predator depends on prey densities only. But a number of respectable researchers [3–5] have claimed that in some environments (when predator needs to search and share foods), the response function of predator may depend on the ratio of prey to predator abundance. For detailed results in predator-prey systems with ratio-dependent response function, we refer to [6–9]. However, when the size of populations is small, the discretization of predator-prey systems is more suitable compared to continuous ones to understand unpredictable dynamic behaviors which exist in the system. Lots of empirical and theoretical works [10–20] have discussed the bifurcations and chaos phenomenon by using numerical simulations or center manifold theory and bifurcation theory.
In recent times, there are a few number of articles discussing the dynamics of ratio-dependent discrete-time predator-prey systems [21–23]. For discrete ratio-dependent predator-prey systems, it is shown that positive equilibrium is globally asymptotically stable [21], the strong and the weak Allee effects are investigated in [22], and it is discussed that periodic solutions exist in [23].
A ratio-dependent predator-prey system takes the form:where and stand for densities of prey and predator, respectively, is predator’s death rate, and is prey’s specific growth rate. is ratio-dependent response function of predator and are rate of capturing and conversion, respectively.
If and ), then system (1) takes in the form [6]:where , , , and are all positive constants. By scaling the variables and parameters as , , and , we write system (2) as with , , and being positive constants. Forward Euler scheme is applied to system (3) to obtain the following system:where is the step size. The motivation of this work is to study systematically the dynamical properties of (4) in detail which includes phase portraits, bifurcation in orbits of periods 9, 19, and 26, chaotic sets, Lyapunov exponents, and fractal dimension.
This paper is organized as follows. Section 2 deals with the existence condition for fixed points of (4) and their stability criterion. In Section 3, we prove that system (4) undergoes a flip bifurcation and a NS bifurcation in the interior of when changes its value in a small neighborhood of a specific parametric space. In Section 4, we provide numerical simulations for one or more control parameters to validate analytical results. In Section 5, we use the method of feedback control to stabilize chaos at the state of unstable trajectories. Finally, a brief discussion is carried out in Section 6.
2. Fixed Points: Existence and Their Stability
The fixed points of (4) are the solution of the following equations:
A simple algebraic computation yields the following result.
Lemma 1. System (4) has two fixed points:
(i) the predator free fixed point , which always exists;
(ii) the interior fixed point , where and This fixed point exists if with
To show the region in the space for which positive fixed point of system (4) exists, we plot in Figure 1. We see that lies inside but not in regions (I)-(II).
Next, we investigate stability of (4) at fixed points. The Jacobian matrix of system (4) at fixed point iswhere
Then the characteristic equation of (6) iswhereIt is noted that the magnitude of eigenvalues of Jacobian matrix evaluated at fixed point determines the local stability of that fixed point. The following propositions represent the stability conditions of fixed points by Jury’s criterion [24].
Proposition 2. For all permissible parameters values, there are four different topological types of and it is a(i)sink if and ;(ii)source if and ( and );(iii)nonhyperbolic if or ( and );(iv)saddle if otherwise.
Condition (iii) can be written as or
It follows that if parameters belong to , then two eigenvalues of are and . Thus a flip bifurcation may occur if parameters belong to . This observation illustrates that predator goes to extinct and prey passes through flip bifurcation leading to chaos as bifurcation parameter varies.
Proposition 3. If with , then positive fixed point of (4) exists and it is a (i)sink if either of the following conditions holds:(i.1) and ;(i.2) and ;(ii)source if either of the following conditions holds:(ii.1) and ;(ii.2) and ;(iii)nonhyperbolic if either of the following conditions holds:(iii.1) and ;(iii.2) and ;(iv)saddle if otherwise, where
From Proposition 3, it follows that if term (iii.1) holds, two eigenvalues of are and . We rewrite the term (iii.1) as follows: or Therefore, a flip bifurcation may emerge from fixed point if parameters belong to or .
Also two eigenvalues of are complex conjugate having magnitude one if (iii.2) holds. We rewrite the term (iii.2) as follows:
and if the parameters lie in , there exists a NS bifurcation emerging from .
3. Analysis of Bifurcation
In this section, attention is paid to recapitulate bifurcations (flip and/or Neimark-Sacker) of system (4) around fixed points using the theories in center manifold and of bifurcation [25–28]. The parameter is chosen as a bifurcation parameter.
3.1. Flip Bifurcation at Fixed Point
We consider system (4) at fixed point and take parameter arbitrarily from (similarly for the case of ).
Let , or , for the case of Then the eigenvalues of positive fixed point are The condition holds if
Using the transformation , and writing , we shift the fixed point of system (4) to the origin. After Taylor expansion, system (4) reduces towhere and
Therefore, we get the following symmetric multilinear vector functions of : where
Let be two eigenvectors (left and right) of associated with eigenvalue , respectively. Then and . After some algebraic calculation, we get
We set , where Then by the standard scalar product in defined by , it is obvious to see that . The direction of the flip bifurcation can be obtained by the sign of , the coefficient of critical normal form [25], and is computed via
The above discussion leads to the following result.
Theorem 4. If (17) holds and varies around the set , then system (4) passes through a flip bifurcation at fixed point . Moreover, the period-2 orbits that bifurcate from are stable (resp., unstable) if (resp., ).
3.2. Flip and Neimark-Sacker Bifurcation at Fixed Point
We first discuss flip bifurcation of system (4) at fixed point . We take parameter arbitrarily from (similarly for the case of ); then it follows that ; that is,Letor
Then the eigenvalues of positive fixed point areThe condition leads to
Using the transformation , and writing , we shift the fixed point of system (4) to the origin. After Taylor expansion, system (4) reduces to
where and
It follows thatand .
Therefore, we get the following symmetric multilinear vector functions of : and
Let be two eigenvectors (left and right) of associated with eigenvalue , respectively. Then and . After some algebraic calculation, we get
We set , where
Then by the standard scalar product in defined by , it is obvious to see that . The direction of the flip bifurcation can be obtained by the sign of , the coefficient of critical normal form [25], and is computed via
The above discussion leads to the following result.
Theorem 5. If (29) holds, , and varies around the set , then system (4) passes through a flip bifurcation at fixed point . Moreover, the period-2 orbits that bifurcate from are stable (resp., unstable) if (resp., ).
We next discuss a Neimark-Sacker bifurcation of system (4).
We first take parameter arbitrarily from and consider system (4) at fixed point .
Since the parameters belong to , the complex eigenvalues of (6) arewhere the condition leads to . Letthen and
Moreover, if , then which leads to
Let be eigenvectors of and corresponding to the eigenvalues and such that
Then by algebraic calculation, we obtain
We set , where
Then it is obvious that , where for . We can decompose vector as , for close to and . Obviously, . Thus, we obtain the following transformed form of system (30) for near : where with and is a smooth complex-valued function. Taylor expression of with respect to yields
According to multilinear symmetric vector functions, we can express the coefficients as follows:
We calculate the coefficient by the following formula to determine the direction in which invariant closed curve appears. where .
It is clear that two conditions (41) and (42) known as transversal and nondegenerate hold well for system (4). Therefore, we have the following result.
Theorem 6. If and changes its value around , then system (4) passes through a Neimark-Sacker bifurcation at fixed point . Moreover, if (resp., ), then there exists an attracting (resp., repelling) invariant closed curve which bifurcates from .
4. Numerical Simulations
Here, diagrams for bifurcation, phase portraits, Lyapunov exponents, and fractal dimension of system (4) will be presented by performing numerical simulation to validate our analytical results. The parameter is chosen as a bifurcation parameter to see the effect of ratio-dependent functional response in the system. We restrict our attention in the case of NS bifurcation only.
Now we consider bifurcation parameters in the following cases.
For case (i): by taking parameters as in Table 1, we see that, at the fixed point , a flip bifurcation occurs at and it shows . This confirms Theorem 4.
Figure 2(a) displays the bifurcation diagram of system (4). This illustrates that stability of fixed point happens for , bifurcation occurs at , and there appears period doubling bifurcation if . Figure 2(b) results in the maximum Lyapunov exponents relating bifurcation in Figure 2(a) confirming the parametric space in which the behaviors of system (4) are chaotic or stable period window.
(a)
(b)
For case (ii): by taking parameters as in Table 1, we see that, at the fixed point , a NS bifurcation emerges at and it shows . This confirms Theorem 6.
Figures 3(a) and 3(b) display the bifurcation diagrams of system (4). This illustrates that stability of fixed point happens for , bifurcation occurs at , and there appears attracting invariant circle if . Figure 3(c) results in the maximum Lyapunov exponents relating bifurcation in Figures 3(a) and 3(b) confirming the parametric space in which the behaviors of system (4) are chaotic or stable period window. Figure 3(d) is local amplification of Figure 3(a) for .
(a)
(b)
(c)
(d)
The phase portraits of bifurcation diagrams related to Figures 3(a) and 3(b) for various values of are plotted in Figure 4, which clearly illustrates the act of smooth invariant circle and how it bifurcates from the stable fixed point. We observe that as increases, the invariant circle suddenly disappears and periods 9, 26, and 19 orbits appear at , respectively. We also see that for , there exists a fully developed chaos in system (4).
For case (iii): the dynamics of map (4) can change when more parameters vary. The diagrams of bifurcation of map (4) for control parameters , , respectively, and the fixing remaining parameters as in case (ii) are disposed in Figures 5(a) and 5(b). The 3D maximum Lyapunov exponent for two control parameters and is plotted in Figure 5(c) and its 2D projection onto is shown in Figure 5(d). It is easy to find values of control parameters for which the dynamics of map (4) is in chaotic, periodic, or nonchaotic status. Obviously, there are chaotic dynamics for and and nonchaotic dynamics for and (see Figure 4), which are compatible with the signs of Lyapunov exponents in Figure 5(c). This exhibits that as the parameter increases, the map (4) shows more chaotic dynamics.
(a)
(b)
(c)
(d)
4.1. Initial Perturbation
Sensitivity to initial values of a system demonstrates that at the beginning the two trajectories are arbitrarily closely overlapped but a significant difference builds up for future trajectories. For numerical simulation, two perturbed trajectories for state variable of system (4) against the number of iterations with initial points and are generated and plotted in Figure 6(a). The difference (error) between two trajectories is computed by and presented in Figure 6(b), where and are solution trajectories of system (4) without and with initial perturbation. We take 1000 iterations for each initial value. From Figure 6, it is clear that the error is negligible initially but it increases for trajectories in future time. Therefore, a slightly initial perturbation may lead to complex dynamic behavior of the system.
(a)
(b)
4.2. Fractal Dimension
The fractal dimension (FD) [29, 30] which characterizes strange attractors of a map is defined by where are Lyapunov exponents and is the largest integer such that and . For our two-dimensional map (4), the fractal dimension takes the form
Taking parameter values given as in case (i) and by computer simulation, two Lyapunov exponents are and for and and for . So the fractal dimensions of map (4) are The strange attractors given in Figure 4 are consistent with the delta values’ corresponding fractal dimension and this illustrates that the ratio-dependent discrete predator-prey system (4) has a very complex dynamic behavior as the parameter increases. We refer to Figure 7 to understand the system dynamics using fractal dimension.
(a)
(b)
5. Chaos Control
We shall apply a state feedback control method [24, 31] in order to stabilize chaos of system (4) at the state of unstable trajectories. By adding a feedback control law as the control force to system (4), the controlled form of (4) becomeswhere the feedback gains are denoted by and and represent positive fixed point of (4).
The Jacobian matrix of the controlled system (53) is given bywhere , , , , and , , , and are determined by (7). The characteristic equation of matrix iswhere and . Let be the roots of (56), that is, eigenvalues of (55). Then
The solution of the equations determines the lines of marginal stability. These conditions confirm that . Assume that ; then from (58) we have ; that is,
Now, suppose that ; then from (57) and (58), we get ; that is,
Next, assume that ; then from (57) and (58), we obtain ; that is, Then the lines , , and (see Figure 8(a)) in the plane determine a triangular region which keeps stable eigenvalues.
(a)
(b)
(c)
In order to observe how feedback method works and controls chaos at unstable state, we have presented some numerical simulations. We set and the rest as in case(i) as fixed parameters. The initial value is , and the feedback gains are . Figures 8(b) and 8(c) illustrate that, at the fixed point , the chaotic trajectory is stabilized.
6. Discussions
We investigated the dynamic complexities of discrete predator-prey system with ratio-dependent functional response (4) in the closed first quadrant of . We established the existence conditions for a flip bifurcation and a NS bifurcation of system (4) by using center manifold theorem and bifurcation theory. Specifically, we showed that system (4) displays the change of complex dynamical behaviors via flip bifurcation and NS bifurcation including an invariant cycle, orbits of periods 9, 19, and 26, and attracting chaotic sets when parameter changes its value around sets and , respectively. This indicates that, at the state of chaos, the system is unstable. In particular, it is observed that if the dynamics of prey is chaotic, then the dynamics of predator will ultimately tend to zero or tend to unstable periodic orbit. This results in the interaction between predator and prey which can approach oscillatory balance behavior. Moreover, the variation of two control parameters which displayed captivating chaotic dynamics exists in the system and one can directly observe the chaotic phenomenon from two-dimensional parameter-spaces. Finally, we applied a strategy of feedback control to stabilize chaos at the state of unstable trajectories resulting in predictable stable orbits.
Conflicts of Interest
The author declares that there are no conflicts of interest.