Spatiotemporal dynamics of a diffusive predator-prey system incorporating social behavior

1.   Introduction

The prey predator model (denoted P-P model) is considered one of many powerful tools for predicting the evolution of species in nature ^[17,18,26]. In recent years, a well-known prey behavior dented by 'herd behavior' or 'social behavior' that been successfully modeled and analyzed. When the prey makes a group defense (which know by prey's social behavior), the predators cannot reach the prey located in the center of the herd, so the predator will hunt only the prey located in the frontier of the herd ^[31]. Therefore, the interaction happens only on the borders of the group. There are numerous approaches to determining the number of prey on the outer line of the group. The simplest one is to consider that the prey makes a group o a square shape, then the number of the prey will be proportional to the square root of the prey density ^[33]. For modeling this specific behavior we presume that the resources make a group with a square shape. Hence, the density of the resources in the frontier of the group is four times the square root of the number of the resources on the herd frontier. Hence, inspired by the manner of proposing the classical Holling I interaction functional we can deduce that the functional response that describes this conduct is $F(M, N) = \gamma N \sqrt{M}$, where $\gamma$ is the hunting rate, and $M$ is the density o the resources, $N$ is the density of the consumers. The manner of building this interaction functional is discussed in details through the paper ^[1]. Further, by changing the group structure or shape the density of the resources will be changed accordingly. As epitome if we consider that the resources make a group with a sphere shape it will influence the number of the resources in the outer bound where it will be $4^{\frac{1}{3}}\times 3^{\frac{2}{3}}\times \pi ^{\frac{1}{3}}\times M^{\frac{2}{3}}$. For generalizing the previous results, it is considered that in ^[3,37] that the functional response $F$ can be generalized by the one $F_2(M, N) = \gamma N M^{\alpha}$, where the parameter $0 < \alpha < 1$ indicates the structure of the resources group. In the case of the square or circle group shape this rate will become $\alpha = \frac{1}{2}$, and in the case of $\alpha = \frac{2}{3}$ we got the case of the sphere of cube group shape, which means that this functional response generalizes all the previous cases of the prey herd structure. Dealing with the resources that exhibit grouping behavior is not always easy for the consumers, where the from one herd to another, and from predator to another, handling with the prey in the outer bound of the group changes and takes different time (different handling time for the predator to handle with a prey) which is been investigated in the first time by Holling to propose the Holling $II$ interaction functional ^[17]. This point of view is applied for modeling the intermingling between the resources and the consumers in this case through the paper ^[6], where a new functional response is obtained $F_3(M, N) = \frac{\gamma N M^\alpha}{1+\gamma t_h M^\alpha}$, which summarizes all the intermingling functions and cases on interaction. This intermingling function was the subject of an investigation on many occasions for the purpose of modeling many behaviors in nature we mention a few ^[10,13]. In fact, using the functional response $F_3$ as intermingling function next to the logistic increasing of the prey and linear mortality of resources (see as an example the paper ^[6]), we get a Gause-type model ^[4], which means that the mathematical investigation is trivial and can be distinguished easily from the paper ^[4]. Investigating the herd behavior in mathematical models is the subject of the recent activities, we cite the researches ^{[2,5,7,8,9,10,12,13,32]}, and for more reading about different mathematical modeling of some natural phenomenon we cite the papers ^{[11,14,16,21,22,23,24,25,30]}. In this research, we will use a different approach where we will incorporate the Leslie-Gower intermingling functional with the functional $F_3(M, N)$ functional responses. Incorporating Leslie-Gower's functional response form with resources social behavior is a recent step and attracts any researchers, we cite for instance the papers ^{[15,18,27,28,35,36]}, hence, it is the subject of interest in this research.

Modeling the intermingling resource-consumer in the case of resources social depends on the spatial positioning, so, it is wise to consider a spatiotemporal model ^[34]. Our purpose in this research is to study the influence of the Leslie-Gower forme on the evolution of the two species, where we will use a comparative analysis for achieving this goal. In these regards, we consider the first model that models the intermingling resource-consumer in the case of the resources social conduct and no Leslie-Gower form. The investigated model is:

with $M(x, \tau)$ (resp. $N(x, \tau)$) is the number of the resources (rep. consumer) at $t$ and position $x$, $\beta M(x, \tau)\left(1-\frac{M(x, \tau)}{L}\right)$ is the logistic increasing of the resources with increasing rate $\beta$ and the crying capacity of the environment $L$ for the resources, $\mu$ in the mortality coefficient for the consumer, $e$ is the conversion rate of the resources into consumer, $\delta$ (resp. $\eta$) is the diffusion rate for the resources (resp. consumer). To mention the Neumann boundary conditions highlights that neither the resources or the consumers cannot move cross the borders. Our main contribution consists to cooperate the Leslie-Gower with the prey social behavior interaction function and determine its effect on the temporal behavior of the solutions. The investigated model is given as:

Our purpose is to investigate with the influence of the new approximation provided in (1.2) in modeling the interaction resources-consumer in the presence of the resources grouping behavior. In fact, we will use a comparative analysis between the two models (1.1) and (1.2), where we will show that the system (1.2) have a very rich dynamics. For achieving these aims we use the sections:

In Sec. 2 we analyze (1.1), where it is proved that it can undergo Hopf bifurcation (H-bifurcation), and cannot have Turing instability (T-instability), which means that it is not possible to have Turing-Hopf bifurcation (T-H bifurcation). The third section is used to analyze the system (1.2), where we will show that the system (1.2) undergoes many types of bifurcation as T-bifurcation, H-bifurcation, and T-H bifurcation. However, the system (1.1) can undergo only Hopf bifurcation (there is no T-bifurcation, and then there is no T-H bifurcation). The normal form of T-H bifurcation is utilized for study of steady states solution near the T-H bifurcation. The obtained mathematical results are confirmed using numerical simulation.

2.   Long time behavior of the model (1.1)

We know that the equilibrium states for the system (1.1) are solutions of the following system

Clearly, (2.1) has three equilibria $E_{0}(0, 0), \ E_{1}(L, 0)$, and the unique positive equilibrium point $E_{\ast}(M_{\ast}, N_{\ast})$ where

which exists if the following conditions hold

The linearized system of (1.1) evaluated at $(M_{\ast}, N_{\ast})$ is

where $D\Delta = diag(\delta\frac{\partial ^{2}}{\partial x^{2}}, \eta\frac{\partial ^{2}}{\partial x^{2}})$ and

Then the characteristic equation of system (2.3) takes the following form

where

and

with

Using the fact that

it follows that (2.6) becomes

and

Putting

Clearly, for $k = 0$ we have $Det_{0} > 0$. Then, if $L > L_{H}$ we have $Tr_{0} < 0$ which means that $E^{\ast}$ is locally stable and if $L < L_{H}$ we get $Tr_{0} < 0$, then $E^{\ast}$ is unstable. Now, taking $L$ as the bifurcation parameter. Thus, we have

Theorem 2.1. Presume that $(H_{1})$ holds and we put

Then, there exists an integer $k^{\ast}$ for which the model (1.1) undergoes a H-bifurcation at $E^{\ast}$ when $L = L_{k}$ for $0\leq k\leq k^{\ast}$. Further, if $k = 0$ the periodic solution is homogeneous (spatially), and if $k = 1, ..., k^*$ the periodic solution is nonhomogeneous (spatially) when $k = 1, ..., k^*$, where $[.]$ is the integer part function.

Proof. Denotes

Clearly, for an integer $k < \bar{k}$ we have $L_{k} > 0$. We consider that $Det_{k}$ defined by (3.10) is a function of $k$. Clearly, $Det_{0} > 0$, therefore there exists $\tilde{k} > 0$ (it can be $+\infty$ if $Det_{k} > 0$ for all integer $k = 0, 1, 2, ...$) a positive integer which represents the first integer that satisfy $Det_{\tilde{k}} > 0$ and $Det_{\tilde{k}+1}\leq0$. Taking $k^{\ast} = \min \{\bar{k}, \tilde{k}\}$. In this case, we guarantees that for $L = L_{k}$, $k = 0, 1, ..., k^{\ast}$ we have $Tr_{k} = 0$ and $Det_{k} > 0$ for all $k = 0, 1, ..., k^{\ast}$ which implies that (2.5) has purely imaginary roots. Letting

be the solution of Eq (2.5) verifying

Then, we get

It follows that the transversal condition is justified at each $L_{k}$, $k = 0, 1, ..., k^{\ast}$, that is to say that system (1.1) undergoes H-bifurcation at $L = L_{k}$. □

Now, we focus on proving that the system (1.1) cannot exhibit the existence of Turing instability. This phenomena is know by the diffusion driven instability, where considering the presence of spatial diffusion (with distinct diffusion rates) can destabilize a stable equilibrium, which means that if an equilibrium is stable in the absence of the spatial diffusion can become instable in the presence of the diffusion (with distinct dispersal rates).

Lemma 2.2. The model (1.1) cannot have T-instability at $E^{\ast}$.

Proof. Before proceeding to prove Lemma 3.6, it is necessary to assume that the conditions that guarantees the existence and the stability of $E^*$ in the absence of diffusion holds which consists to suppose that $(H_{1})$ holds and $L < L_{H}$. Immediately, we find that

Consequently, it follows by $Det_{k} > 0, \ \forall k\in \mathbb{N}$ the non occurrence of T-instability. □

3.   Long time behavior of the model (1.2)

In this section, we consider the diffusive predator-prey model with Leslie-Gower term (1.2). We discuss the existence of Hopf bifurcation, after that we derive the condition for Turing pattern which leads to the occurrence of T-H bifurcation. First of all, we prove that system (1.2) has unique positive solution

3.1.   Existence and boundedness

Here, we investigate the existence and positivity of solution for (1.2).

Theorem 3.1. Assume that $\beta, L, \gamma, \tau_{h}, \delta$ and $\eta$ are all positive, if $M_{0}(x, \tau)\geq 0$ and $N_{0}(x, \tau)\geq 0$ for $(x, \tau) \in [0, l\pi] \times [0, +\infty)$. Hence, (1.2) has a unique positive solution verifying

such that $(M^{**}(\tau), N^{**}(\tau))$ is the unique solution of

Proof. Putting

Since, $f_{N}\leq 0$ and $g_{M}\geq 0$ for $(M, N)\in \mathbb{R}_{+}^{2} = \{(M, N)\vert M\geq 0\}$ and from ^[29] yields that $f, g$ are mixed quasi-monotone functionals in $\mathbb{R}$. Now, letting

From

and $0\leq M_{0}(x, \tau)\leq M_{0}^{**}, 0\leq N_{0}(x, \tau)\leq > N_{0}^{**}, \ \text{for} \ (x, \tau) \in [0, l\pi] \times [0, +\infty],$ we can conclude that $(\tilde{M}, \tilde{N})$ and $(\hat{M}, \hat{N})$ are the upper and lower solutions of (1.2), respectively. From ^[29], we get the existence of a unique globally defined solution $(M(x, \tau), N(x, \tau))$ of system (1.2) satisfying $0\leq M(x, \tau)\leq M^{**}(\tau) \ \text{and} \ 0\leq N(x, \tau)\leq N^{**}(\tau).$

Applying the strong maximum principle to (1.2), we obtain that $M(x, \tau) > 0, \ N(x, \tau) > 0$ for $(x, \tau) \in (0, l\pi) \times (0, +\infty)$. □

Now, we mainly focus on checking the existence of feasible steady states. We can easily verify that $E_{1} = (L, 0)$ are always equilibrium for (1.2). Next, we investigate the existence of a positive steady state of our proposed system (1.2) which is denoted by $E_{2} = (M^{\ast}, N^{\ast})$, then we have the following theorem

Theorem 3.2. Suppose that all parameters are positive, then the system (1.2) has a unique positive steady state $E_{2} = (M^{\ast}, N^{\ast})$ with $M^{\ast} = N^{\ast}$ and $0 < M^{\ast} < L$ (please see Figure 1).

Proof. Obviously, $E_{2}$ is the solution of:

Using the second equation of (3.2), we can notice that $M^{\ast} = N^{\ast}$. Substituting $M^{\ast} = N^{\ast}$ into the first equation of (3.2), gives $h(M^{\ast}) = 0$, with

and

Obviously, $h(0) = -\beta < 0$ and $h(L) = \frac{\gamma L^{\alpha}}{1+\gamma \tau_{h}L^{\alpha}} > 0$. Hence, $h(M)$ is an increasing function for $M\in [0, L]$, and bisect the horizontal axis at $M^{\ast}$, where $0 < M^{\ast} < L$. The proof of Theorem 3.2 is completed. □

Now, we define the Sobolev space

and its complexification, $X_{\mathbb{C}} : = X \oplus iX = \left\{x_{1}+ix_{2}\vert x_{1}, x_{1} \in X\right\},$ with the the inner product $\langle., .\rangle$ as

The linearization of (1.2) at $(M, N)$ is

where $D\Delta = diag(\delta\frac{\partial ^{2}}{\partial x^{2}}, \eta\frac{\partial ^{2}}{\partial x^{2}})$ and

Clearly, the problem

has the eigenvalues $\xi_{k} = \left(\frac{k}{l}\right)^{2}(k = 0, 1, 2, ...)$ where the corresponding normalized eigenfunctions defining in the Sobolev space $X$ are

Now, let

be the nontrivial solution of the system (3.4) yields the existence of $k\in \mathbb{N}_{0}$ such that $\lambda$ satisfy

where, $I$ is $2\times 2$ identity matrix, and $C_{k} = -\left(\frac{k}{l}\right)^{2}diag(\delta, \eta)$. After a straightforward calculation we obtain the characteristic equation of (3.4) as

where

and

For the semi trivial steady state $E_{1}$ we obtain

Obviously, we have $Det_{0}\vert_{E_{1}} = -\sigma \beta < 0$ in the absence of the spatial diffusion, which means that $E_{1}$ is always unstable.

Now, let's concentrate on analyzing the stability and the bifurcation properties of $E_2$. From Theorem 3.2, (1.2) has a unique equilibrium denoted by $E_{2} = (M^{\ast}, N^{\ast})$ where, $M^{\ast} = N^{\ast}$ and $\ 0 < M^{\ast} < L$. Submitting $(M^{\ast}, N^{\ast})$ into (3.9) and (3.10) and using the fact that

then, we obtain

where

and

Now, putting

Remark 3.3. If $\frac{L}{2} < M^{\ast} < L, \ i.e. \ \beta^{\ast} < 0$, we obtain that $Tr_{k}(\sigma) < 0$ and $Det_{k}(\sigma) > 0$, so $E_{2}$ is always locally asymptotically stable (with or without diffusion). Thus, there is no Hopf bifurcation.

Now, assume that the following condition holds

So, we have the following theorem

Theorem 3.4. If $\eta = \delta = 0$, then:

(i) The positive steady state $E_{2}$ of system (1.2) is locally asymptotically stable if $\sigma > \sigma_{0}$ and unstable if $\sigma < \sigma_{0}$ (see Figure 2).

(ii) The local system of the diffusion system (1.2) undergoes Hopf bifurcation for $\sigma = \sigma_{0}$, where $\sigma_{0}$ is defined in (3.15) (please see Figure 3).

Proof. Under the condition $(\mathbf{H_{2}})$ and from (3.13), if $k = 0$ one can immediately get that $E_{2}$ is locally asymptotically stable when $\sigma > \sigma_{0}$ Figure 2 and unstable when $\sigma < \sigma_{0}$. For $\sigma = \sigma_{0}$, Eq (3.8) has a pair of purely imaginary roots $\pm i\sqrt{\sigma C_{\ast}}$. Let $\lambda(\sigma) = \upsilon(\sigma)\pm i\omega(\sigma)$ be the root of (3.8) satisfying $\upsilon(\sigma_{0}) = 0, \ \omega(\sigma_{0}) = \pm i\sqrt{\sigma C_{\ast}}$. A simple calculate gives following transversality condition:

Therefore, we conclude that the non diffusive system associated with (1.2) undergoes H-bifurcation at $E_{2}$ when $\sigma = \sigma_{0}$ (see Figure 3). □

3.2.   H-bifurcation

Now, concentrating on studying the occurrence of time-periodic solutions for (1.2) generated by H-bifurcation. throughout the rest of this paper, we assume that $(\mathbf{H_{2}})$ holds i.e. $\sigma_{0} > 0$ and taking $\sigma$ as the bifurcation parameter. Recall that H-bifurcation appears if

Obviously

where

are the critical points for H-bifurcation. The values of the H-bifurcation are highlighted as follows.

Theorem 3.5.

(i) The diffusive model (1.2) undergoes H-bifurcation at $E_{2}$ when $\sigma = \sigma_{k}$, for $k = 0, 1, ..., k^{\ast}-1$, where $\sigma_{k}$ and $k^{\ast}$ are defined in (3.18) and (3.19), respectively (see Figure 4). Further, for $k = 0$ we get a homogeneous periodic solution and a non homogeneous periodic solution for $k = 1, 2, ..., k^{\ast}-1$.

(ii) The eventual H-bifurcation points $(\sigma_{k})_{0\leq k\leq k^{\ast}-1}$ satisfying the following relationships

Proof.

${\textbf{(i)}}$ From the definition of the integer $k^{\ast}$, we can easily affirm that when $\sigma = \sigma_{k}$, $Tr_{k}(\sigma_{k}) = 0$ and $Det_{k}(\sigma_{k}) > 0$ for $k = 0, 1, ..., k^{\ast}-1$, which follows purely imaginary roots of Eq (3.8). Letting

be the roots of Eq (3.8) which verifies

It follows that if $\sigma$ is in the neighborhood $\sigma_{k}$, the solutions of the characteristic equation (3.8) take the following form

with

and we have

This yields to the verification of the transversality condition for each $\sigma_{k}$ where $k = 0, 1, ..., k^{\ast}-1$.

${\textbf{(ii)}}$ Now, we are in the position to prove the second affirmation. After a simple calculate we find

This implies that $(\sigma_{k})$ is strictly decreasing sequence for all $k = 0, 1, ..., k^{\ast}-1$. The proof of Theorem 3.5 is completed.

□

3.3.   T-instability

In the following, we mainly prove that under certain sufficient condition, the system (1.2) may exhibits T-instability unlike the situation in the system without Leslie-Gower term (please see Sec. 2). Mentioning that T-instability holds when the equilibrium point is locally stable in the non diffusive system and becomes unstable in the case of diffusive system (i.e. $Det_{k} < 0$ for some value of integer $k$). Notice that $E_{2}$ is locally stable if the condition $(\mathbf{H_{2}})$ is satisfied and $\sigma > \sigma_{0}$. In this case we have

In order to study the occurrence of the Turing instability, we define the functional $\Theta$ as

which considered as a quadratic polynomial in $\left(\tfrac{k}{l}\right)^{2}$ and $\sigma_{0}$ is defined in (3.15).

Lemma 3.6. If

the system (1.2) cannot undergo T-instability.

Proof. Clearly, under the condition (3.20) we have $\Theta\left(\left(\tfrac{k}{l}\right)^{2}\right) > 0$ which means that system (1.2) has no diffusion driven instability. □

In the next, we presume that

Obviously, if $\Theta\left(\left(\tfrac{k}{l}\right)^{2}\right) < 0$, then Eq (3.8) has one of the two roots is positive. If $F(\eta, \delta) = \eta \sigma_{0}-\sigma \delta > 0$, $\Theta\left(\left(\frac{k}{l}\right)\right)$ has a minimum at

Evaluating $\Theta$ at this minimum, we get

where $C^*$ is defined by (3.14). Next we show that under $\mathbf{(H_3)}$, $\min\limits_{\left(\frac{k}{l}\right)^{2}}\Theta\left(\left(\frac{k}{l}\right)^{2}\right) < 0$ for some values $\eta/\delta > 0$, which it known by the condition of T-instability. Defining the ratio $\xi = \eta/\delta$ and

Then, $\Pi(\eta, \delta) = 0$ and $F(\eta, \delta) = 0$ are equivalent to

and

Notice that

which means that Eq (3.22) has two has two positive real roots

Easily, one can see that $0 < \xi_{2} < \xi^{\ast} < \xi_{1}$ and if $\eta/\delta > \xi_{1}$, we have $\min\limits_{\left(\frac{k}{l}\right)^{2}}\Theta\left(\left(\frac{k}{l}\right)^{2}\right) < 0$ and $F(\eta, \delta) > 0$. Here, the positive equilibrium $E_{2}$ becomes unstable, which means that T-instability occurs.

Now, defining the set

Hence, we get the results

Theorem 3.7. Presume that $(\mathbf{H_{2}})$ holds and $\sigma > \sigma_{0}$ (for having the stability of the positive equilibrium). Then there exists an unbounded set $R_{T}$ such that for any $(\eta, \delta)\in R_{T}$, the equilibrium $E_{2}$ becomes unstable, that is, Turing instability (for illustrations we refer Figure 5).

3.4.   T-H bifurcation

In this subsection, our aim is to investigate the occurrence of T-H bifurcation. This type of bifurcation happen if there are two integers $k_1$ and $k_2$ where for $k = k_1$, (1.2) has H-bifurcation and for $k = k_2$ the system (1.2) undergoes T-bifurcation, this kind of bifurcation is a bi-dimensional bifurcation which means that we need to choose two bifurcation parameters. Hence, we choose $\sigma$ and $\delta$ as bifurcation parameters. Assuming that $\sigma_{0} > 0$, it is well known that $Tr_{0}(\sigma) = 0$ and $Det_{0}(\sigma) < 0$ are necessary conditions for Hopf bifurcation to occur. From (3.13) and (3.15), $Tr_{0} = 0$ is equivalent to

which is the line of H-bifurcation in $\delta-\sigma$ plan, where the frequency of the oscillations is

Besides, T-instability occurs when $\frac{\eta}{\delta} > \xi_{1}$, where $\xi$ is given in (3.24). It follows that the critical value of the T-bifurcation for the parameter $\sigma$ is

Now, we prove the existence of intersection point between the H-bifurcation curve $\sigma^{H}$ and T-instability curve $\sigma^{T}$ in $\delta-\sigma$ plane. Defining the following function

Clearly, $h$ is monotonously decreases with the increasing of $x$. In addition, we have $\underset{x \rightarrow 0^{+}}{\lim}h(x) = +\infty$. Therefore, we conclude that the H-bifurcation line $\sigma^{H}$ cuts the T-bifurcation curve $\sigma^{T}$ at $(\delta^{T-H}, \sigma^{T-H}) = (\hat{\delta}, \sigma_{0})$ where

Now, we will examine the transversality condition. Fixing $\delta$, and taking $\sigma$ as parameter, letting $\lambda(\sigma)$ the roots of (3.8), hence:

then, we get

Theorem 3.8. Presume that the conditions $(\mathbf{H_{2}})$ and $(\mathbf{H_{3}})$ are verified, then:

${\textbf{(i)}}$ The H-bifurcation line $\sigma^{H}$ cuts the T-bifurcation curve $\sigma^{T}$ in $\delta-\sigma$-parameter space at the unique point $(\hat{\delta}, \sigma_{0})$, where $\sigma_{0}$ and $\hat{\delta}$ are defined in (3.15) and (3.27) (for illustrations we refer Figure 6).

${\textbf{(ii)}}$ At $(\delta, \sigma) = (\hat{\delta}, \sigma_{0})$ the characteristic equation (3.8) has a simple zero root.

In order to illustrate numerically the obtained result in Theorem 3.8, T-bifurcation curve and H-bifurcation curves are plotted in $\delta-\sigma$ plane (please see Figure 5). we fix the parameters $\beta = 0.9, \ L = 10, \ \gamma = 0.2, \ \tau_{h} = 0.5, \ \alpha = 2/3, \ \sigma = 0.5$ and $\eta = 1$. It follows that $M^{\ast} = N^{\ast} = 4.97 < L = 10, \ \sigma_{0} = 0.4478 < \sigma = 0.5, \ \alpha^{\ast} = 0.432 < \alpha = 0.667, \ C_{\ast} = 5.6316$ and $\hat{\delta} = 0.049$. From Figure 6, we can see that the intersection point (1) divide the $\delta-\sigma$ plan into four regions. In $D_{1}$, $E_2$ is stable. $D_{2}$ represents the T-bifurcation region. $D_{3}$ is the domain in which the pure H-bifurcation occurs. In $D_{4}$, both T-instability and H-bifurcation occur. In this situation, the diffusive system (1.2) produces a complex spatiotemporal patterns, where the two instabilities T- bifurcation and H-bifurcation coincide.

4.   Normal forms for T-H bifurcation

Here, we mainly focus on the calculate of the normal form of T-H bifurcation in order to determine the properties and the spatiotemporal dynamics of (1.2) at $E_2 = (M^{\ast}, N^{\ast})$ near the T-H bifurcation point $(\hat{\delta}, \sigma_{0})$. At first, we set $\mu = (\mu_{1}, \mu_{2})\in \mathbb{R}^{2})$ where $\mu_{1} = \delta-\hat{\delta}, \ \mu_{2} = \sigma-\sigma_{0}$. Then, we we apply the translation $\overline{M} = M-M^{\ast}, \ \overline{N} = N-N^{\ast}$ to (1.2) and introduce a new parameter $\mu = (\mu_{1}, \mu_{2})\in \mathbb{R}^{2})$. we denote $\overline{M}$ by $M$ and $\overline{N}$ by $N$. Therefore, the diffusive system (1.2) becomes

For system (4.1) and according to ^[19], also we get

and

where

It follows that

and

where

The coefficients $a_{ij}$ and $b_{ij}$ are given as

Now, we get the corresponding characteristic matrices as

Clearly, $\lambda = \pm i\omega$ with $\omega = \sqrt{Det_{0}}$, are eigenvalues of $\mathbb{D}_{0}(\lambda)$, and $\lambda = 0$ is a simple eigenvalues for $\mathbb{D}_{\hat{\delta}}(\lambda)$, while other eigenvalues have negative real parts. From Theorem 3.8 and by using a simple calculate we can obtain

and

Then, by the procedure developed in ^[19,20], the normal form restricted on central manifold at T-H bifurcation singularity is

where the calculation of $m_{1}(\mu), \ m_{200}, \ m_{011}, \ m_{300}, \ m_{111}, \ n_{2}(\mu), \ n_{110}, \ n_{210}, \ n_{021}$ are given in "Appendix".

By using the new parameter transformation $Z_{1} = r, \ Z_{2} = \rho\cos \vartheta-i\rho\sin \vartheta$, then we get

5.   Numerical simulation

Here, we will provide some figure for illustrating the obtained results. In fact we investigate the following cases:

Figure 7: In this figure we set $\beta = 0.5$, $L = 15$, $\gamma = 1.5$, $\alpha = 0.7$, $\tau_h = 1.7$, $\sigma = 3.1$, $\delta = 0.02$, $\eta = 0.04$, $l = 1$ and the initial data $M(x, 0) = 0.5+0.1\cos(3x)$, $N(x, 0) = 0.5+0.1\cos(2x)$. Here, we obtain the stability of the nonhomogeneous steady state.

Figure 8: In this figure we consider the following values $\beta = 3.1$, $L = 50$, $\gamma = 1.5$, $\alpha = 0.66$, $\tau_h = 0.01$, $\sigma = 1.01$, $\delta = 0.01$, $\eta = 0.04$, $l = 1$ and the initial data $M(x, 0) = 0.5+0.1\cos(3x)$, $N(x, 0) = 0.5+0.1\cos(2x)$. Here, we arrive to the stability of the nonhomogeneous steady state.

Figure 9: In this graphical representation we set $\beta = 1.1$, $L = 50$, $\gamma = 1.5$, $\alpha = 0.66$, $\tau_h = 0.01$, $\sigma = 0.51$, $\delta = 0.01$, $\eta = 0.04$, $l = 1$ and the data $M(x, 0) = 0.5+0.1\cos(3x)$, $N(x, 0) = 0.5+0.1\cos(2x)$. We arrive at stability of the nonhomogeneous steady state.

Figure 10: Here we choose the set of values $\beta = 1.1$, $L = 50$, $\gamma = 1.5$, $\alpha = 0.66$, $\tau_h = 0.01$, $\sigma = 0.1$, $\delta = 0.01$, $\eta = 0.04$, $l = 1$ and the initial data $M(x, 0) = 0.5+0.1\cos(3x)$, $N(x, 0) = 0.5+0.1\cos(2x)$. Here, we arrive to stability of nonhomogeneous periodic solutions.

Figure 11: Here we choose the set of values $\beta = 1.51$, $L = 10$, $\gamma = 1.5$, $\alpha = 0.66$, $\tau_h = 0.01$, $\sigma = 0.051$, $\delta = 0.01$, $\eta = 0.04$, $l = 1$ and the data $M(x, 0) = 0.5+0.1\cos(3x)$, $N(x, 0) = 0.5+0.1\cos(2x)$. Here, we arrive to stability of nonhomogeneous periodic solutions.

Figure 12: Here, we choose the set of values values $\beta = 1.51$, $L = 10$, $\gamma = 1.5$, $\alpha = 0.66$, $\tau_h = 0.01$, $\sigma = 0.08$, $\delta = 0.09$, $\eta = 0.01$, $l = 1$ and the data $M(x, 0) = 0.5+0.1\cos(3x)$, $N(x, 0) = 0.5+0.1\cos(2x)$. Here, we arrive to stability of nonhomogeneous periodic solutions.

6.   Discussion

In this research, we investigated a spatiotemporal P-P model with Leslie-Gower for modeling the saturation of the predator increasing in terms of the density of the prey. The reason behind considering such approximation is to highlight that the evolution of the consumers is affected directly by the density of the resources. The similarity points between the two models, in the absence of the Leslie-Gower scheme interaction functional (1.1), and the presence of this last (1.2), as the occurrence of H-bifurcation in the absence and the presence of diffusion in two studied models. The disagreement between the two considered models consists of the presence of T-instability for (1.1) and the existence of this last in the diffusive model (1.2). Our study was focused on distinguishing the influence of this case of interaction on the value of T-H bifurcation. As it is been highlighted in Figure 4, $\sigma$ (increasing rate for predator) generated by considering the Leslie-Gower scheme interaction functional has a big effect on the existence of T-patterns, and hence it influences the existence of T-H bifurcation. No one can neglect the role of herd behavior in modeling much natural behavior, and the considered model can fit many cases in different species as fish population (sardines), gnus and buffalos which intermingle with the lions, hyenas, which highlights the importance of considering such as approximation.

In fact, there are many scenarios that can behold as the persistence of the two categories with nonhomogeneous patches as it is been shown in Figures 7 and 8, or in nonhomogeneous and periodic patters as the Figures 10–12. These scenarios generated by the presence of the Leslie-Gower scheme functional response (more precisely the T-H bifurcation), which shows the huge importance of considering such as approximation.

Acknowledgments

This study is supported via funding from Prince Sattam bin Abdulaziz University project number (PSAU/2023/R/1444). The author S. Djilali is partially supported by the PRFU project No.C00L03UN130120200004, DGRSDT, Algeria.

Conflict of interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Appendix

Calculations of $m_{1}(\mu), \ m_{200}, \ m_{011}, \ m_{300}, \ m_{111}, \ n_{2}(\mu), \ n_{110}, \ n_{210}, \ n_{021}$. Here, we are in the position to give the expressions of $m_{1}(\mu), \ m_{200}, \ m_{011}, \ m_{300}, \ m_{111}, \ n_{2}(\mu), \ n_{110}, \ n_{210}, \ n_{021}$. We will put only the formulas of these parameters. For more details about the method of calculation we refer the authors to ^[19,20].

where

and