Approximate solution of nonlinear Black–Scholes equation via a fully discretized fourth-order method

1.   Introduction

1.1.   Motivation and literature survey

Ecological balance is defined as a "state of dynamic equilibrium within an organismal community in which genetic, species, and ecosystem diversities remain relatively stable". An ecological imbalance results from a breakdown of an ecosystem's natural equilibrium. Natural habitat degradation, climate change, global warming, and biodiversity loss are all results of an ecological imbalance. Therefore, ecological balance is necessary to preserve the diverse and abundant organismal variety. The survival, reproduction, and fitness of an individual are determined by responses to environmental factors. When one creature uses a habitat or food supply in a certain way, it modifies the temporal and spatial dynamics of the habitat structure and resource distribution for another organism. To maintain the balance and stability, scientists and ecologist study different factors involved in ecology where prey-predator is most essential to consider ^[1]. Ecological relationships between predators and prey are crucial. Predators pursue and consume prey, which has an impact on both population's dynamics. The population dynamics of predators and prey are shown mathematically by the Lotka-Volterra model. The model may be used to estimate how predator and prey populations can evolve under various situations by taking variables like birth rates, death rates, and the availability of resources into consideration. For instance, if there are many preys and few predators, the predator population will probably grow since more food is available. However, if the quantity of predators is too high, it may diminish the number of prey, resulting in a drop in both populations. Ecologists may learn about the numerous links that exist between various species within an ecosystem and create plans for biodiversity preservation and resource management by understanding predator-prey interactions. Prey-predator models have important mathematical implications because they clarify ecosystem dynamics and the interactions of various species. These models enable us to estimate how these populations will evolve over time by using mathematical equations to represent the interactions between predators and prey. For instance, differential equalizations are used in the popular prey-predator Lotka-Volterra model to explain how the populations of predators and prey fluctuate over time. These equations may be used to estimate how predator and prey populations would fluctuate under various circumstances which mentioned above. By analyzing these models, ecologists may learn about the complex interactions between various species in an ecosystem and create plans for biodiversity preservation and resource management ^[2]. Therefore, prey-predator models are building blocks of bio and ecosystems. In an ecosystem, millions of species interact with each other. They compete, change, and disperse simply to seek recourses. In their struggle for existence in this universe, they compete with each other and may attain different attributes or shapes of population interactions. Depending upon their specific settings of applications, they can take the forms of resource-consumer, parasite-host, plant-herbivore, etc. Among prey-predator interactions, predators are those species that kill and eat another species, to be eaten called prey. In this interaction, when the prey population decreases, the predator population also decreases; moreover, when the prey population increases, the predator population also increases. Many models have been proposed to comprehend the process of competition among populations of two species. Among these models, the most well-known model is the following predator-prey model which was proposed by Volterra in 1931 ^[3]:

where the number of prey and predator are denoted by $x$ and $y$ at time $t$, respectively, whereas all involved model's parameters are positive. This model represents the direct relationship between the prey and predator populations. If predators are absent, then the number of prey grows exponentially; when the prey population decreases, the predator population also decreases exponentially. Moreover, in the model, $dxy$ and $bxy$ denote the predator-prey confrontation, respectively, which are useful to predators and harmful to prey. It is noted here that due to the harvesting effect, model (1.1) becomes ^[4]:

Reaction-diffusion, which is another factor in the interaction of prey-predator, was considered by Lazaar et al.; they studied the local and global dynamics of prey-predator by considering the prey refuge ^[5]. Saadeh et al. ^[6] considered the following model to describe different dynamical properties which are effected by commensurate and incommensurate orders:

where $a$, $b$, $c$, and $\delta$ are positive parameters. Elettreby et al. ^[7] studied the dynamics of a discrete prey-predator model with mixed functional response. Many other mathematical models which were considered in continues form showed different results to contribute to the stability of ecosystem. For instance, Chen et al. ^[8] examined the hopf bifurcation of a species interaction model. Chen and Wu ^[9] investigated the dynamics of a diffusive predator-prey model with a harvesting policy and a network connection. Chen and Srivastava ^[10] investigated the bifurcating solution of a predator-prey model. Chen and Wu ^[11] studied bifurcation in the Previte-Hoffman model. Chen and Wu ^[12] studied the dynamics of a diffusive predator-prey model. Britton ^[13] suggested the following Leslie's model:

where all involved model's parameters are positive. The British mathematician Patrick Leslie created the Leslie model in the 1940s for examining the dynamics of fish populations. Leslie was fascinated by how populations of fish and other species changed over time, and he realized that a mathematical model may assist him in forecasting these changes. According to the Leslie model, the population may be split into several age groups, and the birth and death rates of people in each group can be used to forecast how the population will change over time. Since then, a variety of creatures have been studied using the model, making it a crucial tool in ecology and conservation biology. The following assumptions constitute the foundation of the Leslie model ^[14,15]: There are discrete age groups within the population; the birth rate is consistent throughout all age groups; the mortality rate is consistent among all age groups; and there is no inward or outward population migration.

1.2.   Problem statement

Predator-prey models and other population dynamics are studied using the Leslie model; it considers elements including population interactions, mortality rates, and birth rates. The model is helpful for understanding the effects of environmental variables, researching long-term ecosystem dynamics, and anticipating changes through time. Moreover, it is well known that when the populations have non-overlapping generations, the discrete models defined by maps are more logical than the continuous models. Additionally, the discrete-time models offer richer dynamics and more capable computational results for numerical simulations than the continuous ones ^{[16,17,18,19]}. Therefore, in the present study, we will study the hybrid control and bifurcations of the following discrete Leslie's model:

which is a discrete analogue of (1.4) by the Euler forward formula.

1.3.   Layout of the paper

The layout of rest of the paper is as follows. The stability analysis at fixed points is given in Section 2. In Section 3, bifurcations at equilibria are given, whereas Section 4 describes the control of the N-S bifurcation. In Section 5, the obtained results are verified numerically, whereas a concise summary is provided in Section 6.

2.   Stability analysis at equilibria

In this study, we explore stability analysis at equilibria of the discrete Leslie's model (1.5) by the stability theory ^{[20,21,22,23,24,25]}. In this respect, we first find equilibria of the discrete Leslie's model (1.5) by algebraic techniques. By definition of a fixed point, if $E_{xy}(x, y)$ is the equilibrium point of the model (1.5), then

Since $E_{x0}\left(\frac{a}{b}, 0\right)$ satisfied (2.1), then for all model's parameters, $E_{x0}\left(\frac{a}{b}, 0\right)$ is the boundary fixed point of Leslie's model (1.5). For the interior equilibrium point, from (2.1) one has

From $1^{st}$ equation of (2.2), one has

From $2^{nd}$ equation of (2.2) and (2.3), one gets:

Utilizing (2.4) in (2.3), one gets:

In view of (2.4) and (2.5), one can say that for all model's parameters, the Leslie's model (1.5) has an interior fixed point $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$. The variational matrix $\Lambda\biggr |_{E_{xy}(x, y)}$ at $E_{xy}(x, y)$ under the map is as follows:

where

is

Now, in subsequent sections, we will explore the stability analysis at $E_{x0}\left(\frac{a}{b}, 0\right)$ and $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$.

2.1.   Stability analysis at $E_{x0}\left(\frac{a}{b}, 0\right)$

It is stated that at $E_{x0}\left(\frac{a}{b}, 0\right)$, (2.8) implies

whose characteristic roots are

From (2.10), the dynamics of the Leslie's model (1.5) at $E_{x0}\left(\frac{a}{b}, 0\right)$ can be summarized as

Lemma 2.1. $({\rm{i}})$ $E_{x0}\left(\frac{a}{b}, 0\right)$ is never sink;

$({\rm{ii}})$ $E_{x0}\left(\frac{a}{b}, 0\right)$ is a source if

$({\rm{iii}})$ $E_{x0}\left(\frac{a}{b}, 0\right)$ is a saddle if

$({\rm{iv}})$ $E_{x0}\left(\frac{a}{b}, 0\right)$ is non-hyperbolic if

Remark 1. The topological classifications at $E_{x0}\left(\frac{a}{b}, 0\right)$ of Leslie's model (1.5) are given in Figure 1.

2.2.   Stability analysis at $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$

At $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$, (2.8) takes the following form:

whose characteristic equation is

where

Lastly, roots of (2.15) are

where

Hereafter, two Lemmas for $\Delta < 0$ ($\Delta\ge0$ are presented to show the complete topological classifications at $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ of the discrete model (1.5).

Lemma 2.2. $({\rm{i}})$ $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ is a stable focus if

$({\rm{ii}})$ $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ is an unstable focus if

$({\rm{iii}})$ $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ is non-hyperbolic if

Lemma 2.3. $({\rm{i}})$ $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ is a stable node if

$({\rm{ii}})$ $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ is an unstable node if

$({\rm{iii}})$ $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ is non-hyperbolic if

Remark 2. The topological classifications at $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ of Leslie's model (1.5) if $\Delta < 0$ and $\Delta > 0$, respectively, are given in Figures 2 and 3.

3.   Bifurcation analysis

The possible bifurcation analysis are explored in this section at fixed points $E_{x0}\left(\frac{a}{b}, 0\right)$ and $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$. First, we will study the bifurcation analysis at $E_{x0}\left(\frac{a}{b}, 0\right)$. Recall that $\bigwedge \biggr|_{E_{x0}\left(\frac{a}{b}, 0\right)}$ at $E_{x0}\left(\frac{a}{b}, 0\right)$ has eigenvalues, which are depicted in (2.10), and if (2.13) holds, then $\lambda_2\biggr|_{(2.13)} = -1$ but $\lambda_1\biggr|_{(2.13)} = \frac{2d+a}{a}\not = -1 \ \mathrm{or}\ 1$. This grantees the fact that the discrete model (1.5) may undergo the flip bifurcation if $\bigcap: = (a, b, c, d, h, \alpha)$ passes through the curve:

Therefore, in the following theorem, we are going to study the detailed indicated bifurcation analysis at $E_{x0}\left(\frac{a}{b}, 0\right)$.

Theorem 3.1. At $E_{x0}\left(\frac{a}{b}, 0\right)$, the Leslie's model (1.5) does not undergo a flip bifurcation if $\bigcap\in FB\biggr|_{E_{x0}(\frac{a}{b}, 0)}$.

Proof. Since w.r.t $y = 0$, the discrete Leslie's model (1.5) is invariant, and therefore, one restricts (1.5) to line $y = 0$ to explore the indicated bifurcation where it becomes

From (3.2), one has

From (3.1), one denotes $h = h^* = \frac{2}{a}$ and $x = x^* = \frac{a}{b}$. Now, from (3.3) one obtains the following

and

Equation (3.6) clearly indicates that the Leslie's model (1.5) does not undergo a flip bifurcation if $\bigcap\in FB\biggr|_{E_{x0}(\frac{a}{b}, 0)}$, and therefore as a result, the predator population goes to an extension while the prey undergoes a flip bifurcation to chaos. □

Now, the bifurcation analysis at $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ of the Leslie's model (1.5) is investigated. Recall that if $\Delta < 0$, then $\bigwedge \biggr|_{E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)}$ about $E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ has complex conjugate pairs. Moreover, if (2.21) holds, then $\left|\lambda_{1, 2}\right|_{ (2.21)} = 1$. This gives the fact that model (1.5) may undergo a N-S bifurcation if $\bigcap$ crosses the curve:

However, the following theorem gives the detailed analysis of N-S bifurcation at $E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ by the bifurcation theory ^{[26,27,28,29,30,31]}.

Theorem 3.2. If $\bigcap\in NSB\biggr|_{E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)}$ then at $E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$, the N-S bifurcation must exists.

Proof. If parameter $a$ is designated as a bifurcation parameter, then (1.5) becomes

whose interior fixed point is $E^+_{xy}\left(\frac{\left(a^*+\epsilon\right)\alpha}{cd+b\alpha}, \frac{\left(a^*+\epsilon\right)d}{cd+b\alpha}\right)$. Moreover, at $E^+_{xy}\left(\frac{\left(a^*+\epsilon\right)\alpha}{cd+b\alpha}, \frac{\left(a^*+\epsilon\right)d}{cd+b\alpha}\right)$, (2.15) becomes

where

Roots of (3.9) are

From (3.11), one computes

Additionally, one needs to determine that $\lambda_{1, 2}^\upsilon\neq 1, \upsilon = 1, \cdots, 4$ which is equivalent to $\mathcal{T(\epsilon)} \neq -2, 0, 1, 2$ if $\epsilon = 0$. Since if (2.21) holds then $\mathcal{D}(0) = 1$, and hence $\mathcal{T}(0)\not = 2, -2$. Therefore, it further requires that $\mathcal{T}(0)\not = 1, 0$, that is, $\alpha\not = \frac{chd^2\left(hd-2\right)}{b\left(2hd-2-h^2d^2\right)}, \ \frac{chd^2\left(hd-1\right)}{b\left(h d-1-h^2 d^2\right)}$. Now, by the following transformation

$E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ transforms into $E_{trival}$ where $x^* = \frac{a\alpha}{cd+b\alpha}$ and $y^* = \frac{ad}{cd+b\alpha}$. In view of (3.13), (3.8) takes the following form:

Now, for $\epsilon = 0$, we will explore the normal form of (3.14). For this system (3.14), after expanding up to $3^{rd}$-order at $E_{trival}$, becomes

where

Now, using the following translation:

the system (3.14) implies

where

with

and

From (3.19), one gets:

Finally, for the occurrence of the indicated bifurcation, the following quantity is required to be non-zero:

where

From (3.24), the computation yields

From the above manipulation, we therefore say that if $\widehat{\psi}\not = 0$, then the Leslie's model (1.5) undergoes a N-S bifurcation. Additionally, the repelling (attracting) invariant closed curve bifurcates from $E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ if $\widehat{\psi} > 0$ ($\widehat{\psi} < 0)$, respectively. □

4.   Control of N-S bifurcation

In prey-predator Leslie models, controlling the Neimark-Sacker bifurcation is essential since it enables us to stabilize the dynamics of the populations. A sort of bifurcation called a Neimark-Sacker bifurcation occurs when a stable limit cycle, which is a population dynamics repeating pattern, loses stability and gives rise to a new, unstable limit cycle. Chaos in dynamics can result from this, which can be challenging to forecast and manage. In prey-predator Leslie models, we may prevent chaotic dynamics and stabilize the populations by managing the N-S bifurcation. As a result, we are able to better anticipate and lessen the consequences of environmental variables like habitat loss, climate change, and other disturbances. This is crucial for managing natural resources and protecting bio-diversity. In prey-predator Leslie models, the N-S bifurcation needs to be controlled in order to comprehend long-term ecosystem dynamics and create management plans for natural resources. We may learn more about the intricate interactions between many species in an ecosystem, establish plans for biodiversity preservation, and manage resources by researching the interactions between predator and prey populations. Furthermore, by limiting the N-S bifurcation, we can forecast how ecological changes in populations will impact it over time and create plans to reduce the impact of disturbances. Therefore, hereafter we will control the N-S bifurcation by the hybrid control method motivated from existing literatures ^[32,33,34]. If $\beta$ is a bifurcation parameter, then the Leslie's model (1.5) becomes

where $E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ is the interior fixed point. Moreover, the auxiliary equation of $\Lambda\biggr|_{E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)}$ at $E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ for controlled system, which is depicted in (4.1), is

where

The roots of (4.2) with $\beta$ are

Recall that if $\left(2-\beta hd-\frac{\beta abh\alpha}{cd+b\alpha}\right)^2-4\left(\frac{{a\left(\beta^{2} bdh^{2}\alpha +{\beta^{2}ch^{2}d^{2} }-{\beta bh\alpha }\right)}}{cd+b\alpha}+{1- \beta hd}\right) < 0$, then the auxiliary equation, which is depicted in (4.2), has two conjugate roots with modulo $1$. Additionally, it is easy to establish that $E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ of the system (4.1) is non-hyperbolic, stable focus, and unstable if following conditions hold, respectively:

and

Recall that if $a$ is a bifurcation parameter, then the system (4.1) becomes

where characteristics roots are

Furthermore, for the system (4.1) to undergo N-S bifurcation, the following non-degenerate condition(s) should hold:

Furthermore, it is also require that $\lambda_{1, 2}^\upsilon\neq 1, \upsilon = 1, \cdots, 4$ which corresponds to $\mathcal{T(\beta)} \neq -2, 0, 1, 2$ if $\epsilon = 0$. This compute to get $\alpha \not = \frac{chd^2\left(\beta^2 hd-2\right)}{b\left(3\beta hd-2 - \beta^2h^2 d^2 -hd\right) }, \ \frac{chd^2\left(\beta^2hd -1\right)}{b\left(2\beta h d- 1-\beta^2h^2d^2 - h d\right)}$. Now, by (3.13), the system (4.1) becomes

with $x^* = \frac{(a^*+\epsilon)\alpha}{cd+b\alpha}$ and $y^* = \frac{(a^*+\epsilon)d}{cd+b\alpha}$. Now, the system (4.11) becomes

where

Now, using the following translation,

system (4.11) gives

where

with

and

From (4.16), one gets:

From (4.16), the computation yields

Finally, from view (4.20) and (3.23), one can summarize that if $\widehat{\psi}\not = 0$ as $\bigcap\in NSB\biggr|_{E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)}$, then the controlled Leslie's system (4.1) undergoes a N-S bifurcation. Additionally, the dynamical classifications of the controlled Leslie's system (4.1) $E^{+}_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ are given in Figure 4.

5.   Numerical simulations

If $c = 0.8, \alpha = 2.5, d = 1.4, b = 0.5, h = 0.9$, then from obtained parametric condition, which is shown in (2.21), one has $a = 1.9110701532081558$. Therefore, if $a < 1.9110701532081558$, then $E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ is a stable focus; change the dynamics if $a = 1.9110701532081558$ and becomes unstable if $a > 1.9110701532081558$. An illustration, let $a = 1.7 < 1.9110701532081558$, then from Figure 5a it is obtained that $E^+_{xy}\left(1.8565400843881856, 1.0396624472573839\right)$ is a stable focus. Moreover, if $a = 1.78$, $1.8$, $1.83$, $1.87$, $1.9 < 1.9110701532081558$ then Figure 5b–5f demonstrate that corresponding equilibrium is a stable focus. On the other hand, if $a = 1.92 > 1.9110701532081558$, then Figure 6a implies that $E^+_{xy}\left(2.0253164556962022, 1.1282700421940928\right)$ is an unstable focus and as a result, subsequent computations demonstrate that the model undergoes supercritical N-S bifurcation. If $a = 1.92 > 1.9110701532081558$, then from (2.17), (3.12) and (3.25), we have

and

Using (5.3) and (5.1) in (3.23), one has $\widehat{\psi} = -0.9756591603967086 < 0$. Hence, if $a = 1.92 > 1.9110701532081558$, then Leslie's model (1.5) undergoes a supercritical N-S bifurcation. Similarly, if $a = 1.9234$, $1.94$, $1.97$, $1.98$, $2.1 > 1.9110701532081558$, then the corresponding numerical values of $\widehat{\psi}$ are given in Table 1. These calculations demonstrate that Leslie's model (1.5) undergoes a supercritical N-S bifurcation for $a = 1.9234, 1.94, 1.97, 1.98, 2.1 > 1.9110701532081558$ (see Figure 6b–6f). Finally, N-S bifurcation diagrams are presented in Figure 7 alongside the associated Maximum Lyapunov exponent, demonstrating the accuracy of the results in Sections 2 and 3.

Finally, simulations are given for the controlled Leslie's model (4.1). If $b = 0.5, \alpha = 2.5, h = 0.9, c = 0.8, d = 1.4, a = 2.2$ and $\beta = 0.85$, then from (4.4), (4.10) and (4.20), we have

and

Using (5.6) and (5.4) in (3.23), we get $\widehat{\psi} = -0.9756591603967086 < 0$. This implies that the controlled model (4.1) undergoes a supercritical hopf bifurcation (see Figure 8). To conclude, N-S bifurcation diagrams are presented in Figure 9 alongside the associated Maximum Lyapunov exponent, demonstrating the accuracy of the results in Section 4.

6.   Conclusions

In the current study, we investigated the hybrid control and bifurcation analysis of a discrete Leslie's model (1.5). More precisely, for all model's parameters, it is shown that Leslie's model (1.5) has boundary and interior equilibria $E_{x0}\left(\frac{a}{b}, 0\right)$, $E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$, respectively. Additionally, we investigated the topological classifications at equilibria $E_{x0}\left(\frac{a}{b}, 0\right)$ and $E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$ whose main finding are given in Table 2. Moreover, we studied the existence of possible bifurcations at equilibria $E_{x0}\left(\frac{a}{b}, 0\right)$ and $E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$. For the fixed point $E_{x0}\left(\frac{a}{b}, 0\right)$, it is proved that Leslie's model (1.5) does not undergo a flip bifurcation if $\bigcap\in FB\biggr|_{E_{x0}(\frac{a}{b}, 0)}$; however, at $E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)$, Leslie's model (1.5) undergoes only a N-S bifurcation if $\bigcap\in NSB\biggr|_{E^+_{xy}\left(\frac{a\alpha}{cd+b\alpha}, \frac{ad}{cd+b\alpha}\right)}$. In Leslie models, the development of a Neimark-Sacker bifurcation can result in unexpected and chaotic population dynamics, making it challenging to predict the ecosystem's long-term dynamics. However, by comprehending the fundamental causes of N-S bifurcations, we may create plans for population stabilization and environmental impact reduction. Changes in birth and death rates, interactions between predator and prey populations, and environmental issues like habitat loss and temperature change are some of the underlying reasons that can cause N-S bifurcations in Leslie models. These elements have the potential to disrupt population dynamics, resulting in unpredictably erratic population growth. We can create plans for protecting bio-diversity and managing natural resources by researching these aspects. Furthermore, we have studied bifurcation analysis by bifurcation theory. Additionally, we have controlled the N-S bifurcation by employing a hybrid control method. Finally, numerical simulations are also presented to validate the results.

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.