Computer Modeling in Engineering & Sciences

DOI: 10.32604/cmes.2022.019440

ARTICLE

Stability Analysis of Predator-Prey System with Consuming Resource and Disease in Predator Species

Asad Ejaz1, Yasir Nawaz1, Muhammad Shoaib Arif1,3,*, Daoud S. Mashat2 and Kamaleldin Abodayeh3

1Department of Mathematics, Air University, PAF Complex E-9, Islamabad, 44000, Pakistan
2Department of Mathematics, King Abdulaziz University, Jeddah, 22254, Saudi Arabia
3Department of Mathematics and Sciences, College of Humanities and Sciences, Prince Sultan University, Riyadh, 11586, Saudi Arabia
*Corresponding Author: Muhammad Shoaib Arif. Email: marif@psu.edu.sa
Received: 24 September 2021; Accepted: 13 January 2022

1  Introduction

The pronounced work of Lotka et al. [1,2] in the field of predation prey relationship studies has made this area significant for research purposes. Great attention by scientific society has brought this field to the next level, as evident in population dynamics and mathematical ecology. The pioneering work of Kermack and McKendrick [3] has opened a new gospel in epidemiology. Following the work of [1–3], Maiti et al. [4] have studied the ratio-dependent predator-prey model and explored this model’s stochastic and deterministic features [4]. Peculiarly observed four variable balancing points and static factors of the system to understand how the delay affects the system’s movements. Moreover, Hopf bifurcation [5] is discussed by using bifurcation theory.

Wang et al. constructed a model whose stages were constructed on predation prey strategy to study a delay, stage-structured model. This study showed the delay when the prey population was at its immature stage and drew some essential ecological conclusions by analyzing the model in terms of global stability. Considering the study made by Zang et al. [6,7] and his coworkers examined the eco-epidemiological model and alleged that sickness is spreading among prey populations. This study took standard incidents and observed the behavior of the proposed system in all potential equilibrium to formulate conclusions. Additionally, Haque and coworkers examined Hopf bifurcation theory. Following Holling-Tanner functional response, Haque and venturino did a predation prey study; however, in [8], Yi et al. have investigated the concept of control and configurability in the context of predation prey paradigm of behavior. Some authors have discussed the chaotic and non-chaotic behavior of the dynamical system and presented optimal chaos control to overcome chaotic behavior [9]. Hyper chaotic behavior of the non-autonomous cardiac conduction system and overcome techniques are discussed in [10]. Some authors investigated the dynamical behavior of tri-atomic molecules using the fractional calculus approach [11].

Ecology and epidemiology are two distinct fields with fundamental effects on scientific work and society. However, there are numerous parallels between the two fields. Nowadays, scientists have merged both fields to form a new area known as eco-epidemiology that has become increasingly important, as indicated in the studies [12–15]. Proposes and analyses a predator-prey system where predators infect prey due to the polluted environment’s influence [14]. Sinha et al. [14] discovered that diseased prey is more venerable than predators living in a polluted environment. Many writers in the literature suggested eco-epidemiological models that focused solely on disease in the target population [16]. A Lotka-Volterra type predator-prey system was examined by Johri and colleagues in which an illness infected the prey. They investigated the local and global stability and its plausibility; however, Mukhopadhyaya et al. [17] investigated how diffusion and delay affected the system’s dynamics by considering a prey-predator model with Holling type II functional response. Following their efforts, Greenhalgh et al. [18] formulated the influence of diffusivity on the stability and persistence focused discussion of the model using a predator-prey model with disease in the prey population. To compare with Hu et al. [19] three-dimensional prey-predator delayed model, Chattopadhyay et al. [20] created highly organized two and three-dimensional models that articulated the local stability results and compared their findings to Hu and Li’s findings. A fractional chaotic system with order two and three nonlinearities is discussed in [21]; moreover, NSFD is suggested for the model. A fractional predator-prey model is investigated in [22] using some uncertain initial conditions, and the authenticity of the proposed model is proved via numerical simulation. In [23], the authors investigated the two-strain epidemic model. They have analyzed the model having two vaccinations. Liu et al. [24] analyzed semi-trivial, infection-free, and predator-free equilibrium situations using the impulsive predator-prey model; however, Naji et al. [25,26] examined eco-epidemiology’s local and global features of the delayed autonomous predator-prey model. Xiao et al. [27] has introduced the application of partial differential equation using the predator-prey model with infected prey population known as “prey infection system.” The researcher realizes conversing prey to predator is constant (i.e., independent of gestation period) if k=k0e−d(τ)¯ Where dτ¯ are mortality and delay, respectively. Xiao et al. [27] Elaborated that the delay could have a stable and DE stable effect on the system, narrating the positive and negative equilibrium. Linear models for series RL and RC are studied using truncated M-derivative operators [28]. Local asymptotic stability is discussed for transmission dynamics of SARS-Cov-2 [29]. Numerical simulations are provided using RK-2, RK-4, and Euler methods. Some authors used the technique of optimal control to investigate the stochastic and deterministic Hepatitis B Virus model [30].

Additional work on ecology and epidemiology can be seen in [31–33]. Musa et al. [29] examined the predator-prey system in which the predator population was infected with the pathogen. Their findings showed that infection among predators acts as a counterbalance to the predator of enrichment, and parasites removing themselves from food webs could have disastrous penalties. Biological control and resource management of these parasites were also part of their research. Haque [32] Proposed that the introduction of the pathogen in the predatory population acts as a biological control that maintains the predator’s population. Pal et al. [33] investigated a model in which predators were infected with the pathogen and determined their instability; however, several studies were based on the pathogenicity of both predator and prey, focusing on local or global equilibrium. Pattern formation like Turing spots and spiral-like structures in the predator-prey system are discussed [34]. Han et al. [35], using SIS or SIR, proposed four modifications to the prey-predator model, and the threshold of alterations caused by parasite infections was studied. They looked at models with both standard and mass action incidences, as well as a combination of the two. The global stability for the COVID-19 model is discussed in detail in [36], using actual data of Pakistan. Caputo type differential operator is employed to analyze the linear mass springer damper system [37]. Hsieh et al. [38] introduced a model that examined the pathogenicity of both predator and prey hence measuring a large number of threshold parameters for local analysis of distinct equilibria of the proposed system. They confer that the different combinations of variables may lead the infection to be eradicated from the environment due to the killing of infected prey. Hence, they demonstrated that illness in both populations could either promote or impede coexistence. The precise influence of infection on coexistence must be determined in each situation. After this work is completed, it will be able to propose more applications for the current methodology [38–44]. Some researchers have provided the numerical solution of the predator-prey model through numerical schemes like Adomian Decomposition, El Kalla, etc. The comparison between exact and numerical solutions is also provided in plots [45]. Qureshi et al. [46] studied the blood ethanol, concentration model. For fractionalizing the model, three different fractional operators are used.

From the available literature, it is found that some authors have discussed in detail the stability of the predator-prey system with common food resources [47]. The authors have considered the stability of the diffusive and kinetic system. We got the main motivation from the literature for the present work to the best of our knowledge. Nobody has discussed the problem of multi-species coexistence, keeping in view the common food resources and infection in one of the two species (predator). Therefore, we formulated a predator-prey model with a common food resource for both the predator and prey under the assumption that a fatal disease exists in the predator. We have proved by constructing the Lyapunov function that for this complex system, both the species predator (Susceptible and Infected) and prey can coexist in the ecosystem with and without diffusion.

Moreover, we have constructed a third-order numerical scheme for the model, and a faster scheme convergence than BTCS and Crank Nicolson is obtained. To the best of our knowledge, such a system considering disease and common food resources is not much discussed in the literature. We are hopeful that this particular research has an attraction for the eco-epidemiologists because they can better understand the complex behavior of species in the ecosystem by looking at the stability part in Sections 3 and 4. This particular study focuses on the coexistence of both species in the presence of disease in the predator for both diffusive and non-diffusive cases.

The present investigation provides a multi-species eco-epidemiological model discussing the existence of all species, which is the key result in such models. The findings in the present study show the survival of all the species in the same habitat. Some phase portraits are presented to show the behavior of species. Additionally, a faster convergence scheme is developed, which gives better results as compared to some existing schemes. The proposed scheme can help study such models.

The current study focuses on the formulation and analysis of the epidemiological predator-prey model in the sequence described here. Section 2 deals with the formulation of the multi-species model under the assumption that predator species have some fatal disease that can spread from predator to predator only. Moreover, the two classes of predators compete for the common food resources and the predation option. Sections 3 and 4 deal with the stability of the model for both non-diffusive and diffusive cases, respectively. For the sake of simulations, we constructed a third-order numerical scheme for the solution of the model in Section 5. The stability of the proposed scheme is also proved in the same section using the Von Neumann stability technique. In Section 6, numerical simulations are performed, and a comparison of the constructed scheme with Backward in Time and Central in Space and Crank Nicolson is displayed in the form of plots. Section 7 presents the discussion on the plots and findings of the present study, whereas the conclusion is provided in Section 8. Finally, some future directions are suggested in the last part of this research.

2  Model Formulation

The ecosystem is a complex natural structure where interacting species affect the whole environment and the existing species. Many species depend on the same resources, which gives rise to competition among them. Moreover, the interacting species are at the danger of infectious diseases in the real environment. Many researchers presented mathematical models to discuss the effects of infection on population density. The infections are so critical that they disturb the whole ecosystem. This problem motivated several investigators to put their model to study the population dynamics effectively. Here we propose a model to study the predator-prey system. The assumption of the models are as under:

H1: The predator species have two classes; susceptible and infected class.

H2: The disease can only spread from predator to predator, and the disease is mortal.

H3: The transmission is non-vertical.

H4: The predator classes depend on the same food resource f with different consumption rates a and b, respectively; additionally, both classes have predation activity as well.

In the light of the above assumptions, we take US(t),UI(t) as the predators (susceptible and infected) population density and V(t) as prey population density at any time t. Here we assume that U=US+UI is the total population density of predator species.

The above assumptions lead to the following model:

dUS(t)dt=US(afaUS+bUI+α1V−βUI−m1),t>0(1)

dUI(t)dt=UI(bfaUS+bUI+α2V+βUS−(µ+m1)),t>0(2)

dV(t)dt=V(cg−α1US−α2UI−n1),t>0(3)

where

Φ0={US(0)=US,0≥0,UI(0)=UI,0≥0,V(0)=V0≥0}(4)

is the set of initial conditions. Here m1 and n1 are the natural mortality rates of the predator and prey population, respectively, whereas µ is the death rate of an infected predator due to disease, f stands for the same food resources for both predator classes; and a,b are the consumption rates for susceptible and infected predators. The parameters α1 and α2 respectively represent the predations rates for both predator classes, whereas β stands for the transmission rate of infection from infected predators to susceptible predators. The parameter g displays the food resources for the prey population, and c is the consumption rate. Moreover, due to infection, the mortality rate of infected predators is high as compared to the susceptible class.

We can rewrite the above system of equations as

dUSdt=US(afaUS+bUI+α1V−βUI−m1),t>0(5)

dUIdt=UI(bfaUS+bUI+α2V+βUS−m2),t>0(6)

dVdt=V(cg−α1US−α2UI−n1),t>0(7)

Here the parameter m2 stands for the death rate of infected predators. It is obvious from the above discussion that m2 > m1.

3  Stability for the Predator-Prey System of Equations

In the present section, we compute the equilibria of the system. The equilibrium points are simply the solutions of the following system of equations:

0=US(afaUS+bUI+α1V−βUI−m1)(8)

0=UI(bfaUS+bUI+α2V+βUS−m2)(9)

0=V(cg−α1US−α2UI−n1)(10)

One can solve the above system of equations to find the equilibrium points which are

E1=(0,fm2,0),E2=(fm1,0,0),E3=(US,3,UI,3,V3)

E4=(US,4,UI,4,V4),E5=(US,5,UI,5,V5),E6=(US∗,UI∗,V∗)

where

US,3=0,UI,3=cg−n1α2,V3=cgm2−fα2−n1m2α2(cg−n1)

US,4=am22−bfβ−bm1m2β(am2−bm1),UI,4=afβ−am1m2+bm2β(am2−bm1),V4=0

US,5=cg−n1α1,UI,5=0,V5=cgm1−fα1−n1m1α1(cg−n)

US∗=((α1f+m1(cg−n1)α2−(cg−n1)(α1m1−β(cg−n1)))b−afα22(aα2−bα1)(−α2m1+m2α1−β(cg−n1))(11)

UI∗=((α2f+m2(cg−n1)α1−(cg−n1)(α2m1+β(cg−n1)))a−bfα12(aα2−bα1)(−α2m1+m2α1−β(cg−n1))(12)

V∗=((b(cg−n1)m1−a(cg−n1)m2+f(aα2−bα1)β+(−α2m1+m2α1)(−bm1+m2a)(aα2−bα1)(−α2m1+m2α1−β(cg−n1))(13)

In the present section, we study the stability of the fixed point, which deals with the coexistence of all the species. The existence condition for equilibrium point is

ϕ∗={(US∗,UI∗,V∗),US∗>0,UI∗>0,V∗>0}(14)

Theorem 1: Let a,b and f be positive and (11) holds then (US∗,UI∗,V∗) of (5)–(7) is stable in

ω={(US,UI,V),US>0,UI>0,V>0}

Proof: We construct the Lyapunov function to determine the stability of the coexistence equilibrium point (US∗,UI∗,V∗). Consider the following:

L(US,UI,V)=∫US∗USη−US∗ηdη+∫UI∗UIκ−UI∗κdκ+∫V∗Vξ−V∗ξdξ(15)

By taking derivative with respect to “t”, we have

dLdt=US−US∗USdUSdt+UI−UI∗UIdUIdt+V−V∗VdVdt

dLdt=(US−US∗)(afaUS+bUI+α1V−βUI−m1)+(UI−UI∗)(bfaUS+bUI+α2V+βUS−m2)+(V−V∗)(cg−α1US−α2UI−n1)

After some simplification, we get

dLdt=(US−US∗)afaUS∗+bUI∗−aUS−bUI(aUS+bUI)(aUS∗+bUI∗)+(UI−UI∗)bfaUS∗+bUI∗−aUS−bUI(aUS+bUI)(aUS∗+bUI∗)

dLdt=−f(aUS+bUI−(aUS∗+bUI∗))2(aUS+bUI)(aUS∗+bUI∗)≤0(16)

The above result proves the Lyapunov stability for the system of Eqs. (5)–(7).

4  Stability for Diffusive System

The dynamical population model, random walk influences the structure of the habitat. This movement can be varied from one place to the other depending on the factors like food resources and protection in the environment. In addition, species can vary with respect to density in the area they live. The following section presents the effect of diffusion on the coexistence of species.

4.1 Global Stability for the Self-Diffusive System

The self-diffusion system for (5)–(7) can be written as

dUS(x,t)dt=d11△US+US(afaUS+bUI+α1V−βUI−m1),xϵΩ,t>0(17)

dUI(x,t)dt=d22△UI+UI(bfaUS+bUI+α2V+βUS−m2),xϵΩ,t>0,(18)

dV(x,t)dt=d33△V+V(cg−α1US−α2UI−n1),xϵΩ,t>0(19)

where ∂US∂ϑ=∂UI∂ϑ=∂V∂ϑ=0 , xϵ∂Ω,t>0

Here US(x,t), UI(x,t) and V(x,t) respectively represent the population density of susceptible predator, infected predator, and prey at location x ϵ Ω and time t. Here Ω represents the bounded domain having smooth boundary ∂Ω, and ϑ is the outward unit normal on the smooth boundary. The Neumann boundary conditions show the zero population flux across the boundary.

The initial condition is as under

Φ0={US(x,0)=US,0(x)≥0,UI(x,0)=UI,0(x)≥0,V(x,0)=V0(x)≥0}(20)

Here d11,d22,d33 are non-negative self-diffusion parameters. According to Hollis results, the solution to the above system has global existence and boundedness [16]. The self-diffusive system has an equilibrium point E6=(US∗,UI∗,V∗) under the assumption (11).

Theorem 2 Let Ω be bounded domain with smooth boundary ∂Ω. Suppose a,b,d11,d22,d33 and f be positive and (11) holds then (US∗,UI∗,V∗) of (14)–(16) is stable in

ω={(US,UI,V),US>0,UI>0,V>0}

Proof. It is obvious that {(US,UI,V),US=0}, {(US,UI,V),UI=0} and {(US,UI,V),V=0} are the invariant manifolds. We establish the Lyapunov function to determine the global stability for the self-diffusive system as under

L(U,VS,VI)=∫Ω∫US∗USη−US∗ηdηdX+∫Ω∫UI∗UIκ−UI∗κdκdX+∫Ω∫V∗Vξ−V∗ξdξdX(21)

Now taking derivative with respect to “t” on both sides, we get

dLdt=∫ΩUS−US∗US∂USdtdX+∫ΩUI−UI∗UI∂UIdtdX+∫ΩV−V∗V∂VdtdX

=∫Ω(US−US∗US)(d11△US+US(afaUS+bUI+α1V−βUI−m1))

+∫Ω(UI−UI∗UI)(d22△UI+UI(bfaUS+bUI+α2V+βUS−m2))

+∫Ω(V−V∗V)(d33△V+V(cg−α1US−α2UI−n1),

=M1+M2+M3+M4(22)

where

M1=−d11US∗∫Ω1US2|∇US2|dX≤0(23)

M2=−d22UI∗∫Ω1UI2|∇UI2|dX≤0(24)

M3=−d33V∗∫Ω1V2|∇V2|dX≤0(25)

M4=∫Ω(US−US∗US)[US(afaUS+bUI+α1V−βUI−m1)]dX+∫Ω(UI−UI∗UI)[UI(bfaUS+bUI+α2V+βUS−m2)]dX+∫Ω(V−V∗V)[V(cg−α1US−α2UI−n1)]dX≤0(26)

M4=−∫Ωf(aUS+bUI−(aUS∗+bUI∗))2(aUS+bUI)(aUS∗+bUI∗)≤0(27)

From

dLdt≤0(28)

Hence (US∗,UI∗,V∗) is globally stable.

5  Construction of Numerical Scheme

An implicit scheme on three-time levels is constructed. The construction of the scheme is started with the difference equation given as

uin+1=auin+uin−14+Δt{b(∂u∂t)in+1+c(∂u∂t)in+e(∂u∂t)in−1}(29)

The Taylor series expansions for uin+1,uin−1,(∂u∂t)in+1 and (∂u∂t)in−1 are given as

uin+1=uin+Δt(∂u∂t)in+(Δt)22(∂2u∂t2)in+(Δt)36(∂3u∂t3)in+O((Δt)4)(30)

uin−1=uin−Δt(∂u∂t)in+(Δt)22(∂2u∂t2)in−(Δt)36(∂3u∂t3)in+O((Δt)4)(31)

(∂u∂t)in+1=(∂u∂t)in+Δt(∂2u∂t2)in+(Δt)22(∂3u∂t3)in+O((Δt)3)(32)

(∂u∂t)in−1=(∂u∂t)in−Δt(∂2u∂t2)in+(Δt)22(∂3u∂t3)in+O((Δt)3)(33)

Substituting Taylor series expansions (30)–(33) into Eq. (29), we obtain the following:

uin+Δt(∂u∂t)in+(Δt)22(∂2u∂t2)in+(Δt)36(∂3u∂t3)in+O((Δt)4)=auIn+14(uin−Δt(∂u∂t)in+(Δt)22(∂2u∂t2)in−(Δt)36(∂3u∂t3)in+O((Δt)4)+Δt{b(∂u∂t)in+bΔt(∂2u∂t2)in+b(Δt)22(∂3u∂t3)in+c(∂u∂t)in+e((∂u∂t)in−Δt(∂2u∂t2)in+(Δt)22(∂3u∂t3)in)}(34)

Comparison of coefficients of uin,Δt(∂u∂t)in,(Δt)2(∂2u∂t2)in and (Δt)3(∂3u∂t3)in on both sides of Eq. (34) leads to the following:

1=a+14(35)

1=−14+b+c+e(36)

12=18+b−e(37)

16=−124+b2+e2(38)

Solving Eqs. (35)–(38) gives the values of an unknown parameters a,b,c and e as

a=34,b=1948,c=56,e=148(39)

5.1 Stability Analysis

The stability of the numerical scheme for the diffusive system is checked by applying the Von Neumann stability criterion. For this criterion, the difference equation obtained by applying a numerical scheme is transformed to another equation containing amplitudes and phase angle of the waves. After that, stability condition is imposed on the obtained equation, which is based on the quotient of amplitudes of the wave at ``(n+1)th'' and nth levels. For starting the stability procedure, consider the following diffusive linear predator-prey equation having the form:

∂u∂t=d1∂2u∂2x+α¯u(40)

The discretization of Eq. (40) using the presented scheme is given as

uin+1=auin+14uin−1+Δt{bd1(ui+1n+1−2uin+1+ui−1n+1(Δx)2)+bα¯uin+1+cd1(ui+1n−2uin+ui−1n(Δx)2+α¯uin)+cα¯uin+ed1(ui+1n−1−2uin−1+ui−1n−1(Δx)2+eα¯uin−1)}(41)

According to Von Neumann stability criterion, consider the following transformations:

uin+1=En+1eiIψ,ui±1n=Ene(i±1)Iψ,ui±1n−1=En+1e(i±1)Iψ where I=−1(42)

Substituting transformations (42) into Eq. (41), the following is obtained:

En+1eiIψ=aEneiIψ+14En−1eiIψ+Δt{bd1(e(i+1)Iψ−2eiIψ+e(i−1)Iψ(Δx)2)En+1+α¯En+1eiIψ+cd1(e(i+1)Iψ−2eiIψ+e(i−1)Iψ(Δx)2)En+cα¯EneiIψ+ed1(e(i+1)Iψ−2eiIψ+e(i−1)Iψ(Δx)2)En−1+eα¯En−1eiIψ}(43)

Dividing both sides of Eq. (43) by eiIψ yields

En+1=aEn+14En−1+Δt{bd1(2cosψ−2(Δx)2)En+1+α¯En+1+cd1(2cosψ−2(Δx)2)En+cα¯En+ed1(2cosψ−2(Δx)2)En−1+eα¯En−1}(44)

Collecting coefficients of En+1 on the left-hand side of Eq. (44) provides

(1−2dd1−Δtα¯)En+1=aEn+14En−1+2cdd1(cosψ−1)En+Δtα¯En+1+2edd1(cosψ−1)En−1+eΔtα¯En−1(45)

where d=Δt(Δx)2.

Eq. (45) can be expressed as

En+1=AEn+BEn−1(46)

where

A=a+2cdd1(cosψ−1)+cΔtα¯1−2dd1−Δtα¯ and B=b+2edd1(cosψ−1)+eΔtα1−2dd1−Δtα¯(47)

Since the presented numerical scheme is constructed on three-time levels, one more equation is required to get an amplification matrix. An additional equation is expressed as

En=1.En+0.En−1(48)

System of Eqs. (46) & (48) can be expressed in the following matrix-vector equation:

[En+1En]=[AB10][EnEn−1](49)

The amplification factor is a matrix, and the stability conditions can be imposed on the Eigenvalue of the amplification matrix and these conditions are expressed as

|A+A2+4B2|≤1and|A−A2+4B2|≤1(50)

The scheme will be stable if it satisfies inequalities (50) for positive eigenvalues of the amplification matrix, and if the eigenvalues of an amplification matrix are negative, then the stability condition is expressed as

|B|2≤1(51)

where A & B are expressed in (47).

6  Numerical Simulation

The presented third-order multistep finite difference scheme can efficiently be applied for solving the eco-epidemiological population model. Here, this scheme will be applied for the numerical simulation of the model presented in (17)–(19). The model presentation under the scheme applied is given as

US,in+1=aUS,in+14US,in−1+Δt[b{d11(US,i+1n+1−2US,in+1+US,i−1n+1)/(Δx)2+US,in+1(afaUS,in+1+bUI,in+1+α1Vin+1−βUI,in+1−m1)}+c{d11(US,i+1n−2US,in+US,i−1n)/(Δx)2+US,in(afaUS,in+bUI,in+α1Vin−βUI,in−m1)}+e{d11(US,i+1n−1−2US,in−1+US,i−1n−1)/(Δx)2+US,in−1(afaUS,in−1+bUI,in−1+α1Vin−1−βUI,in−1−m1)}](52)

UI,in+1=aUI,in+14UI,in−1+Δt[b{d22(UI,i+1n+1−2UI,in+1+UI,i−1n+1)/(Δx)2+UI,in+1(bfaUS,in+1+bUI,in+1+α2Vin+1+βUS,in+1−m2)}+c{d22(UI,i+1n−2UI,in+UI,i−1n)/(Δx)2+UI,in(bfaUS,in+bUI,in+α2Vin+βUS,in−m2)}+e{d22(UI,i+1n−1−2UI,in−1+UI,i−1n−1)/(Δx)2+UI,in−1(bfaUS,in−1+bUI,in−1+α2Vin−1+βUS,in−1−m2)}](53)

Vin+1=aVin+14Vin−1+Δt[b{d33(Vi+1n+1−2Vin+1+Vi−1n+1)/(Δx)2+Vin+1(cg−α1US,in+1−α2UI,in+1−n1)}+c{d33(Vi+1n−2Vin+Vi−1n)/(Δx)2+Vin(cg−α1US,in−α2UI,in−n1)}+e{d33(Vi+1n−1−2Vin−1+Vi−1n−1)/(Δx)2+Vin−1(cg−α1US,in−1−α2UI,in−1−n1)}](54)

The above Eqs. (52)–(54) give the discretization of the diffusive system.

7  Discussion

We presented an eco-epidemiological model for the predator-prey system in the present study. Two classes of predator species (susceptible and infected) with disease transmission from predator to predator are considered. The transmission of infection is supposed to be non-vertical in the present research. We have proved the coexistence of predator and prey in the same habitat. Figs. 1 & 2 show the convergence of three numerical schemes for the modified model of the eco-epidemiological diffusive mathematical model using different parameter values. It is clear from the plot that faster convergence is obtained using the proposed numerical scheme. The maximum number of norms for each Eqs. (52)–(54) is determined over the number of iterations used by each scheme, and Figs. 1 & 2 illustrate the convergence speed. Since the proposed scheme is third-order accurate, it converges faster than existing first and second-order schemes in time. It is to be noted that the second-order central space discretization is used for three schemes. Figs. 3 & 4 demonstrate the comparison of solutions obtained by utilizing three numerical schemes over time considering two different data sets, while Figs. 5 & 6 show the comparison of schemes over the spatial coordinate x. The phase portraits for the ODEs model are given in Figs. 7–9 that express the relationships between Us,UI & V. The initial conditions for drawing phase portraits are obtained in the interval (0,1). In Figs. 1 & 2, Nx,Nt, tf & L respectively show the number of grid points, the number of time levels, final time, and length of the domain.

Figure 1: Convergence of three schemes using a=0.1,b=0.2,f=0.1,α1=0.3,β=0.2,m1=0.1,α2=0.2,m2=0.1,c=0.1,g=0.1,n1=0.1,Nx=50,Nt=300,tf=10,L=7,d1=d2=d3=d4=0.4

Figure 2: Convergence of three schemes

Figure 3: Comparison of three schemes over time using a=0.1,b=0.01,f=0.2,α1=0.1,β=0.2,m1=0.3,α2=0.5,m2=0.3,c=0.3,g=0.1,n1=0.3,Nx=50,Nt=500,tf=15,L=7,d1=d2=d3=d4=0.4