Abstract

The paper investigates the dynamical behaviors of a two-species discrete predator-prey system with Crowley–Martin functional response incorporating prey refuge proportional to prey density. The existence of equilibrium points, stability of three fixed points, period-doubling bifurcation, Neimark–Sacker bifurcation, Marottos chaos, and Control Chaos are analyzed for the discrete-time domain. The time graphs, phase portraits, and bifurcation diagrams are obtained for different parameters of the model. Numerical simulations and graphics show that the discrete model exhibits rich dynamics, which also present that the system is a chaotic and complex one. This paper attempts to present a feedback control method which can stabilize chaotic orbits at an unstable equilibrium point.

1. Introduction

The dynamical behavior of the prey-predator system in ecology has been the subject of study for many researchers and hence they have provided a substantial contribution to the growth of the population models [18] in the continuous-time domain. In case of populations that have overlapping generations and the birth processes occuring continuously, the ordinary differential equations are used for modeling of predator-prey interaction [930]. Many species, such as monocarpic plants, and semelparous animals have discrete non overlapping generations, and their births occur in regular breeding seasons. Their interactions are characterized by difference equations or formulated as discrete-time mappings. Discrete-time predator-prey models are capable of exhibiting much more complicated dynamics in comparison to the corresponding continuous-time models [3145].

The dynamical study of the prey-predator model with different functional responses is an effective approach for time domain. Upadhyay and Naji [46] presented three-species food chain model with Crowley–Martin type functional response for the dynamical study. Zhou [47] presented qualitative analysis of the predator-prey model for Crowley–Martin functional responses. Tripathi et al. [48] applied Crowley–Martin functional response on a delayed density dependent predator-prey model. Li et al. [49] presented a predator-prey model with Crowley–Martin functional response for positive periodic solutions.

A certain degree of protection can be offered by nature to a set number of prey populations by providing refuges. The effects of prey refuges on the population dynamics are intricate and complex in nature; however, it can be considered to be constituted by two components for modeling purposes. The first component of refuge effects, which governs positive growth of prey and negatively that of predators, consists of a reduction in prey mortality owing to decreasing predation success. The second one may arise as the trade-offs and byproducts arising from the prey’s hiding behavior, which could cause advantages or even be detrimental for all the interacting populations. González-Olivares and Ramos-Jiliberto [8] presented prey refuges in a simple prey-predator system. Huang et al. [9] presented stability analysis of a prey-predator system incorporating a prey refuge. Ma et al. [12] discussed the effect of prey refuges on a prey-predator system. Chen et al. [13] studied the prey refuge on a Leslie–Gower predator-prey model.

Extensive studies have been carried out by many researchers on refuge concept [930], and most of these papers consider refuge as a constant number of prey or proportional to prey. These earlier works have been taken into account that the use of refuges by a fraction of prey, or a constant number prey, is influential in existing a stabilizing effect in the dynamics of the interacting populations. This paper will test the above statement assuming that the quantity of prey in refuges is proportional to the prey. We also conduct an analysis of the dynamic properties of such a system via modification of the well-known Lotka–Volterra predator-prey model with prey self-limitation incorporating Crowley–Martin functional response [46, 47]. In this paper, we consider the refuge term proportional to prey density.

The rest of the paper is organized as follows. In Section 2, a discrete-time prey-predator model under nonoverlapping generation with refuge is formulated. Analysis of the local stability of the proposed model is dealt with in Section 3. The bifurcation analysis of the proposed model is studied through period-doubling bifurcation and Neimark–Sacker bifurcation in Section 4. Sections 5 and 6 give Marottos chaos and control chaos procedures, respectively. The proposed model is supported with the help of numerical simulations in Section 7. Section 8 makes the conclusion to this paper.

2. Description of the Proposed Predator-Prey Model

The proposed model belongs to a generalized predator-prey model which incorporates logistic growth having no stage structure for neither the prey nor the predator, which is represented by the following differential equations:where x and y denote the density of prey and predator populations at any time t, respectively, and are all positive constants and have their biological meanings accordingly. The intrinsic per capita growth rate of prey population is denoted by r, k is the environmental carrying capacity of prey population, c is the maximal per capita consumption rate of predators, d is the efficiency with which predators convert consumed prey into new predators, and f is the per capita death rate of predators. The term represents the functional response of the predator population.

According to Maynard Smith [4], there is a quantity of prey population which incorporates refuges. Thereafter, the functional response is modified to incorporate prey refuges in the above model, and the term becomes the functional response in the form of Based on a laboratory study of dragonfly, Crowley and Martin present functional response , where describe the handling time and is the magnitude of interference among predators. In the proposed model, we consider refuge , and to derive a discrete-time model let and , where and are the densities of the prey and predator populations in discrete time . Let and , then the generation of the prey and predator populations by replacing t by n is obtained as follows:

Let us consider and then from (2) we obtain the following proposed discrete-time predator-prey system as follows:where , and β are all positive constants. By biological meaning of the model variables, we only consider the system in the region in the plane.

3. Fixed Points and Stability Analysis of the Prey-Predator System

Fixed points of system (3) are determined by solving the following nonlinear system of equations:

By simple calculation, we get three nonnegative fixed points as follows: , , , and , where is positive root of the equation , where , , , and , and .

Lemma 1. Equation has a unique positive root if , three positive roots if , and a unique positive root if .

Now, the local behavior of system (3) is studied for each equilibrium point of the prey-predator system. The stability of system (3) is carried out by computing the Jacobian matrix corresponding to each equilibrium point. The Jacobian matrix J for system (3) iswherewhere

The characteristic equation of matrix J is , where Trace of matrix  Determinant of matrix 

Hence, model (3) is a dissipative dynamical system, if , is a conservative dynamical one, if and only if , and is an undissipated dynamical system otherwise.

Proposition 1. The equilibrium point is Sink if , Saddle if , and Nonhyperbolic if .

Proposition 2. The equilibrium point is Sink if and , Source if and , Saddle if and or and , and Nonhyperbolic if or

Proposition 3. Interior equilibrium point is Sink if one of the following conditions holds: and and

Proposition 4. Interior equilibrium point is Source if one of the following conditions holds: and and

Proposition 5. Interior equilibrium point is Saddle if one of the following conditions holds: and and

Proposition 6. Interior equilibrium point is Nonhyperbolic if one of the following conditions holds: and and and and The prey-predator model (3) undergoes period-doubling bifurcation at , when parameters vary in a small neighborhood of or . Model (3) undergoes Neimark–Sacker bifurcation at , when parameters vary in a small neighborhood of or .

4. Bifurcation Analysis

This section carries out an investigation of the conditions for the existence of Neimark–Sacker bifurcation (NSB) and period-doubling bifurcation (PDB) at the positive fixed point of proposed system (3). There is an emergence of different kinds of bifurcations from the fixed point in dynamical systems, when a particular parameter passes through its critical value. Many dynamical properties of a system can be discussed owing to emergence of NSB and PDB. We discuss NSB and PDB for the positive fixed point of the prey-predator system (3) taking b as a bifurcation parameter.

4.1. Neimark–Sacker Bifurcation

This section presents NSB of the proposed prey-predator model (3) at for the parameters are located in the following set:

Similar arguments can be applied to the other case:

In analyzing the NSB, b is used as the bifurcation parameter. Furthermore, () is the perturbation of and we consider a perturbation of the model as follows:

Let , then the equilibrium point is transformed into the origin, and further expanding f and as a Taylor series at to the third order, the proposed model (10) becomeswhere

The characteristic equation associated with the linearization of model (11) at is given by . Roots of the characteristic equation are

From , when , we have and .

In addition, when this is equivalent to .

Now to study the normal form, let and . We define , and using the transformation model (11) becomeswhere the functions and denote the terms in model (14) in variables with the order at least two.

In order to undergo Neimark–Sacker Bifurcation, we require the following discriminatory quantity to be nonzero:where

From the above analysis we have the following result.

Theorem 1. If , then model (3) undergoes NSB at when the parameter varies in a small neighborhood of the origin. Moreover, if then an attracting (repelling) inariant closed curve bifurcates from for .

4.2. Period-Doubling Bifurcation

It is clear that one of the eigenvalues of the positive fixed point is and the other is neither 1 nor if parameters of the model are located in the following set:

Here, we discuss PDB of the proposed model (3) at when parameters vary in a small neighborhood of . Similar arguments can be applied to the other cases. In analyzing the PDB, b is used as the bifurcation parameter. Furthermore, () is the perturbation of and we consider a perturbation of the model as follows:

Let then the equilibrium point is transformed into the origin, and further expanding f and as a Taylor series at to the third order, model (18) becomeswhere , , , , , , , , , , , , , , , , , , , and .

We define , and it is obvious that T is invertible. Using the transformation , then model (19) becomeswhere the functions and denote the terms in model (20) in variables with the order at least two.

From the center manifold theorem, we know that there exists a center manifold of model (20) at in a small neighborhood of , which can be approximately described as follows:where

Model (20) which is restricted to center manifold has the following form:where

For flip bifurcation, we require the two discriminatory quantities and to be nonzero:

Finally, from the above analysis, we have the following result.

Theorem 2. If and , then model (3) undergoes PDB at when the parameter b varies. Moreover, if then the period-2 points that bifurcation from are stable (unstable).

5. Existence of Marottos Chaos

This section presents chaotic nature of system (3) in the sense of Marrotto.

Definition 1. Let the function be differentiable in . The point is an expanding fixed point of F in if and all eigenvalues of exceed 1 in norm for all .

Definition 2. Assume that Z is an expanding fixed point of F in for some . Then, Z is said to be a snap-back repeller of F if there exists a point with , , and for some positive integer M.

Theorem 3. If , then is an expanding fixed point of .

Proof. For map,The eigenvalues corresponding with the fixed point are given bywhereTherefore, the fixed point has a pair of complex eigenvalues, and the norm of them exceeds unity if and Thus, we can state the following theorem.

Theorem 4. is a snap-back repeller in .

Proof. According to the definition of a snap-back repeller, one needs to find one point such that , , and , for some positive integer M, where Map F is defined by (3).
To proceed, notice thatNow, a map has been constructed to map the point to the fixed point after two iterations if there are solutions different from for equations (29) and (30). The solutions different from for equation (30) satisfy the following equation:Substituting and into equation (29) and solving , we haveBy simple calculations, we obtainwhereThe solutions of equations (31) and (32) will farther subject to and then is a snap-back repeller in .
Thus, the following theorem is established.

Theorem 5. If , the solutions and of equations (31) and (32) satisfy in addition , , , and then is a snap-back repeller of Map (3), and hence Map (3) is chaotic in the sense of Marotto.

6. Chaos Control of the Proposed System

Several methods can be used for obtaining chaos control in discrete-time models. A few of these methods are the state feedback method, pole-placement technique, and hybrid control method. We give a feedback control method to stabilize chaotic orbits at an unstable positive fixed point of the prey-predator system (3), with the following controlled form:with the following feedback control law as the control force:where and are the feedback gain and is a positive fixed point of the model.

The Jacobian matrix J for system (35) at iswhere

The corresponding characteristic equation of matrix J is

Let and are the eigenvalues of the above characteristic equation:

The lines of marginal stability are determined by solving the equations and . These conditions guarantee that the eigenvalues and have modulus less than 1.

Suppose , from (40), we have line as follows:

Suppose , from (40), we have line and as follows:

The stable eigenvalues lie within a triangular region by line , and .

7. Numerical Simulations

This section presents the qualitative dynamical behaviors of the system by presenting Lyapunov exponent, bifurcation diagram, and phase plane for specific parameter values. We present the numerical simulation to observe rich dynamics of our model, so we consider the hypothetical parameters values of the proposed discrete prey-predator system. The value of the parameters are consider as , , , , , and with initial population for simulation and graphical presentation of different dynamical observations. Figures 13 are iterative plots for different refuge parameters, and we observed that refuge has a stabilizing effect. Figures 49 are bifurcation diagrams and corresponding local amplification diagram of the prey and predator for varying parameters and β, respectively.

The phase portraits associated with Figure 5 are displayed in Figures 13. A chaotic attractor for and 0.01. A five-pieces chaotic attractor for and 0.03. Period-10 orbit for . Period-68 orbit for . Attracting invariant cycle for , and 0.152. Fixed point for and 0.160.

The orbit diagram of species for a, which is a scaling parameter of the intrinsic per capita growth rate of prey population. Figure 4 shows that the system trajectory evolves from a fixed point to NSB and finally into a chaotic attractor. There is one visible periodic windows in the bifurcation diagram for both prey and predator. We also give local magnification of the corresponding bifurcation figure. Figure 10 is the corresponding Lyapunov exponents of Figure 4.

Figure 5 shows the orbit diagram of prey and predator population for the refuge parameter with other fixed parameter values. The system trajectory evolves into a periodic orbit from the chaotic attractor and finally into a fixed point. There is only one obvious periodic window in all the bifurcation process for both prey and predator. We also give local magnification of the corresponding bifurcation figure. Figure 11 is the corresponding Lyapunov exponent of Figure 5.

Figure 6 shows the orbit diagram of the population for the parameter of per capita consumption rate of predators. The system trajectory evolves from a fixed point to NSB and finally into a chaotic attractor. We also give local magnification of the corresponding bifurcation figure. Figure 12 is the corresponding Lyapunov exponent of Figure 6.

Figure 7 shows the orbit diagram of species for the parameter of efficiency with which predators convert the consumed prey into a new predator. The system trajectory evolves from a fixed point to NSB and finally into a chaotic attractor. There is one big visible periodic window in the bifurcation diagram for both prey and predator. We also give local magnification of the corresponding bifurcation figure. Figure 13 is the corresponding Lyapunov exponent of Figure 7.

Figure 8 shows the orbit diagram of prey and predator population for the parameter α with other fixed parameter values. The system trajectory evolves into a fixed point from the periodic orbit. We also give local magnification of the corresponding bifurcation figure. Figure 14 is the corresponding Lyapunov exponent of Figure 8.

Figure 9 shows the orbit diagram of prey and predator population for the parameter β with other fixed parameter values. The system trajectory evolves into a periodic orbit from the chaotic attractor and finally into a fixed point. There is only one obvious periodic window in all the bifurcation process for both prey and predator. We also give local magnification of the corresponding bifurcation figure. Figure 15 is the corresponding Lyapunov exponent of Figure 9.

7.1. Chaos Control

We observe chaotic behavior of prey-predator, as shown in Figure 16(a). In this case fixed point is unstable. In the feedback control method for feedback gain and we observe the fixed point is stable, as shown in Figure 16(b).

7.2. Fractal Dimension

The fractal dimension is defined by using Lyapunov exponents as follows.

Let f be a map on . Consider an orbit with Lyapunov exponents and let p denote the largest integer such that .

Fractal dimension and Lyapunov exponents for different values of b for the value of the parameters as , , , , and with initial population are presented in Table 1.

8. Conclusion

In this paper, we have studied the effects of refuges on a predator-prey interaction, by using the analytical and graphical approach. The refuges are considered as prey refuge proportional to prey density. We evaluate the effects with regard to the local stability of the interior equilibrium point and the long-term dynamics of the interacting populations. The results show that the effects of refuges can stabilize the interior equilibrium point of the proposed prey-predator discrete-time domain model. We have studied bifurcations in a discrete predator-prey model with refuge. Furthermore, it is also shown that the model exhibits various bifurcations of codimension 1, including period-doubling bifurcation, and Neimark–Sacker bifurcation as the values of parameters vary. Numerical simulations and graphical presentation show the rich dynamics in bifurcation and chaotic nature of the system model.

Data Availability

No data were used to support this study.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.