Bifurcation structure of nonconstant positive steady states for a diffusive predator-prey model

1.   Introduction

In the middle of the 20th century, the disease of Poliomyelitis (polio) was one of the contagious diseases which has been transmitted from the western regions of the Asian continent. Jonas Salk, in the year 1952 developed a vaccination process for this outbreak. Poliovirus became the cause of this outbreak among most societies. It affects the nerves of the human body and may cause whole-body paralysis within a short duration. The associated disease mainly affects children in the first three years of age but it may also affect anyone at any age. The poliovirus is transmitted to the children's bodies through the mouth, eventually reaching the intestines of the digestive system. The poliovirus can be transferred from one person to the other. The starting signs of the associated infection are fever, headache, fatigue, vomiting, and pain produced in the limbs. There is still no effective treatment for the virus, but its spread can be controlled through the process of vaccination and immunizations which commonly save many children's lives.

Mostly, assumptions regarding dynamic problems of daily life are necessary to manage these problems; therefore, in every mathematical model the infections are only transmitted to the susceptible or healthy class due to the infective population. Several diseases like polio are contagious during the time of incubation. So one can conclude that the interactions between the healthy population classes exposed classes have their importance in the spread of diseases. Thus, scholars have discussed the relationship between the susceptible and the exposed classes. Many scientists have attempted to construct mathematical models of Polio. The authors [1] developed a mathematical system spread of polio after the re-introduction of poliovirus. Garfinkel and Sarewitz et al. [2] briefly studied the conditions that would lead to polio eradication. Furthermore, to prove efficacy of polio vaccines in the control of polio infection, most scholars have constructed various models for the spread of polio [3]–[5]. The authors of [6] specifically considered the SIR problem with pulse vaccination and have shown that pulse vaccination can lead to eradication under specific conditions related to the magnitude of vaccination proportion. In some of the articles like [7], the authors dealt with a two-dimensional SIS polio system with the addition of vaccination terms for bifurcation in the opposite direction. The authors of [8] developed a polio mathematical system to investigate its dynamical behavior in the presence of vaccinations strategy with nonlinear incidence, which also deals with the optimality of the model that converges the vaccine quantity to the threshold value. Agarwal and Bhadauria [9] developed a polio epidemic model with vaccine quantity to investigate the effect of vaccination once a strategy is applied to the susceptible and exposable populations. In this study, the polio model was adapted from [9], and it has been divided into five subgroups, namely susceptible persons $\mathcal{S}$, exposed persons $\mathcal{E}$, acutely infected persons $\mathcal{I}$, chronically infected persons $\mathcal{V}$ and recovered persons R:

with the following initial conditions

The parameters used in model (1.1) are described as follows: A is the rate of immigration in the society, µ is the rate of death that occurs naturally, β is the probability of disease transmission because of the infectious population, rβ is the rate of chance of infection because of the exposed population, r is the rate of reduction in the transmission of infection by exposed class, ν is the rate of individuals moving from the susceptible to vaccinated class, ν₁ is the rate at which exposed individuals are vaccinated, b is the rate at which exposed individuals move to the infective class and α is the rate of death due to the disease.

The utilization of fractional calculus in the fractional-order model, which combines differentiation and integration, offers a more effective approach to comprehending real-world problems as compared to classical derivatives [10]–[16]. The concept of fractional derivatives, originally introduced by Riemann-Liouville based on the power law, has recently faced scrutiny from researchers due to the presentation of experimental results questioning the adequacy of a single fractional operator, such as the Caputo, Caputo-Fabrizio and Atangana-Baleanu operators, as a tool to describe complex phenomena in science and engineering [17]–[21]. These findings suggest that a broader range of fractional operators may be necessary to accurately capture the behavior of diverse systems [22]–[26]. The authors of [27] developed a model by using a new definition of the constant proportional Caputo operator, which describes the generalized memory effects. A new recursive algorithm has been constructed to solve certain initial value problems involving a fractional differential equation [28]. This new approach Adomian decomposition method (ADM) is based on the application of the ADM, and it involves combining the decomposition with a recurrence formula and utilizing the solutions of the generalized Abel equation. According to the researchers, a new fractional derivative can be defined by utilizing the exponential kernel [29], [30]. The utilization of non-singular kernel fractional derivatives in the modeling of epidemics offers valuable insights into the dynamics of infectious disease outbreaks, particularly when applied to trigonometric and exponential functions [31]–[33]. Initialization can pose challenges when working with non-singular kernel fractional derivatives, especially in cases involving non-singular kernels. Regardless of their type, equations in the following form have been observed to present similar challenges, as noted by the authors of [34]:

The condition F(a,a)≠0 introduces distinct limitations on differential equations involving non-singular kernels, implying that if F(x,τ) and y(x) are continuous and F(a,a)≠0, then y(a)≠0. This condition gives rise to peculiar behaviors in these equations, necessitating careful consideration during the development of models in various research fields. As mentioned in [34], this issue remains unsolved. However, the problem associated with the operator introduced in [17] has been solved by various researchers [35], [36]. To overcome the challenges associated with non-singular operators, the authors of [37] proposed a modification to the operator that involves utilizing a Mittag-Leffler kernel. They demonstrated that the resulting fractional differential equations, based on this modified operator, are easier to initialize than those based on non-singular kernels. Additionally, they showcased that the modified operator, known as the modified Atangana-Baleanu-Caputo (mABC) derivative, is capable of solving multiple fractional differential problems, which is not possible for the Atangana-Baleanu-Caputo (ABC) derivative. The mABC derivative possesses an integrable singularity at the origin. In [38], the authors developed a unique numerical method for the mABC derivative, utilizing finite differences; the method streamlines the initialization of the corresponding fractional differential equations. Researchers have also developed modified fractional difference operators by employing the mABC derivative and Mittag-Leffler kernels [39].

We propose a memory-affected mABC model to model the spread of polio. Although the epidemic model described by Eq (1.1) utilizes classical derivatives and should be considered, our memory-affected mABC model provides additional insights and enhances the accuracy of hepatitis B spread predictions. To convert the ordinary system of equations given by Eq (1.1) into the mABC operator, we make the following transformation:

with the following initial conditions

The rest of the paper is arranged as follows. In Section 2, we recall some basic definitions and statements from the literature. In Section 3, we give some formulas and show the non-zero solution to the homogeneous fractional initial value problem. With the use of the Laplace Adomian decomposition method, we find the approximate solution for the mABC derivative in Section 4. In the same section, we also establish the stability and uniqueness results for the considered system. The approximate series solution is graphically presented in Section 5. We conclude our work in Section 6.

2.   Basic results

Here, we recall basic results from the literature on fractional calculus.

Definition 2.1. [37] Let f(t)∈L1(0,T) be a function; then, the mABC derivation is presented as follows:

where Eθ is known as the Mittag-Leffler function for one parameter while Eθ,θ denotes the Mittag-Leffler for two parameters; from the above definition one may prove that $^\mathrm{mABC}\mathbf{D}_t^\vartheta\mathrm{f}=0.

Definition 2.2. [37] Let f(t)∈L1(0,T) be a function, then, the integral of mABC operator is given as

Lemma 2.2.1. For f′∈L1(0,∞) and the order ϑ∈(0,1), we obtain

The Laplace transform for the mABC is defined as

3.   Analysis of the model

We state that there is a non-zero solution to the homogeneous fractional initial value problem. To achieve that, we use the following formulas:

Lemma 3.0.1. [37] Suppose that the fractional initial value problem is as follows:

where 0<ϑ<1.

(1) For Ω=M(ϑ)1−ϑ, the solution is as follows:

(2). For Ω≠M(ϑ)1−ϑ, the solution is as follows:

where, τϑ=Ω(1−ϑ)M(ϑ).

Proof. (1) Given that

for 0<t, we have

which shows that the proof is complete.

Applying Eq.(3.1) and Eq.(3.2) for t>0, we have

Thus, the proof is complete.

Lemma 3.0.2. [37] Assume the following fractional differential equations

where 0<ϑ<1 and Ω≠−M(ϑ)1−ϑ; the solution of the above initial value system is given by

where

and zϑ=M(ϑ+Ω(1−ϑ)).

Proof. Utilizing Eq (3.1) and (3.2) we may verify that

by Eq (2.4), we obtain

Further calculation then yields the following:

By substituting Eq (3.14) into Eq. (3.16), the following results are obtained:

This finishes the proof.      �

Remark 1. Suppose that $\mathcal{G}_i\in{C}[0,T]$, with i=1,2,3,4,; then, we have

subject to the following conditions:

Therefore, $\widetilde{\mathcal{S}}=\mathcal{S}_0, \widetilde{\mathcal{E}}=\mathcal{E}_0, \widetilde{\mathcal{I}}=\mathcal{I}_0$ and $\widetilde{\mathcal{V}}=\mathcal{V}_0$; hence, the soultion given by Eq (3.12) presents the condition of continuity. On the basis of the above conditions the solution exists for the proposed model.

4.   Approximate solution via the Laplace Adomian decomposition method

Applying the Laplace transform to the mABC, as given in [37] to the considered model (1.3) gives

Proceeding we get

By application of the initial conditions we get

Decomposing each quantity as yields

We write the non-linear term as follows:

where

Plugging all of the above values into Eq (4.3) we get

after calculating the above Eq. (4.4) we get

Considering the first three terms and by applying the inverse Laplace transform we get

By putting the values $\mathcal{S}_1,\ \mathcal{E}_1,\ \mathcal{I}_1,\ \mathcal{V}_1$ into the second term of each quantity, we get

Lastly, each quantity can be written as follows:

4.1.   Stability and uniqueness of the solution

Theorem 1. Suppose that a Banach space is denoted by (B,|.|) with a mapping T:B→B which satisfies

for every x,y∈B; also, 0≤Θ and 0≤π<1. Then, T is Picard T-stable.

Theorem 2. Suppose T to be a self map defined as given below:

Thus, the iteration is T-stable in L1(x,y), if we obtained the following results:

Proof. We need to prove that T has a fixed point; for this we use (q,p)∈N×N

Taking the first equation of Eq (4.19) and calculating the norm on both sides, we have

with the help of the triangular inequality, and by solving Eq (4.20), we obtain

By establishing the relation as follows:

replacing the above relation in Eq (4.21), we get the following relation

Furthermore, the convergent sequence $\mathcal{V}_q$ is bounded.

Additionally, one may obtain different constants N₁, N₂, N₃ and N₄ for every t such that

Therefore, considering Eq (4.23) and Eq (4.24), we obtain

here, w₁ and w₂ are the functions of $\mathcal{L}^{-1}\bigg[\frac{(1-\vartheta)(s^\vartheta+\mu_{\vartheta})}{\mathbf{M}(\vartheta)s^\vartheta}\mathcal{L}\bigg]$. Considering the same procedure for the remaining equations, we have

Hence, T has a fixed point. Applying Eq (4.25) and (4.26), we assume that

Therefore, the conditions of Theorem 1 are satisfied; hence, the proof is complete.

Theorem 3. The aforementioned method gives a unique solution for the considered model.   �

5.   Numerical simulation

The graphical representations for the proposed model under the modified operator have been established by using the data taken from [9]. Four studied compartments of the model were tested on different fractional orders by using the obtained numerical scheme. The initial values and corresponding parameter descriptions are given in Table 1.

In Figure 1(a–d), the dynamics of all four compartments are shown for different fractional orders under the modified fractional operator which has extra terms as compared to the ABC operator. The susceptible population gains stability after increasing to some peak value and the stability is best achieved for relatively small fractional orders instead of the higher order. The exposed class declines and then increases to the equilibrium point to achieve stability. During this period, the exposed population is transferred to the infected population, whose density is increasing and then stabilizes, or decreases as a result of individuals recovering from the infection after the vaccination is applied to the population. The vaccination process is terminated once the infection is controlled.

In Figure 2(a–d), the dynamics of all four compartments are shown for different arbitrary orders under the modified fractional operator with extra informative terms. The susceptible class gains stability after declining to some small value and the stability is achieved faster with small fractional orders. The exposed class here declines quickly and then increases to the equilibrium point to achieve the convergence. During this period, the susceptible population is transferred to the exposed class, and then to the infected class whose density increases and then stabilizes, or decreases as a result of recovering from the infection after the vaccination is injected into the total population. The vaccination process is reduced as the transmission rβ is increased and the infection is controlled after more populations.

In Figure 3(a–d), the rate of transmission is decreased which affects the dynamics of all four compartments as shown for different arbitrary orders. The susceptible class increases and gains stability after reaching some high value. The exposed class here declines much more quickly and then converges to its equilibrium point. In this case, after a small increase, the infection is controlled quickly. Here the vaccination rate also increases, i.e. as the vaccination rate increases the infection may be easier to control.

In Figure 4(a–d), the rate of transmission is again decreased which affects the dynamics of all the four compartments, as shown for different arbitrary orders. In this case, the dynamics are about the same as observed for the previous dynamics but the dynamics can be controlled more quickly. Such dynamics provide the total density of each compartment in the form of a continuous spectrum lying between 0 and 1.

6.   Conclusions

This study developed the analysis of the impact of vaccination on the spread of polio in the human population by using a four-compartment system with the generalized novel mABC fractional operator. The dynamical behaviors of all compartments have been successfully examined for different fractional orders to demonstrate the validity of an extra degree of freedom in the selection of the derivative order. The fixed-point theory shows the uniqueness of the solution in the generalized novel modified format. The stability analysis of the model has been conducted by T-Picard-type stability techniques. The numerical solution for the model has been achieved by using the Laplace transform, along with the Adomian technique in the format of the decomposition process. Such a study involved the analysis of the complex geometry in the dynamical system. The study was performed by using various fractional orders and iterations with different transmission rate; also, each curve has been plotted for six different fractional orders and compared with the integer-order result. The study also shows that the impact of vaccination is more significant on the total population. In terms of stability, relatively smaller fractional orders yield high accuracy. In all of the numerical simulations, decreasing the transmission rate led to better control of the polio infection.

Use of AI tools declaration

The authors declare that they have not used artificial intelligence tools in the creation of this article.