Controlling the chaos and bifurcations of a discrete prey-predator model

1.   Introduction

1.1.   Motivation and literature review

The aquaculture or fish forming in seas, oceans, rivers and streams hold massive potential as a tool to provide us healthy and nutritious food. It also presents a good contribution to development in areas of employment like fish farming, rearing, fish culture, handling and fish processing, etc. Due to population growth, it is a great challenge around the globe to find the sources of food for those people. Fish provide a healthy source of protein. It is highly nutritious and provides a cheap and healthy source of protein. Fish is not only a source of food for humans, it also serves as a source of prey for most of the sea life. Thus, naturally, fish plays an important role in aquatic ecosystems. As the demand for fish increases, the harvesting of fish leads to overfishing, which results in the extinction of some species entirely. Overfishing is the act of catching too many fish before they can reproduce, which depletes fish populations over time. Increased demand for seafood, technologies which make fishing easier and a lack of appropriate fishing laws and enforcement cause overfishing. Overfishing reduces biodiversity and may lead to the extinction of some groups of fish or entire fish species. Those animals that rely on fish for food may struggle to find enough sustenance and also become endangered. Increased investment in fishery management and research is needed to understand the fish populations and ecosystems better. If an all-out extinction is to be prevented, meaningful measures must soon be taken. Scientists and researchers believe that this can be achieved by only harvesting the optimum yield. Scientists and researchers develops several numerical models to estimate the innovative trends in population. Within an environment, the population dynamics can be modeled by a system of autonomous differential equations. However, there is no analytical solution for certain differential equations, such as the nonlinear ones. In these situations, the qualitative approach is applied in addition to the numerical strategy.

A generalized numerical method is added for biological growth; it incorporates the acknowledged capabilities encompassed in generalized logistic, Brody, von Bertalanffy, Gompertz, Richards logistic and other models. In addition, a few theoretical mathematical features of interaction between prey and predator have been provided on the presumption that predation has little or no impact on the prey's growth population. Predator population factors have been integrated into development models to take into account the fact that prey populations grow according to logistic and von Bertalanffy theory. It is difficult to accurately predict the prey-predator interaction. Hence, different results might obviously be produced. Studying various angles of these prey-predator interactions by exploring the dynamical characteristics generates a fascinating and significant biological phenomenon. By examining the effects of toxic substances on aquatic ecosystems, Huda et al. ^[1] have explored the behavior of the following prey-predator fishery model:

where $K$ and $r_1$ respectively represent the carrying capacity and logistic growth rate of the prey population; $\gamma_1$ and $\gamma_2$ are the rates of interaction between the prey and predator fish populations; $\alpha$ and $\beta$ respectively denote the prey populations and toxicity coefficients for predatory; $A$ denotes the environmental protection level for prey; $c$ is predator growth rate; $d$ and $q$ respectively denote the natural mortality rate of predatory fish populations and coefficient for predatory population capture activity and the finally; $E$ denotes the harvesting rate of the predator population. Pujaru and Kar ^[2] have investigated the dynamics of the following predator-prey system:

where $r$ and $s$ denote the growth rates of the predator and prey populations, and $L$ and $K$ respectively denote the carrying capacities of the predator and prey populations. Kar ^[3] has examined the behavior of the following fishery model with time delay:

where $x(t)$ and $y(t)$ respectively denote the densities of prey and predator populations; $g(x)$ and $x p(x)$ respectively denote the growth rate of the prey and response function; $\alpha$ and $d$ respectively denote the conversion factor and predator's death rate; $q$ is the predator's catchability coefficient; $q E y$ is the catch rate function and $\tau$ denotes the time delay. Liu et al. ^[4] have studied the dynamical behaviors of the following continuous model:

where $a_i, b_i\ (i = 1, 2)$, $c$, $d$ and $u$ are positive constants, $\frac{y}{u+y}$ denotes the Allee effect and the algicidal bacteria and harvesting functions for algae are respectively denoted by $q_2 m E y$ and $q_1 m E x$. Moreover $E$, $d$ and $\frac{c x y}{d+y^2}$ respectively denote the combined capture effort, half-saturation constant and algal growth. Keong et al. ^[5] have examined the behavior of the following fishery model:

where $\alpha$ denotes the ratio of the prey's growth rate; $\beta$ is the ratio of the product of catchability coefficient of prey and harvesting effort to the growth rate of prey; $\beta$ is the ratio of the product of coefficient of toxicity and growth rate of prey to the square of growth rate of predator by prey and predator's growth rate is $\delta$. Chen et al. ^[6] have explored Hopf bifurcation of a species interaction model. Chen and Wu ^[7] have examined the dynamical behavior of a predator-prey system with a harvesting policy and network connection. For more interesting results in this field, we refer the reader to the work of the most eminent mathematicians ^[8,9,10] and the references cited therein. On the other hand, many investigators have explored the dynamics of discrete models designated by maps or difference equations rather than continuous models described by differential equations because discrete-time models are more reasonable in the case of non-overlapping generations; additionally, such models give more efficient computational results for numerical simulations ^[11,12]. For instance, Elettreby et al. ^[13] have investigated the complex dynamical behavior of the predator-prey model described below:

where all parameters are positive. Santra et al. ^[14] have examined the behavior of the following predator-prey system:

with all involved parameters being positive. Zhang et al. ^[15] have explored the Hopf bifurcation of a delayed biological economic model. Zhang and Zou ^[16] have explored the dynamics of the following model:

where $x$ and $y$ respectively denote the prey and predator densities. Moreover, $\frac {r_0}{K} \left(K-x_t\right) \left(x_t - c\right)$ denotes the per capita growth rate of the prey in the absence of a predator. Chakraborty et al. ^[17] have investigated the dynamical behaviors of the following prey-predator model:

where $h$ and $q$ respectively denote step size and catchability coefficient for the predator population, $a$, $c$, $e$ and $s$ are positive constants and $E$ denotes the harvesting effort. On the other hand, in recent years, many authors have investigated the bifurcation behavior of fractional-order systems. For instance, Mua et al. ^[18] have investigated the bifurcation and hybrid control tactics of a chemical reaction model with delays. Xu et al. ^[19] have explored the dynamics in a discrete predator-prey competitive model with feedback controls. Xu et al. ^[20] have explored the bifurcation mechanism for fractional-order three-triangle multi-delayed neural networks. Xu et al. ^[21] have explored the bifurcation in a fractional-order predator-prey system with mixed delays. For more interesting results in this regard, we refer the reader to ^{[22,23,24,25]}.

1.2.   Mathematical modeling of a discrete prey-predator fishery model with harvesting

This section presents the mathematical model of a discrete predator-prey fishery model with harvesting. In order to formulate the continuous-time fishery model with harvesting, our investigation has been divided into the following three essential parts ^[26]:

● Consider the following model equation in which fish harvesting occurs without a predator:

where $Ex$ denotes harvesting per unit time and $E$ denotes a positive constant that increased the measure of the effort. Additionally, $K$ denotes the natural carrying capacity and $r$ is the linear per capita growth rate. It is important to denote that, if $E = 0$, then model (1.10) reduces to the non-harvesting fish model.

● Now, consider the following model equation in which fish harvesting occurs in the presence of a predator:

where $s$ and $\nu$ respectively denote the positive parameters for the birth rate and death rate of the predator population.

● Finally, if one takes into account that competition within the predator population increases the death rate of the species, the desired fish harvesting model takes the following form:

where $\beta$ is the measure of the predator's rate of prey consumption and $\alpha$ is the measure of the predator's rate of intra-specific competition. On the other hand, the discrete version of the continuous model (1.12), according to the Euler forward formula takes the following form:

which further becomes

where the integral step size is $h$.

1.3.   Main contributions

The purpose of the present study is to explore the dynamical characteristics of a discrete prey-predator fishery model with harvesting, i.e., model (1.14). More precisely, our key contributions in this paper are as follows:

● Demonstration of the local dynamics at equilibria based on the linear stability theory of the discrete fishery model (1.14).

● Examination of the positive invariance of the discrete fishery model (1.14).

● Identification of bifurcation sets and detailed bifurcation analysis at equilibria of the discrete fishery model (1.14).

● Study of chaos via the state feedback control method.

● Validation of theoretical results numerically through the use of Matlab and Mathematica.

1.4.   Structure of the paper

The rest of the paper is organized as follows: In Section 2, a brief local dynamical analysis at fixed points and under the condition of positive invariance are studied, whereas bifurcations and chaos control in the discrete fishery model are examined in Section 3. In Section 4, our main findings are verified numerically. The paper's conclusion, along with future work, is given in Section 5.

2.   Local dynamics

The local dynamics of the discrete fishery model (1.14) is examined in this section. To do this, we first determine the equilibria and the linearized form of fishery model (1.14) as follows:

Theorem 2.1. In $\mathbb{R}_+^2 = \left\{(x, y): x, \ y\ge 0\right\}$, fishery model (1.14) has three fixed points. More precisely the following is true:

$({\rm{i}})$ $\forall\ K, \ r, \ E, \ s, \ \alpha, \ \beta, \ \nu$, $\mathfrak{L}_1 = \left(0, 0\right)$ is a trivial fixed point of the fishery model (1.14).

$({\rm{ii}})$ If $r > E$, then $\mathfrak{L}_2 = \left(\frac{K(r-E)}{r}, 0\right)$ is a semitrivial fixed point of fishery model (1.14).

$({\rm{iii}})$ If $E < \min \left\{\frac{\nu \beta+\alpha r}{\alpha}, \frac{r K s-r \nu}{K s}\right\}$, then $\mathfrak{L}_3 = \left(\frac{K(\nu \beta+\alpha r-\alpha E)}{K s \beta+\alpha r}, \frac{r K s-r \nu-K s E}{K s \beta+\alpha r}\right)$ is an interior fixed point of fishery model (1.14).

Proof. See in Appendix A. □

Now, at $\mathfrak{L}$, the linearized form of (1.14) under the map

is

where

and

Now, it should be noted that, at $\mathfrak{L}_{1}$, (2.3) becomes

with

Theorem 2.2. For $\mathfrak{L}_{1}$, the following characteristics hold:

$({\rm{i}})$ If

and

then $\mathfrak{L}_{1}$ of fishery model (1.14) is a sink.

$({\rm{ii}})$ If (2.8) holds and

then $\mathfrak{L}_{1}$ of fishery model (1.14) is a source.

$({\rm{iii}})$ If (2.8) holds and

or

then $\mathfrak{L}_{1}$ of fishery model (1.14) is a saddle.

$({\rm{iv}})$ If

or

then $\mathfrak{L}_{1}$ of fishery model (1.14) is non-hyperbolic.

Proof. See in Appendix B. □

Now, at $\mathfrak{L}_{2}$, (2.3) gives

with

Theorem 2.3. For $\mathfrak{L}_{2}$, the following characteristics hold:

$({\rm{i}})$ If

and

then $\mathfrak{L}_{2}$ of fishery model (1.14) is a sink.

$({\rm{ii}})$ If (2.17) holds and

then $\mathfrak{L}_{2}$ of fishery model (1.14) is a source.

$({\rm{iii}})$ If (2.17) holds and

or

then $\mathfrak{L}_{2}$ of fishery model (1.14) is a saddle.

$({\rm{iv}})$ If (2.17) holds and

or

then $\mathfrak{L}_{2}$ of fishery model (1.14) is non-hyperbolic.

Proof. See in Appendix C. □

Finally, at $\mathfrak{L}_{3}$, (2.3) gives

with

where

Additionally, the roots of (2.24) are

where

Theorem 2.4. Let $\Delta < 0$; then, for $\mathfrak{L}_{3}$, the following holds:

$({\rm{i}})$ If

then $\mathfrak{L}_{3}$ of fishery model (1.14) is a stable focus.

$({\rm{ii}})$ If

then $\mathfrak{L}_{3}$ of fishery model (1.14) is an unstable focus.

$({\rm{iii}})$ If

then $\mathfrak{L}_{3}$ of fishery model (1.14) is non-hyperbolic.

Proof. See in Appendix D. □

Theorem 2.5. Let $\Delta > 0$; then, for $\mathfrak{L}_{3}$, the following holds:

$({\rm{i}})$ If

then $\mathfrak{L}_{3}$ of fishery model (1.14) is a stable node.

$({\rm{ii}})$ If

then $\mathfrak{L}_{3}$ of fishery model (1.14) is an unstable node.

$({\rm{iii}})$ If

or

then $\mathfrak{L}_{3}$ of fishery model (1.14) is non-hyperbolic.

Proof. See in Appendix E. □

2.1.   Invariant set

In order to examine the positively invariant set, it is noted that the model (1.14) can also be expressed as

Definition 2.6. The set $S$ is called invariant with respect to the operator $\Omega$ if $\Omega (S)\subset S$.

Regarding the positively invariant set of the discrete model (1.14), one has the following proposition.

Proposition 2.7. The set

is positively invariant with respect to $\Omega$.

In the next section, we will analyze bifurcations at $\mathfrak{L}_{1, 2, 3}$ according to the bifurcation theory ^[27,28].

3.   Bifurcations and controlling chaos

It is noted that, if (2.12) holds, then from (2.6), one gets

but,

which suggests that fishery model (1.14) may undergo flip bifurcation if $(\alpha, \beta, \nu, r, K, E, s)$ belongs to the following set:

But, by calculation, flip bifurcation cannot exist; so, $\mathfrak{L}_{1}$ is degenerate for higher codimensions. Further, if (2.21) holds, then it is also easy to verify from (2.15) that

but, $\lambda_1|_{(2.21)} = -1$, which implies that fishery model (1.14) may undergo flip bifurcation if $(\alpha, \beta, \nu, r, K, E, s)$ belongs to the following set:

The following theorem gives the proof of the fact that, if $(\alpha, \beta, \nu, r, K, E, s)\in \mathcal{FB}|_{\mathfrak{L}_{2}}$, then fishery model (1.14) do not undergo flip bifurcation.

Theorem 3.1. At $\mathfrak{L}_{2}$, if $(\alpha, \beta, \nu, r, K, E, s)\in \mathcal{FB}|_{\mathfrak{L}_{2}}$, then no flip bifurcation occurs.

Proof. Since model (1.14) is invariant with respect to $y = 0$, one can restrict it on $y = 0$ in order to study flip bifurcation, where (1.14) takes the following form

From (3.3), one can write

Finally, if

and

then, from (3.4), we get

and

From (3.7), it should be noted that no flip bifurcation at $\mathfrak{L}_{2}$ occurs if $(\alpha, \beta, \nu, r, K, E, s)\in \mathcal{FB}|_{\mathfrak{L}_{2}}$. □

Further, if (2.30) holds, then $V|_{\mathfrak{L}_{3}}$ at $\mathfrak{L}_{3}$ has complex eigenvalues where $\left|\lambda_{1, 2}\right|_{ (2.30)} = 1$, which implies that fishery model (1.14) may undergo Neimark-Sacker bifurcation if $(\alpha, \beta, \nu, r, K, E, s)$ belongs in the following set:

Theorem 3.2. At $\mathfrak{L}_{3}$, Neimark-Sacker bifurcation exists if $(\alpha, \beta, \nu, r, K, E, s)\in \mathcal{NB}|_{\mathfrak{L}_{3}}$.

Proof. If $(\alpha, \beta, \nu, r, K, E, s)\in \mathcal{NB}|_{\mathfrak{L}_{3}}$ with $h$ as the bifurcation parameter, then model (1.14) takes the following form:

Moreover, the roots of $V|_{\mathfrak{L}_{3}}$ at $\mathfrak{L}_{3}$ for the $\epsilon$-involving model (3.9) is

where

From (3.10) and (3.11), we have

and

Now, occurrence of Neimark-Sacker bifurcation for fishery model (3.9) at $\mathfrak{L}_{3}$ requires that $\lambda_{1, 2}^{\tau}\neq 1, \ \tau = 1, \cdots, 4$ if $\epsilon = 0$, which is equivalent to $p(0)\neq -2, 0, 1, 2$. However, if (2.30) holds, then from (3.11), we get that $q(0) = 1$. So, $p(0)\not = -2, 2$, and, hence, one only requires that $p(0)\not = 0, 1$, which gives

Now, $\mathfrak{L}_{3}$ of the discrete fishery model (3.9) is transformed to the origin by using the following transformations:

where

and

From (3.15) and (3.9), we get

Next, if $\epsilon = 0$, then we study the normal form of (3.16). For this, system (3.16) becomes

with

Now, we can transform the linear part of (3.17) into canonical form by applying the following transformation:

with

From (3.19) and (3.17), we get

where

and

From (3.22), we get

Finally, regarding the occurrence of Neimark-Sacker bifurcation, it is required that the following quantity be non-zero ^[29,30,31]:

where

So, the computation yields

Finally, from (3.25) and (3.27), if we get that $\Gamma\not = 0$ as $(\alpha, \beta, \nu, r, K, E, s)\in \mathcal{NB}|_{\mathfrak{L}_{3}}$, then, at $\mathfrak{L}_{3}$, the discrete fishery model (1.14) undergo Neimark-Sacker bifurcation. In addition, supercritical (resp. subcritical) Neimark-Sacker bifurcation occurs if $\Gamma < 0$ (resp. $\Gamma > 0$). □

Remark 1. The Neimark-Sacker bifurcation has important biological implications in prey-predator models. In these models, the Neimark-Sacker bifurcation indicates the beginning of periodic decreases and increases in the populations of the predator and prey species. It has similar cycles: as the prey population grows, more predators become available, which, in turn, leads to a decline in the prey population, and so on. Because of the complex and reciprocal links that exist between various species, this pattern emerges.

Finally, if (2.33) holds, then the eigenvalues of $V|_{\mathfrak{L}_{3}}$ at $\mathfrak{L}_{3}$ satisfy that $\lambda_1|_{(2.33)} = -1$, but

which confirms that the discrete fishery model (1.14) may undergo flip bifurcation if $(r, s, K, E, \alpha, \beta, \nu)$ belongs in the following set:

But, the following theorem shows that flip bifurcation must exist if $(r, s, K, E, \alpha, \beta, \nu)\in \mathcal{FB}|_{\mathfrak{L}_{3}}$.

Theorem 3.3. At $\mathfrak{L}_{3}$, flip bifurcation exists if $(r, s, K, E, \alpha, \beta, \nu)\in \mathcal{FB}|_{\mathfrak{L}_{3}}$.

Proof. Since $(r, s, K, E, \alpha, \beta, \nu)\in \mathcal{FB}|_{\mathfrak{L}_{3}}$, the discrete fishery model (1.14) takes the form of (3.9) if $h$ is in a neighborhood of $h^*$. Additionally, fishery model (1.14) becomes

where

by (3.15). Moreover, system (3.30) becomes

where

by

Now, in a small neighborhood of $\epsilon$, the center manifold $M^cO$ at $O$ for the map (3.32) is studied. Precisely, $M^cO$ at $O$ can be expressed as

where, by computation, one gets

So, the map (3.32), restricted to $M^cO$, can be expressed as

where

Finally, the following non-zero discriminatory quantities imply that map (3.37) undergoes flip bifurcation ^[32,33,34]:

After extensive manipulation, one gets

and

where

Therefore, from (3.41), along with (3.42), if $\epsilon_2\not = 0$ as $(\alpha, \beta, \nu, r, s, K, E)\in \mathcal{FB}|_{\mathfrak{L}_{3}}$, then, at $\mathfrak{L}_{3}$, fishery model (1.14) undergoes flip bifurcation. Additionally, by the bifurcation theory ^[35,36,37], if $\epsilon_2 > 0$ (resp. $\epsilon_2 < 0)$, then period-$2$ points bifurcated from $\mathfrak{L}_{3}$ are stable (resp. unstable). □

Remark 2. We are choosing the step size, represented by the symbol $h$, as the bifurcation parameter, and it is significant. This is because $h$ represents the time interval between successive iterations in a discrete-time model. By changing the value of $h$, we can simulate different rates at which populations of prey and predators change over time. Biologically, $h$ reflects the reproductive and growth rates of these populations. A smaller value of $h$ means that we have a finer time resolution, allowing for more frequent population size updates. This can help to capture rapid population changes, such as fast reproduction or predation events, and to provide a more detailed representation of the dynamics. A larger value of $h$ means that we have a coarser time resolution, with fewer updates of the population sizes. This can be useful in scenarios in which the population changes occur over longer time scales, such as seasonal variations or slower ecological processes. A larger value of $h$ can help to capture the overall trends and dynamics of the populations. By varying the step size as the bifurcation parameter, we can explore how different rates of population change affect the stability and bifurcation behavior of the system. This can provide insights into the resilience of the ecological system and its response to changes in reproductive rates, predation rates or other important ecological processes. In summary, choosing the step size $h$ as the bifurcation parameter in a discretized prey-predator system with harvesting allows us to analyze how different rates of population change influence the system's dynamics. This can have biological significance as a result of capturing the effects of reproductive and growth rates on the stability and behavior of the ecological system.

Hereafter, chaos control is studied by adding $u_t$ as the control force to fishery model (1.14), as motivated by the existing literature ^[38,39]:

where

$k_1, \ k_2$ denotes feedback gains and

Moreover, $V^C|_{\mathfrak{L}_{3}}$ at $\mathfrak{L}_{3}$ for control system (3.43) becomes

where

Next, if $\lambda_{1, 2}$ are roots of $V^C|_{\mathfrak{L}_{3}}$ at $\mathfrak{L}_{3}$, then

Theorem 3.4. Controlled system (3.43) is asymptotically stable if both eigenvalues of $V^C|_{\mathfrak{L}_{3}}$ at $\mathfrak{L}_{3}$ satisfy that $|\lambda_{1, 2}| < 1$ where $k_1$ and $k_2$ denote the feedback gains.

Proof. See in Appendix F. □

4.   Numerical simulations

To validate the theoretical findings, the following simulations are presented.

Case A: Let

then, from (2.30), we obtain that

which implies that, if

then $\mathfrak{L}_{3}$ of fishery model (1.14) is a stable focus. Further, if

then Figure 1a shows that

of the discrete fishery model (1.14) is a stable focus. Also, if

then Figure 1b–1f shows that the corresponding equilibrium $\mathfrak{L}_{3}$ of fishery model (1.14) is a stable focus.

On the other hand, if

then Figure 2a shows that the interior fixed point

changes the dynamics and, as a result, stable curves appear. Now, we further numerically prove that, if

then the discrete fishery model (1.14) undergoes supercritical Neimark-Sacker bifurcation, that is, $\Gamma < 0$. So, if $h = 1.56$, then, from (3.13), we obtain

Additionally, from (3.10) and (3.27), one has

and

Applying (4.1) and (4.2) in (3.25), one gets

which confirms that our findings are numerically true; so fishery model (1.14) undergoes supercritical Neimark-Sacker bifurcation. Similarly, if

Figure 2b–2f implies that stable curves also appear, and therefore, the discrete fishery model (1.14) undergoes supercritical Neimark-Sacker bifurcation. If

then $\Gamma < 0$ (see Table 1).

Finally, the maximum Lyapunov exponent and bifurcation diagrams are plotted in Figure 3.

Case B: Let $\alpha = 5.8, \ \beta = 3.0, \ \nu = 0.09, \ r = 15.52, \ s = 4.501, \ K = 0.52$ and $E = 0.001$; then, from (2.33), we obtain

Moreover, for $\alpha = 5.8, \ \beta = 3.0, \ \nu = 0.09, \ r = 15.52, \ s = 4.501, \ K = 0.52$ and $E = 0.001$, we also get

From a theoretical point of view, the interior fixed point $\mathfrak{L}_{3}$ of the discrete two-dimensional fishery model (1.10) is a stable node if

So, if

then Figure 4a implies that the interior fixed point

of fishery model (1.14) is a stable node. Also, if

then Figure 4b–4d also implies that

of fishery model (1.14) is also a stable node.

Further, if

then $\mathfrak{L}_{3}$ becomes an unstable node and, meanwhile, flip bifurcation occurs, that is, if

then, from (3.40), one gets

Moreover, in view of (3.41) and (3.42), one gets

which shows that stable period-$2$ points bifurcate from

Finally, the maximum Lyapunov exponent and bifurcation diagrams are plotted in Figure 5.

Case C: The results of Section 3 will be verified in this case. For this, if $\alpha = 5.8, \ \beta = 3.0, \ \nu = 0.09, \ r = 15.52, \ s = 4.501, \ K = 0.52, \ E = 0.001$ and $h = 0.95$, then, from (F.1)–(F.3), one gets

Therefore, lines (4.3)–(4.5) give the triangular region for (3.43) that satisfies that $|\lambda_{1, 2}| < 1$ (see Figure 6).

5.   Conclusions

In this work, we have examined the local behavior at fixed points, chaos control via the state feedback control method and bifurcations of the fishery model (1.14). Moreover, it is shown that, for all of the model's parameters, the discrete fishery model (1.14) has trivial fixed point; however, it has semitrivial and interior fixed points under certain parametric condition(s). By using linear stability theory, we have examined the local stability at equilibria of the discrete fishery model (1.14). Further, we have also examined the existence of bifurcations at fixed points and proved that, at trivial and semitrivial fixed points, the discrete fishery model (1.14) does not undergo flip bifurcation, whereas, at interior fixed point, model (1.14) undergoes both Neimark-Sacker and flip bifurcations, and we have studied these bifurcations by applying bifurcation theory accordingly. Further, chaos control has been studied by using the state feedback control method. Finally, the main results have been confirmed numerically. Our numerical simulations shows that, at interior fixed point, the discrete fishery model (1.14) undergoes supercritical Neimark-Sacker bifurcation and biologically, the occurrence of supercritical Neimark-Sacker bifurcation means that there exists a periodic or quasiperiodic oscillations between prey and predator populations.

Future work: Bifurcation analysis of the (ⅰ) fractional-order discrete-time prey-predator model (1.14) and (ⅱ) discrete-time prey-predator model (1.14) with time delay, are our next goals of study.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Conflict of interest

The authors declare that they have no conflicts of interest regarding the publication of this paper.

Appendix

A.   Existence of equilibria

If $\mathfrak{L} = \left(x, y\right)$ is a fixed point of fishery model (1.14), then

After straightforward calculation, from (A.1), one gets

which implies that, if

then

is an interior fixed point of fishery model (1.14). On the other hand, for

Equation (A.1) is satisfied in an identical manner; so, $\mathfrak{L}_{1, 2}$ respectively denote trivial and semitrivial fixed points of fishery model (1.14).

B.   Local dynamics at $\mathfrak{L}_{1}$

Recall from (2.6) that, at $\mathfrak{L}_{1}$, the characteristic roots of $V|_{\mathfrak{L}_{1}}$ are

Therefore, by the stability theory $\mathfrak{L}_{1}$, of the fishery model (1.14) is a sink if

and

which imply that if

and $E > r$ then $\mathfrak{L}_{1}$ is a sink. In similar way, one can obtain that $\mathfrak{L}_{1}$ is a source if

a saddle if

or

and non-hyperbolic if

C.   Local dynamics at $\mathfrak{L}_{2}$

From (2.15), if

and

i.e.,

and

then $\mathfrak{L}_{2}$ of fishery model (1.14) is a sink. Similarly, it is also easy to prove that $\mathfrak{L}_{2}$ of fishery model (1.14) is a source if (2.17) holds and

a saddle if (2.17) holds and

or

and, finally, non-hyperbolic if (2.17) holds and

D.   Local dynamics at $\mathfrak{L}_{3}$ if $\Delta < 0$

If $\Delta < 0$, then the characteristic roots of $V|_{\mathfrak{L}_{3}}$ at $\mathfrak{L}_{3}$ are

which yield that, if

i.e.,

then $\mathfrak{L}_{3}$ of fishery model (1.14) is a stable focus. In a similar way, one can obtain that $\mathfrak{L}_{3}$ of fishery model (1.14) is an unstable focus if

and non-hyperbolic if

E.   Local dynamics at $\mathfrak{L}_{3}$ if $\Delta > 0$

If $\Delta > 0$, then the roots are

By the linear stability theory, $\mathfrak{L}_{3}$ of the discrete of fishery model (1.14) is a stable node if $|\lambda_{1, 2}| < 1$, i.e.,

In a similar way, one can also establish that $\mathfrak{L}_{3}$ of fishery model (1.14) is an unstable node if

and non-hyperbolic if

or

F.   Existence of triangular region for controlled system (3.43)

It is noted that marginal stability can arise from the restrictions $\lambda_1 = \pm 1$ and $\lambda_1\lambda_2 = 1$, which imply that $|\lambda_{1, 2}| < 1$. If $\lambda_1\lambda_2 = 1$, then, from (3.47), we get

If $\lambda_1 = 1$, then, from (3.46) and (3.47), we get

Finally, if $\lambda_1 = -1$ then from (3.46) and (3.47) we get

Therefore, from (F.1)–(F.3), lines $L_1$–$L_3$ in the $(k_1, k_2)$-plane from the triangular region, which implies that the characteristic roots of controlled system (3.43) satisfy that $|\lambda_{1, 2}| < 1$.