Abstract

We have discussed the dynamical behaviour of a single-species population model in a polluted environment which describes the effect of toxicants on ecological system. Boundedness, positivity, and stability analysis of the model at various equilibrium points is discussed thoroughly. We have also studied the effect of single discrete delay as well as double discrete delays on the population model. Existence conditions of the Hopf bifurcation for single time delay are investigated. The length of delay preserving the stability is also estimated. The direction and the stability criteria of the bifurcating periodic solutions are determined by using the normal form theory and the center manifold theorem. The stability of the model with double time delays is investigated by using the Nyquist criteria. By choosing one of the delays as a bifurcation parameter, the model is found to undergo a Hopf bifurcation. Some numerical simulations for justifying the theoretical results are also illustrated by using MATLAB, which shows the reliability of our model from the practical point of view.

1. Introduction

In the world today, the pollution of the environment is a threatening problem due to the rapid development of industrialization. The presence of toxicants in the environment decreases the growth rate of species and its carrying capacity. In recent years, the most serious problem to society is the change in the environment caused by pollution, affecting the long-term survival of species, human life style, and biodiversity of the habitat. A large quantity of the toxicants and contaminants enter into the ecosystem one after another which seriously threaten the survival of the exposed population including humans. Therefore, the study of the effects of toxicants on the population and the assessment of the risk to populations are becoming quite important.

The problem of estimating qualitatively the effect of a toxicant on a population by mathematical models is a very effective way. The study of deterministic dynamic population models with toxicant effect was proposed by Hallam and his colleagues in 1980s [13]. This model was revisited by many researchers [4, 5]. He and Ma in [6, 7] discussed the survival of a single-species in a polluted environment, considering the organism’s uptake of toxicant from the environment and egestion of toxicant into the environment. In [4] Ma studied a Leslie resource-consumer model and obtained the threshold between the persistence and extinction of the consumer. In [8], Buonomo et al. studied the effect of variation of the population on the toxicant concentration in the organism and environment. Samanta and Maiti [9] studied a dynamical model of a single-species system in a polluted environment under two cases: constant exogenous input of toxicant and rapidly fluctuating random exogenous input of toxicant into the environment by means of ordinary and stochastic differential equations. There are many other works on the effects of a single toxicant on various ecosystems by using mathematical models [1015].

Time delays have been incorporated in ecological models by many authors [1627]. Introducing such delays in ecological models makes them more realistic and reasonable. Delay may have very complicated impact on the dynamics of a system since delay can cause the change of stability and can bifurcate periodic solutions.

In this paper, we have developed a single-species population model in polluted environment. In Section 2, we present a brief sketch of the construction of the model. Boundedness and positivity analysis of the solutions is shown which implies that the system is ecologically well behaved. In the next section, we have discussed the existence and stability analysis of various equilibrium points under zero-exogenous input and nonzero constant exogenous input. Next, we obtain necessary conditions for the existence of interior equilibrium and local and global stability of the system at . The analysis of single delayed model is described in Section 4. This analysis shows that there is a critical value of delay below which the system is stable and above which the system becomes unstable at the interior equilibrium . The system undergoes a Hopf bifurcation around at . In Section 5, we have estimated the length of delay to preserve the stability around . The direction and the stability criteria of the bifurcation periodic solutions by using the normal form theory and the center manifold theorem are discussed in Section 6. In the next section, the analysis of double delayed model is discussed. The stability of the model with double delays is investigated by the Nyquist criteria. By choosing one of the delays as a bifurcation parameter, the model is found to undergo a Hopf bifurcation. Some of the important analytic results are numerically verified by using MATLAB in Section 8. Finally, Section 9 contains the general discussions of the paper and ecological implications of our mathematical findings.

2. Basic Mathematical Model

In this paper, we analyze a model which describes the effect of toxicant on a single-species population. Here we assume the following:: concentration of the population biomass,: concentration of the toxicant in the population,: concentration of the toxicant in the environment.

The model satisfies the following assumptions. (A1)There is a given toxicant in the environment and the living organisms absorb into their bodies part of this toxicant so that the dynamics of the population is affected by the toxicant. (A2)For the growth rate of population, we assume that the birth rate is and the death rate is , where , , , and are assumed to be positive constants.

We consider the model with initial data , , and . Here,: depletion rate of toxicant in the environment due to its intake made by the population,: depletion rate of toxicant in the population due to egestion,: depletion rate of toxicant in the population due to metabolization process,: depletion rate of toxicant in the environment,: exogenous toxicant input rate which is assumed to be a smooth bounded nonnegative function of .

We can see that if , then will be going to extinct. So we suppose that The model we have just specified has nine parameters, which make the analysis difficult. To reduce the number of parameters, we make the following transformations: , , and . So, system (1) becomes where , , , , and .

Therefore,  from (2), we  have

Theorem 1. Each component of the solution of system (3) subject to initial conditions (4) is positive and bounded.

Proof. From the first equation of the system (3) we get
From the second equation of (3), we get
Similarly, from the third equation of (3), we get
Therefore, , , and , for all .
From the first equation of (3), we have
By a standard comparison theorem and by (5), we have,
From the third equation of (3), we have
So, we conclude that each component of the solution of system (3) subject to (4) is positive and bounded for all . This completes the proof.

3. Stability Behaviour of the Model

3.1. Case I: Zero Exogenous Input

Theorem 2. If , then the model (3) has nonnegative equilibria which is unstable and which is locally asymptotically stable. The interior equilibrium is not feasible.

Proof. The variational matrix of system (3) at is
The eigen values are , , and . So, it is obvious that, is unstable (hyperbolic saddle).
Now the varitiational matrix of system (3) at is
The characteristic equation of is , where and , since .
The eigen values are and .
Since , and , therefore the signs of the real parts of , are negative. Hence is locally asymptotically stable.
It is noted here that the other equilibrium point is not feasible, since are opposite in sign since , by (2).
This completes the proof.

3.2. Case II: Nonzero Exogenous Input

When , the model (3) has two nonnegative equilibria, and . The variational matrix of system (3) at is given by The characteristic equation of is .

So, is locally asymptotically stable if and only if and if , and then becomes unstable.

The interior equilibrium point of system (3) is given by where It can be provided that the unique interior equilibrium point of system (3) exists if and only if the following two conditions are satisfied: Summarizing the above analysis we come to the following theorem.

Theorem 3. If , then the equilibrium point of the system (3) is locally asymptotically stable if and only if and the unique interior equilibrium point of system (3) exists if and only if the following two conditions are satisfied:

3.3. Local Stability of

The variational matrix at is given by where The characteristic equation is where By the Routh-Hurwitz criterion [2830], it follows that all eigenvalues of characteristic equation have negative real part if and only if

Therefore, we come to the following theorem.

Theorem 4. is locally asymptotically stable if and only if the inequalities (24) are satisfied.

Furthermore, the conditions arising through Hopf bifurcation cannot be derived explicitly in terms of system parameters. Let us assume that stands for any parameter involved with the model system (3). If there exists a critical magnitude such that , , , and , then undergoes a Hopf bifurcation at (Liu’s criterion) [2830].

3.4. Global Stability of

is not always globally asymptotically stable. The conditions which guarantee the global stability of are stated in the following theorem.

Theorem 5. Since , , and are bounded, there exist positive constants () such that , , and . If the following inequalities hold: then is globally asymptotically stable.

Proof. We consider the following positive definite function which has a strict minimum at : Differentiating both sides with respect to along the solution of (3), we get (after some simple calculations) where Now, sufficient conditions for to be negative definite are that is,
Hence is a Lyapunov function and so is globally asymptotically stable [2830]. This competes the proof.

4. Model with Discrete Delay

It is already mentioned that time delay is an important factor in biological system. It is also reasonable to assume that the effect of the environmental toxicant on the population growth will not be instantaneous, but mediated by some discrete time lag required for incubation. So, we modify system (3) and (4) as follows: with initial conditions where is a nonnegative continuous functions on . For a biological meaning, we further assume that .

All parameters are the same as in system (3) except that the positive constant represents the activation period or reaction time of the toxicant in the population and is a constant.

The system (31) has the same equilibria as in art.3.2. The main purpose of this section is to study the stability behaviour of in the presence of discrete delay . We linearize system (31) by using the following transformations: Then transformed linear system is given by where Let us take the solution of the system (34) as , , . This leads to the following characteristic equation: where

A necessary condition for a stability change of is that the characteristic equation should have purely imaginary solutions. In order to see the delay induced instability, let us consider as bifurcation parameter and assume a purely imaginary solution of (3) in the form . Therefore, substituting in (37) and separating real and imaginary parts, we get Eliminating by squaring and adding (38), we get the equation for determining as where Substituting in (40), we get a cubic equation given by Now, from (38), we get Thus when , the characteristic equation (37) has a pair of purely imaginary roots . Let be a root of (37) such that the following two conditions hold: Differentiating both sides of (37) with respect to , we get Now, if is the first positive root of (41) and hence is the first positive root of (40), and then Thus, the transversality condition is satisfied and the steady state becomes unstable when . Moreover, a Hopf bifurcation [2830] occurs when passes through the critical value , where is the smallest positive value of given in (42).

We notice that, in (41), is continuous everywhere with and as . Therefore, the cubic equation (41) always has at least one positive root. Consequently the stability criteria of the system for will not necessarily ensure the stability of the system for . In the following theorem, we have given a criterion for switching the stability behaviour of .

Theorem 6 (see [2830]). If exists with the conditions (24) and be a positive root of (41), then there exists a such that(i) is locally asymptotically stable for ,(ii) is unstable for ,(iii)the system (31) undergoes a Hopf-bifurcation around at , and the minimum is taken over all positive such that is a solution of (41).

5. Estimation of the Length of Delay to Preserve Stability

We linearize the system (31) about its interior equilibrium and get where Taking Laplace transform of the system (47), we get where and , , and are the Laplace transform of , and , respectively.

Now, we will use the “Nyquist theorem” [21] which states that if is the arc length along a curve encircling the right half of the plane, then a curve will encircle the origin a number of times equal to the difference between the number of poles and the number of zeros of in the right half of the plane. Using “Nyquist theorem” [21, 31], it can be shown that the conditions for local asymptotic stability of are given by [32] where and is the smallest positive root of (52).

We have already shown that is stable in absence of delay (provided that inequalities given in (24) are satisfied). Hence, by continuity, all eigenvalues will continue to have negative real parts for sufficiently small provided that one can guarantee that no eigenvalue with positive real part bifurcates from infinity as increases from zero. This can be proved by using Butler’s lemma [32].

In this case, conditions (51) and (52) give

Now, (53) and (54), if satisfied simultaneously, are sufficient conditions to guarantee stability. We will utilize them to get an estimate on the length of delay. Our aim is to find an upper bound on , independent of so that (53) holds for all values of , and hence in particular at .

We rewrite (54) as

Maximizing subject to , , we obtain

Hence, if then clearly from (56) and (57), we have .

From the inequality (53), we get

As is locally asymptotically stable for (provided that inequalities given in (24) are satisfied), for sufficiently small , (58) will continue to hold. Substituting (55) into (58) and rearranging, we get Using the bounds we obtain

Now, from (59) and (60), we get where then stability is preserved for .

Summarizing the previous discussions we come to the following theorem.

Theorem 7. The delayed model (31) and (32) will be locally asymptotically stable at with the conditions (24) if the delay lies within the interval where is given by (63).

6. Direction of Hopf Bifurcation

In this section, we will derive the direction of the Hopf bifurcation [17, 21, 33] and sufficient conditions of the stability of bifurcating periodic solution from the interior equilibrium of the system (31) at the critical value . We will employ the approach of the normal form method and center manifold theorem introduced by Hassard et al. [33].

Let , , , , , and for . , and we can rewrite the system (31) and (32) as follows: where , are given, respectively, by where By the Riesz representation theorem, there exists a matrix function whose components are bounded variation function in such that In fact we can choose where is the Dirac delta function.

For , we define Then the system (64) is equivalent to For , where is the 3-dimensional space of row vectors, the adjoint operator of A is defined as For and , we define a bilinear inner product taking the first argument to be conjugate linear and the second to be linear of the form where . It is easy to verify that and are a pair of adjoint operators.

From the discussions in the previous section, we know that are the eigenvalues of . Hence, the eigenvalues of are .

Suppose that is eigenvector of corresponding to , and then we have It follows from (67), (68), and the definition of that where is identity matrix of order , that is, Solving (76), we have , where Now, suppose that is the eigenvector of corresponding to , and then we have It follows from (67), (68), and the definition of that, where is identity matrix of order , that is, Solving (80), we have , where Now, Thus in order to assume , we can choose Also since , we have Therefore, .

Now, we first compute the coordinates to describe the center manifold at . Let be the solution of (64) when . Define On the center manifold , we have where and are local coordinates of center manifold in the direction of and . Note that is real if is real. We only consider real solutions.

From (85), we get For the solution of (64) and since , we have where Since , we have where

Then,

It follows that where Since , then we obtain The comparison of the coefficients of , , , from (89) and (95) produces Clearly, involves the terms , , and we still need to compute , as which is rewritten as where On the other hand, on We derive the following inequalities from (98)–(100): For , it follows from (89), (97), and (98) that Computing this equation with (99), we get For , it follows from (69), (101), and (103) that Solving for , we get And similarly from (69), (101), and (104), it follows that where and are both three dimensional constant vectors and can be determined by setting in . This leads to Comparing the corresponding coefficients of the above equation with (99) for , we get From the definition of and (101) for , we get Notice that Substituting (106) and (111) into (109), we obtain which leads to that is, On simplification, we get where

Similarly substituting (107) and (112) into (110), we get which leads to that is, On simplification, we get where

Thus, and can be determined from (90) and (93). Furthermore, can also be determined. Therefore, based on the previous analysis, can be determined in terms of the parameters and delays involved in the system. Thus, we can compute the following values: where sign of determines the direction of the Hopf bifurcation, sign of determines the stability of the bifurcating periodic solution, and sign of determines the period of the bifurcating periodic solution in the center manifold at the critical value . Therefore, by the results of Hassard et al. [33], we can summarize the properties of the Hopf bifurcation at the critical value in the following theorem.

Theorem 8. In system (31) and (32), the following results hold.(i)If , the Hopf bifurcation is supercritical (subcritical).(ii)If , the bifurcated periodic solutions are stable (unstable).(iii)If , period of the bifurcating periodic solution increases (decreases).

7. Model with Double Delays

In this section we have discussed the dynamical behaviours of the model (3) and (4) with double discrete time delays and . Here, represents the activation period or reaction time of the toxicant in the population biomass; represents the activation period or reaction time of the toxicant in the population from the environment and is a constant. The modified model is as follows: with initial conditions where is a continuous function on and is a continuous function on . For a biological meaning, we further assume that .

Now, we will estimate the length of the delay which preserves the stability of the system (125a) and (125b). The corresponding linearized system of (125a) and (125b) about the interior equilibrium is given by where Taking Laplace transform of the system (126), we get where and are the Laplace transforms of , , and , respectively.

Now, using “Nyquist theorem” [21, 31], it can be shown that the conditions for local asymptotic stability of are given by [32] where and is the smallest positive root of (131).

Nyquist criterion [21, 31] implies that if is a solution of (131), then Next, we wish to find an upper bound for (say, ) independent of , such that (130) and (131) hold for all , and hence in particular for when .

Maximizing (133) with , , , , , we get If is a positive root of (134), then we obtain Then, clearly, .

Equation (130) indicates that the following inequality should hold: where If we can find such that where , then the Nyquist criterion holds. By (138), we can estimate the values of and We choose Now, We can choose From (138), we get Thus So, we come to the following theorem.

Theorem 9. The delayed model (125a) and (125b) will be locally asymptotically stable at with the conditions (24) if

Next, we will study the existence of Hopf bifurcation of system (125a) and (125b) by choosing one of the delays as a bifurcation parameter. Let us take as a bifurcation parameter. The characteristic equation corresponding to system (125a) and (125b) is Hence, is a simple root of (146) at . By separating real and imaginary parts, we obtain

Differentiating equation (146) implicitly with respect to , we obtain

Now, where Cleary, the sign of in (149) is positive. Because of the coherence of the sign of and , we get that if , then .

Therefore, transversality condition holds and Hopf bifurcation occurs at , provided that .

Summarizing the above analysis, we get the following result.

Theorem 10. If there exists such that (147) holds and , that is, , then a Hopf bifurcation occurs at as passes through .

8. Numerical Simulations

Analytical studies can never be completed without numerical verification of the derived results. In this section, we present computer simulation of some important analytical results discussed in the previous sections. Besides verification of our analytical findings, these numerical solutions are very important from a practical point of view.

We take the parameters of the system as given in Table 1. Using those values under zero exogenous input, that is, , it is observed that the equilibrium becomes unstable and interior equilibrium is not feasible. Only the axial equilibrium is locally asymptotically stable, which is shown in Figure 1. Therefore, Theorem 2 is verified. When (i.e., nonzero exogenous input) using the parameter values given in Table 2, it is observed that the equilibrium exists and is locally asymptotically stable with , but the interior equilibrium does not exist, which is shown in Figure 2. On the other hand, when using the parameter values given in Table 1, it is found that , , and becomes unstable and exists, which is shown in Figure 3. Therefore, Theorem 3 is verified.

Using the parameter values given in Table 1 with , it is observed that the conditions of Theorem 4 are satisfied as , , . Therefore, by Theorem 4, is locally asymptotically stable. Figure 3 shows that and approach to their steady-state values , , and , respectively, in finite time. The -plane, -plane, -plane, and -plane projections of the solution are shown in Figures 4, 5, 6, and 7, respectively. Here the conditions of Theorem 5 are also satisfied and consequently is globally asymptotically stable. The phase diagram is shown in Figure 8.

9. Discussions and Conclusions

In this paper, we have developed a single-species population model under the influence of toxicant in the population and in the environment by means of ordinary differential equations in terms of their concentrations. It is shown that the model system (3) with (4) is bounded and the solutions are positive, which in turn implies that the system is ecologically well behaved. The analysis of the existence of various equilibrium points under zero exogenous toxicant input and nonzero exogenous toxicant input leads us to Theorems 2 and 3, respectively. Theorem 2 states that under zero exogenous input the trivial equilibrium is unstable, the interior equilibrium is not feasible, and only the axial equilibrium is locally asymptotically stable. On the other hand, when exogenous input is nonzero, the existence of the interior equilibrium depends on some conditions and the equilibrium is locally asymptotically stable under some condition. The stability criteria for given in Theorems 4 and 5 are the conditions for locally and globally stable coexistence of the population biomass, population toxicant, and environmental toxicant.

It is proved by several researchers that the effect of time delay must be taken into account in the population models to make them ecologically more meaningful and useful. It seems also reasonable that the effect of the environmental toxicant on the organismal activities and the population growth will not be instantaneous. It must be mediated with some discrete time lag required for incubation. From this point of view, we have formulated single and double time delayed model given in (31) and (125a) and (125b), respectively, where the delays may be looked upon as the reaction time of the environment toxicant on the population biomass and the population toxicant.

The analysis of the delayed models (31) and (125a) and (125b) shows more complicated behaviour than their nondelayed counterpart. The rigorous analysis of the single delayed model (31) leads us to Theorem 6 which mentions that the stability criteria in absence of delay are no longer enough to guarantee the stability in the presence of delay, but rather there is a critical value of the gestation delay such that the system is stable for and becomes unstable for at the interior equilibrium . From the next part of our analysis, we obtain the estimated length of delay to preserve stability at and the direction of Hopf bifurcation by using normal form theory and center manifold theorem. The stability of the double delayed model (125a) and (125b) is investigated by the Nyquist criteria which leads us to Theorem 9. By choosing one of the delays as a bifurcation parameter, the model (125a) and (125b) is found to undergo a Hopf bifurcation under some conditions stated in Theorem 10.

Analytical studies can never be completed without numerical verification of the results. So, some numerical simulations for justifying the important theoretical results are also illustrated by using MATLAB, which shows the mathematical and ecological reliability of the proposed model.