NUMERICAL ANALYSIS OF PREDATOR-PREY MODEL IN PRESENCE OF TOXICANT BY A NOVEL APPROACH

It is assumed in this paper that the environmental toxicant affects both prey and predator population. We consider the Holling Type-I predator-prey model and stage-structured predator-prey model in presence of environmental toxicant in the form of Volterra integro-differential equations. We solve the system of equations by a novel approach Trapezoidal Base Homotopy Perturbation Method. Using the dynamical behaviour of the systems as a validating tool, we validate the above mentioned method as an efficient method to solve Volterra integro-differential predator-prey model.


INTRODUCTION
The modeling of predator-prey populations is an attractive research area for scientists from different fields like biologists, ecologists, mathematicians, economists, etc. Because there is the scope of proposing a wide variety of models to analyze different mechanisms such as control of populations, oscillations of population, delay maturation time effects and pattern formation in predator-prey models with diffusion, migration, refuge, harvesting, etc. One of the important attributes in these models is the predator-prey coexistence, which is a key feature to obtain realistic models of predator-prey interactions. Many times eliminating one population from the environment and analyzing the other makes the models easy to handle. But coexistence makes the population model complicated to analyze. Mathematically, the first type of coexistence is established by a stable equilibrium which are the most basic solutions to the models. The second type of coexistence is given by stable periodic solutions and this is the most commonly accepted type in predator-prey systems [16]. Since Lotka (1925) and Volterra (1927) proposed the Lotka-Volterra population model, many mathematical models have been formed based on realistic explicit and implicit biological assumptions [1,2,3,4,5,6].
One of the essential effects on population dynamics which is the need of the hour today to be studied in detail is the effect of environmental pollution. Many mathematical biologists have worked in this area. The problem of estimating qualitatively the effect of a toxicant on a population by mathematical models began in the early 1980s. For a general class of single population models with toxicant stress, Ma et al [7] obtained a survival threshold under the hypothesis that the capacity of the environment is large relative to the population biomass and that the exogenous input of toxicant into the environment is bounded. In 1987, Ma and Hallam [8] studied the two-dimensional non-autonomous Lotka-Volterra model by the average method and obtained sufficient conditions for persistence and for the extinction of the population. The threshold of the survival for a system of two species in a polluted environment was studied by Huaping and Ma [9]. Population toxicant coupling has been applied in several contexts including Lotka-Volterra and chemostat-like environments, resulting in ordinary, Integro-differential, and stochastic models [10,12,13,15]. The Volterra model for a species population growth in a closed structure is as follows: where a > 0 is the birth rate coefficient, b > 0 is the crowding coefficient, c > 0 is the toxicity coefficient, p 0 is the initial population and p = p(t) denotes the population at time t.
Considerable research work has been invested recently by Small [14], TeBeest [17] among others to the development of efficient strategies to determine numerical and analytic solutions of the population growth model 1. Although a closed-form solution has been achieved by them, but it was formally shown that the closed-form solution cannot lead to an insight into the behavior of the population evolution. Analytical approximations and Pade approximants for problem 1 are considered by Wazwaz [18]. In [19], a comparison of the Adomian decomposition method and Sinc-Galerkin method is given and it is shown that the Adomian decomposition method is more efficient for the solution of Volterra's population model 1. In [20] [23] is efficiently used to solve fractional differential equations. This is an attempt to apply a modified method in analogy to Deithlm et.al [23], named as Trapezoidal Based Homotopy Perturbation Method(TBHPM) to solve integro-differential equations. Instead of comparing the results obtained by TBHPM with other established methods, in this paper, we have considered the stability of the system as a tool to validate the efficiency of the method. Also, we expect that the presence of toxicant in the environment (i) affects both the population. (ii) But as the predator population is dependent on prey population for food, therefore due to the toxicant effect as well as lack of prey, the predator population reduces faster than the prey population, (iii) in the case of the SSPP model toxicant effects will be more on the immatured prey than the matured prey. We shall observe that the solutions obtained by TBHPM establish these expectations.

INTEGRO-DIFFERENTIAL EQUATION
Trapezoidal Based Adams-Bashforth Method applied by Diethelm et. al in [23]. In this paper, we have modified this method using the Homotopy perturbation method and call it as Trapezoidal Based Homotopy Perturbation Method(TBHPM) and we have applied this method for a class of integro-differential equation. Consider the nonlinear integrodifferential equation of the following form The initial value problem 1 is equivalent to a Volterra integral equation In the sense that a continuous function solves Eq. 3 if and only if it solves Eq.2.
Firstly the product trapezoidal quadrature formula is applied to replace the integrals of Eq. 3, setting h = T N ,t n = nh, n = 0, 1, 2.....N ∈ Z + . The integral part in the Eq. 3 is again simplified by using Trapezoidal rule. Then Eq. 3 can be discretized as follows: Eq. 3 may be rewritten as The right hand side of system 6 contains term x(t n+1 ), so it is an implicit scheme. We shall evaluate the term x(t n+1 ) in right hand side by HPM and the scheme is named as Trapezoidal rule based Homotopy Perturbation method(TBHPM).

Holling Type-I Predator-Prey Model with Environmental Pollution:
In this section, we investigate the perdator-prey model of the form where u = u(t) and v = v(t) are the population of identical individuals at a time 't' which exhibits crowding and sensitivity to the amount of toxic produced [14]. The meaning of the symbols used in system 7 are: r is the intrinsic growth rate and k is the environmental carrying capacity of prey. a and b represent respectively the predation rate of potential prey and conversion factor of predator after consuming potential prey. d is the death rate of predator and e is the toxic coefficient. If e = 0 we have the well-know logistic equation discussed in details in [15]. The integrals contained in the last terms of model 7 indicates the total amount of toxins produced since time zero. As stability is the key tool to validate the numerical strategy here, the bellow stated remark gives a short note on it.
Remark 3.1.1. The Jacobian matrix is corresponding to Eq. 7 is In absence of toxic effect, the system 7 has three equilibrium states given by abk E 0 and E 1 are always exist,but E 2 will exit only when bk > d and the model 7 posses a coexistence equilibrium at the point E 2 .
(3) Eigen values at E 2 are For stability of this point we must have bk > d. Again the point is a stable node if −4b 2 k 2 + 4bdk + dr > 0 and stable spiral if −4b 2 k 2 + 4bdk + dr < 0. In the following section, we have numerically analysed this behaviour by using TBHPM and shown that the method gives the expected results.
3.1.1. Solution of the Integro-differential Equation 7 by TBHPM. In this section we shall apply the TBHPM described in section 2 to investigate the proposed integro-differential equation 7 numerically and establish the efficiency of the method in solving integro differential equation. Applying TBHPM method to model 7 yields Where, s j,n+1 is defined by Eq. 5.
We represent the terms u(t n+1 ) and v(t n+1 ) in the right hand side by u p (t n+1 ) and v p (t n+1 ), which shall be predicted by using Homotopy Perturbation Method. We may constructed the homotopy as where According to HPM substitute in the above equations We obtain the following equations extinction of predator population is faster than the prey populations. This is because predator are affected by toxicity as well as lack of food source(as prey population is affected due to toxicity). This analysis done by TBHPM for the predator-prey dynamics in presence of toxicant establishes that proposed method is reliable method to solve integro-differential equations.

Stage-Structured Predatpor-Prey(SSPP) Model with Environmental Pollution.
Every species has two or more stages in its entire lifespan. Therefore mathematical biologists have proposed SSPP model incorporating various factors [24]- [34]. With reference to the model worked out in [34], in this section we consider a stage structured predator-prey model incorporating toxicant effect: All the parameters r, δ 1 , d 1 , β , γ 1 , m, δ 2 , d 2 , γ 2 , e 1 , e 2 , δ 3 , d 3 are all positive and there biological meanings are as follows: r is the intrinsic growth rate of immature prey depending on parents, δ 1 , δ 2 and δ 3 represent Intraspecific competition rate between immature prey, matured prey population and predator population respectively. d 1 , d 2 and d 3 are the natural death rate of immature prey, matured prey and predator population respectively. β is grownup rate of immature prey population. γ 1 and γ 2 represent the attack rate of predator to immature and matured prey population respectively. c 1 and c 2 are respectively the conversion factors of predator after consuming immature and matured prey with 0 < c 1 < 1, 0 < c 2 < 1. e represents the toxicant effect. Stability criteria of the system 8 in absence of toxicant effect and in presence prey refuge are discussed in details in [34].  (1) Trivial E 0 = (0, 0, 0) always exists.
Here u * = δ 2 (v * ) 2 +d 2 v * β and v * is the positive root od the equation and v * is a positive root of the equation where, Interior equilibrium point exists if A 4 < 0 with (A 2 > 0 or A 3 < 0) and d 3 < Moreover the predator extinct equilibrium is stable if Based on this information we shall verify the efficiency of the TBHPM in solving the SSPP model.

Solution of SSPP model by TBHPM.
Applying TRHP method to model 8 yields We represent the terms u(t n+1 ),v(t n+1 ) and w(t n+1 ) in the right hand side by u p (t n+1 ),v p (t n+1 ) and w p (t n+1 ) and find them by applying the HPM. To establish the efficiency of the method we consider two types of equilibrium: (i) predator extinct equilibrium and (ii) co-exixtence or interior equilibrium.    From Table 1, we noticed that when e = 0, the population of immature prey, mature prey, and predators continues to reach a point of equilibrium at (3.70681, 2.9111, 0). Table 2 and Table 3 shows that as the toxicity level increases, the population of immature prey decreases rapidly than mature prey. Furthermore, we conclude that in    Table 4 displays the results when e = 0, the population of immature prey, mature prey, and predators begins to achieve an interior equilibrium point at (5.79227, 4.08069, 2.05682), as we saw earlier. Table 5 and Table 6 demonstrate that as the degree of toxicity rises, Immature prey populations decline rather than mature prey populations. Furthermore, we believe that predator species decline more rapidly than other populations in the face of toxicity and diminishing prey populations.

CONCLUSION
In this work, we have investigated two predator-prey models under the influence of environmental pollutants. Taking convergence towards the stable point of equilibrium point as a validity factor, we have shown the efficiency of TBHPM in solving the population model. We have solved Holling Type-I predator-prey model with toxic term and analyzed results graphically.
Whereas in section 3.2 for the SSPP model results are presented in tabular form. Graphical analysis of the solution of the system 7 gives expected results. From the figures, it is seen that even though both predator and prey population are exposed to the same amount of toxicant, predator population decay faster than prey population. This indicates that the predator population faces death due to natural death, toxic effect, and lack of food. In the case of model 8, in absence of toxicity, we observe that the solution converges to predator extinct equilibrium and interior equilibria. Whereas, in presence of toxicity in the environment, immatured prey are affected by its effect more than the matured prey. Also, all the table show it clearly that in presence of toxicity predator population decays faster than other population. Moreover, the toxicity effect is directly proportional to the growth of the predator and prey population. These observations strongly support the proposed TBHPM.

CONFLICT OF INTERESTS
The author(s) declare that there is no conflict of interests.