Discrete Dynamics in Nature and Society
Volume 2015 (2015), Article ID 756315, 13 pages
http://dx.doi.org/10.1155/2015/756315

Research Article

On the Dynamical Behavior of Toxic-Phytoplankton-Zooplankton Model with Delay

Mehbuba Rehim, Weixin Wu, and Ahmadjan Muhammadhaji

College of Mathematics and System Sciences, Xinjiang University, Urumqi, Xinjiang 830046, China

Received 11 October 2014; Accepted 5 January 2015

Academic Editor: Gualberto Solís-Perales

Copyright © 2015 Mehbuba Rehim et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

A toxin producing phytoplankton-zooplankton model with inhibitory exponential substrate and time delay has been formulated and analyzed. Since the liberation of toxic substances by phytoplankton species is not an instantaneous process but is mediated by some time lag required for maturity of the species and the zooplankton mortality due to the toxic phytoplankton bloom occurs after some time laps of the bloom of toxic phytoplankton, we induced a discrete time delay to both of the consume response function and distribution of toxic substance term. Furthermore, based on the fact that the predation rate decreases at large toxic-phytoplankton density, the system is modelled via a Tissiet type functional response. We study the dynamical behaviour and investigate the conditions to guarantee the coexistence of two species. Analytical methods and numerical simulations are used to obtain information about the qualitative behaviour of the models.

1. Introduction

Phytoplankton are one of the most important components of the marine ecosystem. They not only form a basis for all aquatic food chains but also perform a very useful service by producing a huge amount of oxygen for human and other living animals after absorbing carbon-dioxide from surrounding environments [1]. Zooplankton are microscopic animals that eat other plankton and serve as a most favorable food source for fish and other aquatic animals. During the recent years, many authors have studied the system between zooplankton and phytoplankton. Authors in [2] have dealt with a nutrient-plankton model in an aquatic environment in the context of phytoplankton bloom. In [3] the effect of seasonality and periodicity on plankton dynamics is investigated. In [4], two plankton ecosystem models with explicit representation of viruses and virally infected phytoplankton are presented.

The most common features of the phytoplankton population is rapid increase of biomass due to rapid cell proliferation and almost equally rapid decrease in populations, separated by some fixed time period. This type of rapid change in phytoplankton population density is known as “bloom” [5]. Due to the accumulation of high biomass or to the presence of toxicity, some of these blooms, more adequately called “harmful algal blooms” [6], are noxious to marine ecosystems or to human health and can produce great socioeconomic damage. There has been a global increase in harmful plankton blooms in last two decades [7–9].

Because of the difficulty in measuring plankton biomass, mathematical modeling of plankton population is an important alternative method of improving our knowledge of the physical and biological processes relating to plankton ecology [10–15]. In [16], nutrient-plankton-zooplankton interaction models with a toxic substance which inhibits either the growth of phytoplankton, zooplankton, or both trophic levels are proposed and studied. In [17], the authors have constructed a mathematical model for describing the interaction between a nontoxic and toxic phytoplankton with a single nutrient.

Based on the fact that the release of toxin from phytoplankton species is not an instantaneous process but is mediated by some time lag required for maturity of the species and the zooplankton may die after some time lapse of the bloom of toxic phytoplankton (see http://www.mote.org/, http://www.mdsg.umd.edu/), models incorporating time delay in diverse biological models are extensively reviewed by Beretta and Kuang [18], Gopalsamy [19], Cooke and Grossman [20], and Cushing [21]. The discrete time delay has potential to change the qualitative behavior of the dynamical systems [22–29]. Chattopadhyay et al. proposed a delay model incorporating time lag in toxin liberation by phytoplankton to avoid predation by zooplankton [28]. In [28], the authors introduced distribution delay and discrete to toxin liberation term. Due to discrete time delay in toxin liberation, the local existence of periodic solution through Hopf bifurcation has been obtained in [28].

In [5, 10, 28–31], the following plausible toxic-phytoplankton-zooplankton system has been studied:

where is the density of the toxin-producing phytoplankton (TPP) population and is the density of zooplankton population at any instant of time . In model (1), is the intrinsic growth rate and is the environmental carrying capacity of TPP population. The term describes the functional response for the grazing of phytoplankton by zooplankton, and is the maximum uptake rate for zooplankton. denotes the ratio of biomass conversion and is the natural death rate of zooplankton. The function represents the distribution of toxic substance which ultimately contributes to the death of zooplankton populations, and the parameter denotes the rate of toxic substances produced by per unit biomass of phytoplankton.

Model (1) has been studied for the following cases: when is linear but is Holling type II [10, 28, 29] or Holling type III [10]; when and are both Holling type II [5, 10] or Holling type III [10]; when is Holling type II while is Holling type III [10].

In above cases, and are both increasing functions of over the entire interval . However, in some cases, very high substrate concentrations in the lakes actually inhabit the growth of phytoplankton cells. Moreover, with the substrate concentrations increasing unlimitedly, some kind of microorganism will die eventually [32]. To describe the above phenomenon accurately, we consider from a different point of view. Adopting the idea used in [32], we assume that there exists a constant such that is increasing over the interval and is decreasing on the interval . More precisely, we use the so-called Tissiet functional response of the form of (see, e.g., [32]). This type of functional response takes care of the fact that the predation rate decreases at large toxic-phytoplankton density.

The remainder of this paper is organized as follows. The model is described in Section 2. In Section 3, we state and prove the positivity and boundedness of the solutions. Then, in Section 4, equilibria, their existence, and local asymptotic stability are considered. Under the aids of numerical simulation method, we further analyse the model and determine if there is a parameter range for the delay parameter , where more complicated dynamics occur. In Section 5, the permanence of the system is discussed by some analytic techniques on limit sets of differential dynamical systems. Finally, a brief discussion is presented in Section 6.

2. State of the Model

In a real ecological context, the interaction between phytoplankton and zooplankton will not be essentially instantaneous. Instead, the response of zooplankton to contacts with phytoplankton is likely to be delayed due to a gestation period. Another fact, during the interaction between phytoplankton and zooplankton, is that the liberation of toxic substances by phytoplankton must be mediated by some time lag which is required for the maturity of toxic-phytoplankton. Let stand for the time delay in conversion of food to viable biomass for the species, and is the discrete time period required for the maturity of phytoplankton cells to liberate toxic substances. Based on model 3 in [5], we intend to study a model system with the assumption that and are described by same type of function, namely, Tissiet functional response. Then we arrive at the following model:

where the parameter is the intrinsic growth rate and is the environmental carrying capacity of TPP population. The constant is the maximum per capita grazing rate, denotes the ratio of biomass conversion, denotes the rate of toxic substances produced by per unit biomass, is the natural death rate of the zooplankton. All the parameters in system (2) are positive constants with their usual ecological meanings.

Here we observe that if there is no delay (i.e., ) and , then . Hence, throughout our analysis, we assume that .

For the sake of simplicity of mathematical analysis, in this paper, we consider model (2) for the special case: the time delay in conversion of food to viable biomass for the species equal to the discrete time period required for the maturity of phytoplankton cells to liberate toxic substances; that is, .

By performing the following scaling for model (2):

we obtain a dimensionless system in the state variables , , which can be written, by removing the stars, in the form

with initial conditions:

where .

3. Positive and Boundedness

In this section, we consider the positive and boundedness of the solutions of model (4) with initial condition (5). One has the following theorem.

Theorem 1. Let , on and , . Then

(a)all the solutions of (4) with initial condition (5) exist on for some constant and are unique and positive for ,(b), , where ,(c)if , then , where . Further, the subset

is positively invariant with respect to (4).

Proof. By Theorems  2.1 and 2.3 in Hale and Lunel [33], solutions of (4) with initial data exist on for some and are unique. Suppose is a solution of (4) for . Without loss of generality, assume that is the maximum internal of the solution and if the solution exists for any . Integrating the first equation of (4) gives

To prove the for any , use the method of contradiction. Suppose there exists a such that

From the second equation of system (4), we have

This contradicts . Hence, for all . This completes the proof of conclusion (a) in Theorem 1.
It follows from the first equation of (4) that , which implies that . Define , . Then from (4) we obtain

Applying the theorem of differential inequality we obtain that

Therefore, . Thus, there is a constant , such that for all . This completes the proof of conclusion (b) in Theorem 1.
From the first equation of system (4) we get

which implies that .
For any , let be the solution of (4) with the initial function . If there is a such that , then for some and . Hence, it follows from the first equation of (4) that

which is a contradiction to . So, for all .
It is easy to prove that if there is a such that , then for all one has . Which implies that for all . This completes the proof of conclusion (c) in Theorem 1.

Remark. The conclusions (b) and (c) in Theorem 1 indicate that if the ratio of biomass conversion is less than certain value, then phytoplankton population will be persistent.

4. Equilibrium, Stability, and Hopf Bifurcation

4.1. Existence of Equilibria

It is easy to see that model (4) has two boundary equilibria and . To discuss the existence of the positive equilibria, we work on the equation

Denote the left-hand side of the second equation in (14):

Then , ,

Let , we have , where

From this we know that the function is monotone increasing in the interval and monotone decreasing in the interval . Hence, reaches its maximum on the interval . Therefore, we have the following results.

Theorem 2. The following statements hold.

(1)If , then equation has no roots in the interval (Figure 1(a)). In this case, for system (4), there is no positive equilibrium.(2)If and , then equation has only one root, , in the interval , where (Figure 1(c)). In this case, for system (4), there exists a single positive equilibrium .(3)If and , then equation has two distinct roots, and , in the interval , where (Figure 1(d)). In this case, for system (4), there are two distinct positive equilibria and .(4)If , then equation has a unique root, , in the interval , (Figure 1(b)). In this case, for system (4), there is a unique positive equilibrium .(5)If and , then equation has two distinct roots, and , in the interval , where (Figure 1(e)). In this case, for system (4), there is also unique positive equilibrium ,

where , , are determined by the first equation in (14).

(a)

(b)

(c)

(d)

(e)

Figure 1: The graph of function .

4.2. Local Stability and Hopf Bifurcation

Let denote any one of the equilibrium points of system (4). Linearizing (4) about we obtain

In what follows, existence of the interior equilibria, dynamical properties of the zooplankton-free equilibrium, and the interior equilibrium are investigated.

Theorem 3. For any time delay , is locally asymptotically stable if and is unstable if ; the trivial solution of the linearized system of (4) about is stable for .

Proof. From (18) we know that, at any equilibrium point , the characteristic equation is

where

in which .
We first consider local stability of . For , (19) becomes

One of the roots of (21) is , and other roots are given by solution of the following equality:

It follows from [34] that

(i)if , then all roots of (22) have negative parts for any time delay . Hence, is locally asymptotically stable for any time delay ;(ii)if , then (22) has roots which have positive real parts for any time delay . Hence, is unstable for any time delay ;(iii)if , this is a critical case and (22) is equivalent to

From [34], it has that except , any root of (23) has negative real part for any time delay . Hence, the trivial solution of the linearized system of (4) about is stable for any time delay .

Theorem 4. Suppose that the positive equilibria , , of system (4) exist. For system (4), the following statements hold.

(i)If and , then are locally asymptotically stable for , are unstable for , and there is a periodic solution around for . If , then is unstable for any , where is determined by (30), .(ii) is unstable for any .(iii)If , then the trivial solution of the linearized system of (4) about is stable for any . If , then the trivial solution of the linearized system of (4) about is unstable for any ,

where , .

Proof. When , the transcendental equation (19) becomes

where

if , and if , .
By the argument in Section 4.1, as long as exists, it must be , ; thus, ; as long as exists, it must be ; thus, ; as long as exists, it must be ; thus, ; hence, from Routh-Hurwitz theorem, we have that if , , then and are locally asymptotically stable. If , , then and are unstable; if , , then and are nonhyperbolic equilibria; is unstable; is a critical case.
Suppose that , , is a root of (19) for some . Substituting in the characteristic equation (19), and separating real and imaginary parts, then we get

Eliminating from (26), we get a biquadratic equation in as

Its roots are

By (20) we obtain

We consider the stability of ,  .
Since , if ,  , then there is only one root, ,  , such that

That is, (19) has one imaginary solution, ,  . From (26), we obtain the following set of for which there are imaginary root:

where and

If , then , when , , are locally asymptotically stable. Hence, if and , then are locally asymptotically stable for (,  ).
Differentiating (19) with respect to , we obtain

It follows that

From (35), we have

Thus, we have

Hence, there is a Hopf bifurcation at , . Therefore, if and , then is locally asymptotically stable for , is unstable for , and there is a periodic solution around for . If and or , then is unstable for any ,  .
We consider the stability of .
Since , if , then (19) has only one imaginary solution. From proof of , . Hence, is unstable for any . If and , then (19) has only one imaginary solution. If and , then (19) does not have imaginary solution. Wnen ,   is unstable. Therefore, is unstable for any .
We consider the stability of .
Consider the characteristic equation (19). is a root of (19) for all . Suppose is a root of (19), then substituting in the characteristic equation, we have

and separating real and imaginary parts, then we get

From (39), we obtain

(i) Assume . Suppose (19) has a root ,   for all . Then, from (40), we have

The left-hand size of (41) is equal to

Obviously,

and by assumption,

Therefore,

According to (41), this is impossible, since we are assuming . Hence, all roots of (19) have nonpositive real parts, and this implies that the trivial solutions of the linearized system of (4) about are stable.
(ii) Assume . Consider the following real function:

It is clear that and . There exists a such that if , , we also have

Since yields , which gives , hence, there exists a such that when , . Therefore, there exists at least one , , such that , that is, (19) has at least a positive root. This indicates that trivial solution of the linearized system of (18) about is unstable.

4.3. Numerical Simulation

For numerical simulation, let us first consider case (i) in Theorem 4. We choose parameter values as , , , and . For this choice, , . The stability of the interior equilibrium point depends upon the magnitude of delay . From (30) and (31) we obtain , . Consider , . When , Figures 2(a) and 2(b) show that the system is approaching the equilibrium point . When increases to , becomes unstable (Figure 2(c)). Stable oscillations appear when and the Hopf bifurcation periodic solutions are depicted in Figure 3.

(a)

(b)

(c)

Figure 2: ((a)-(b)) Behavior and phase portrait of system (4) with . The positive equilibrium is stable. (c) The oscillations of and versus time when show that the positive equilibrium is unstable. Here and the initial value is .

Figure 3: Behavior and phase portrait of system (4) when . Hopf bifurcation occurs from the interior equilibrium . The initial value is .

For case (ii) in Theorem 4, consider the following set of parametric values: , , , and . Then we have , . By Theorem 4, the system (4) has two equilibria and . From direct calculation, we get that , , , , , , and . For , the numerical simulation shows that equilibrium point is stable (Figure 4). When increases to , the equilibrium becomes unstable (Figure 6). There is a periodic solution around the equilibrium when (Figure 5).

Figure 4: Behavior and phase portrait of system (4) with shows the equilibrium is stable. Here and the initial value is .

Figure 5: Behavior and phase portrait of system (4) with . Hopf bifurcation occurs from the interior equilibrium . The initial value is .

Figure 6: Behavior and phase portrait of system (4) with . The interior equilibrium is unstable. Here and the initial value is .

For case (iii) in Theorem 4, let us take the parametric values , , , and . For these set of values, by direct calculation, we get , , and . By Theorem 1, system (4) has only one equilibrium , and from Theorem 4 we know that the equilibrium is stable (Figure 7(a)). Take parametric values as , , , and . Then we have , , , and . By Theorem 2, system (4) has only one equilibrium . The numerical simulation shows is unstable (Figure 7(b)).

(a)

(b)

Figure 7: Behavior and phase portrait of system (4). (a) The equilibrium is stable for any when . (b) The equilibrium is unstable for any when . The initial value is .

5. Permanence

In this section, we will use the same methods as [35] to prove the permanence of system (4).

Theorem 5. Suppose that , , and . Then, for any time delay , system (4) is permanent.

Proof. From Theorem 1, we see that it is enough to consider the solution with initial condition . Theorem 1 implies that the omega limit set of is nonempty, compact, and invariant and . It follows from definition of permanence and Theorem 1, that we only need to show

Here is some positive constant which dose not depend on the initial function .
Let us first show

In fact, if (49) is not true, then, from , we see that there exists a positive time sequence such that

Note that the solution is bounded on by Theorem 1. It follows from (4) that is uniformly continuous on . Hence, from Ascoli’s theorem there is a subsequence of , still denoted by , such that

holds uniformly on in the wide sense. From Theorem 1, we have that for any and that, for any , the function of is the solution of (4) with the initial function . Here we note that and for any .
We claim that for any . From Theorem 1, we know that if , then the solution of (4) exists and and . Thus, from , we have that for any . Then, it follows from (4) that for any , and for any . Hence,

Since is bounded for , we must have , which implies that for . It follows from (4) and the invariant of that for any . This shows that above claim holds. In particular, we have that

For sufficiently small and sufficiently large , , we have

Hence, for large , from (50) and above inequality, we obtain

which is a contradiction to . This completes the proof of .
Next, let us show

For any initial functions sequence , let be the solution of (4) with the initial function . Let be the omega limit set of . We have that there exists some compact and invariant set such that as . Here, means Hausdor distance [35, 36].
If (48) does not hold, for some initial function sequence such that , we have that there is some such that for some . Now, let be the solution of (4) with the initial function . Then, by the invariant of , we have that for all and . Note and the positivity of all solutions, we easily have that for all . Hence, it follows from (4) that , , and for all . This implies that , , for all . If , we see that the negative semiorbit is unbounded. This is contradiction.
If , we have that , for all . This shows that . Let us show is factually isolated [35, 36]. that is, there exists some neighborhood of in such that is the largest invariant set in . In fact, let us choose

for some sufficiently small positive constant and . We will show that is the largest invariant set in for some .
If not, for any sufficiently small , there exists some invariant set such that is not empty. Let and be the solution of (4) with the initial function . Then, for all .
If , by the invariance of and Theorem 1, we also have the contradiction that or that the negative semiorbit of (4) through is unbounded.
If , from the Theorem 1, we see that for all . Now, let us consider the continuous function:

for some constant . Because of , we have . The time derivative of along the solution satisfies

Since and , we have that . We can choose , such that . From (49), there exists a constant such that for some constant and . Hence, it follows from (59)

for all . Thus, as . This contradicts Theorem 1 and shows that is isolated.
We easily see that the semigroup defined by the solution of (4) satisfies the conditions of Lemma  4.3 in [36] with . Thus, by Lemma  4.3 in [36], we have that there is some , such that . Here, denotes the stable set of .
If , again by the invariant of and Theorem 1, we also have the contradiction that or that the negative semiorbit of (4) through is unbounded.
If , from Theorem 1, we see that , for all . It follows from that , , which contradicts (49). This shows that (48) holds. Thus, (4) is permanent. This completes the proof of Theorem 5.

6. Discussion

In a real ecological context, the interaction between phytoplankton and zooplankton will not be essentially instantaneous. Instead, the response of zooplankton to contact with phytoplankton is likely to be delayed due to a gestation period. Another fact, during the interaction between phytoplankton and zooplankton, is that the liberation of toxic substances by phytoplankton must be mediated by some time lag which is required for the maturity of toxic-phytoplankton. And in some cases, very high substrate concentrations in the lakes actually inhabit the growth of phytoplankton cells. Moreover, with the substrate concentrations increasing unlimitedly, some kind of microorganism will die eventually [32]. Based on the above fact, in this paper, we introduce time delay to a phytoplankton-zooplankton interaction model with exponential substrate uptake and exponential distribution of toxic substance term. The model (4) accounts for some natural phenomenon. By using comparison principle for functional differential equations and traditional analysis technique for transcendental equations [34], we give a detailed analysis on boundedness of solutions of system (4) and local asymptotic stability of the equilibria of system (4). Our results show that time delay is factually harmless for the local asymptotic stability of the zooplankton free equilibrium of (4), but it is not always harmless for the positive equilibrium; that is to say, because of the time delay the positive equilibrium becomes unstable (Theorem 4). Based on some known techniques on limit sets of differential dynamical systems, we show that, for any time delay, the phytoplankton-zooplankton interaction model is permanent if and only if one positive equilibrium exists.

Under the aids of the numerical simulation we further investigate the delayed model system (4). Figures 1–6 illustrate that the delay plays crucial role in determining the asymptotic behavior of solutions of model (4). The stability or oscillatory coexistence depends upon not only the parametric restriction but also upon the gestation delay (liberation delay). Our system exhibits stable behavior when , where is the threshold value of the parameter . This threshold value of is determined by our numerical simulations. This value of is 3.5 for case (i) of Theorem 4 and is 5.1 for case (ii). These stable solutions are shown in Figures 2(a), 2(b), and 4. When crosses , there is a delay induced instability demonstrated in Figures 2(c) and 5 for . Stable oscillations appear when and the Hopf bifurcation periodic solution occurs (see Figures 3 and 5). Thus, there is a range of gestation delay (liberation delay) which initially imparts stability, then induces instability, and ultimately leads to periodic behavior.

Conflict of Interests

All authors declare that there is no conflict of interests regarding the publication of this paper.

Acknowledgments

The authors are grateful to the anonymous reviewers for their careful reading and valuable suggestions that greatly improved the quality of this paper. This work is supported by the National Natural Science Foundation of China (Grants no. 11261058, 11261056, 11271312, and 11361059).

References

J. Duinker and G. Wefer, “Das CO₂-problem und die Rolle des Ozeans,” Naturwissenschaften, vol. 81, no. 6, pp. 237–242, 1994.

B. Mukhopadhyay and R. Bhattacharyya, “Modelling phytoplankton allelopathy in a nutrient-plankton model with spatial heterogeneity,” Ecological Modelling, vol. 198, no. 1-2, pp. 163–173, 2006.

M. Gao, H. Shi, and Z. Li, “Chaos in a seasonally and periodically forced phytoplankton-zooplankton system,” Nonlinear Analysis: Real World Applications, vol. 10, no. 3, pp. 1643–1650, 2009.

C. J. Rhodes, J. E. Truscott, and A. P. Martin, “Viral infection as a regulator of oceanic phytoplankton populations,” Journal of Marine Systems, vol. 74, no. 1-2, pp. 216–226, 2008.

T. Saha and M. Bandyopadhyay, “Dynamical analysis of toxin producing phytoplankton-zooplankton interactions,” Nonlinear Analysis. Real World Applications, vol. 10, no. 1, pp. 314–332, 2009.

T. J. Smayda, “What is a bloom? A commentary,” Limnology and Oceanography, vol. 42, no. 5, part 2, pp. 1132–1336, 1997.

D. M. Anderson, “Toxic algae blooms and red tides: a global perspective,” in Red Tides, Biology, Environmental Science and Toxicology, pp. 11–21, Elsevier, New York, NY, USA, 1998.

G. M. Hallegraeff, “A review of harmful algal blooms and their apparent global increase,” Phycologia, vol. 32, no. 2, pp. 79–99, 1993.

T. Smayda, “Novel and nuisance phytoplankton blooms in the sea: evidence for a global epidemic,” in Toxin Marine Phytoplankton, pp. 29–40, Elsevier, New York, NY, USA, 1990.

J. Chattopadhayay, R. R. Sarkar, and S. Mandal, “Toxin-producing plankton may act as a biological control for planktonic blooms—field study and mathematical modelling,” Journal of Theoretical Biology, vol. 215, no. 3, pp. 333–344, 2002.

P. J. Wallhead, A. P. Martin, and M. A. Srokosz, “Spatially implicit plankton population models: transient spatial variability,” Journal of Theoretical Biology, vol. 253, no. 3, pp. 405–423, 2008.

R. R. Sarkar and J. Chattopadhayay, “The role of environmental stochasticity in a toxic phytoplankton-non-toxic phytoplankton-zooplankton system,” Environmetrics, vol. 14, no. 8, pp. 775–792, 2003.

J. Chattopadhyay, R. R. Sarkar, and A. El Abdllaoui, “A delay differential equation model on harmful algal blooms in the presence of toxic substances,” IMA Journal of Mathematics Applied in Medicine and Biology, vol. 19, no. 2, pp. 137–161, 2002.

A. Y. Morozov, S. V. Petrovskii, and N. P. Nezlin, “Towards resolving the paradox of enrichment: the impact of zooplankton vertical migrations on plankton systems stability,” Journal of Theoretical Biology, vol. 248, no. 3, pp. 501–511, 2007.

B. Mukhopadhyay and R. Bhattacharyya, “Role of gestation delay in a plankton-fish model under stochastic fluctuations,” Mathematical Biosciences, vol. 215, no. 1, pp. 26–34, 2008.

S. R.-J. Jang, J. Baglama, and J. Rick, “Nutrient-phytoplankton-zooplankton models with a toxin,” Mathematical and Computer Modelling, vol. 43, no. 1-2, pp. 105–118, 2006.

S. Roy, “The coevolution of two phytoplankton species on a single resource: allelopathy as a pseudo-mixotrophy,” Theoretical Population Biology, vol. 75, no. 1, pp. 68–75, 2009.

E. Beretta and Y. Kuang, “Geometric stability switch criteria in delay differential systems with delay dependent parameters,” SIAM Journal on Mathematical Analysis, vol. 33, no. 5, pp. 1144–1165, 2002.

K. Gopalsamy, “Harmless delays in model systems,” Bulletin of Mathematical Biology, vol. 45, no. 3, pp. 295–309, 1983.

K. L. Cooke and Z. Grossman, “Discrete delay, distributed delay and stability switches,” Journal of Mathematical Analysis and Applications, vol. 86, no. 2, pp. 592–627, 1982.

J. M. Cushing, Integrodifferential Equations and Delay Model in Population Dynamics, Springer, New York, NY, USA, 1977.

A. Hastings, “Delays in recruitment at different trophic levels: effects on stability,” Journal of Mathematical Biology, vol. 21, no. 1, pp. 35–44, 1984.

A. Martin and S. Ruan, “Predator-prey models with delay and prey harvesting,” Journal of Mathematical Biology, vol. 43, no. 3, pp. 247–267, 2001.

S. G. Ruan, “The effect of delays on stability and persistence in plankton models,” Nonlinear Analysis: Theory, Methods & Applications, vol. 24, no. 4, pp. 575–585, 1995.

Y. Xiao and L. Chen, “Modeling and analysis of a predator-prey model with disease in the prey,” Mathematical Biosciences, vol. 171, no. 1, pp. 59–82, 2001.

R. Mehbuba and I. Mudassar, “Dynamical analysis of a delay model of phytoplankton-zooplankton interaction,” Applied Mathematical Modelling, vol. 36, no. 2, pp. 638–647, 2012.

Y. Lv, J. Cao, J. Song, R. Yuan, and Y. Pei, “Global stability and Hopf-bifurcation in a zooplankton-phytoplankton model,” Nonlinear Dynamics, vol. 76, no. 1, pp. 345–366, 2014.

J. Chattopadhyay, R. R. Sarkar, and A. El Abdllaoui, “A delay differential equation model on harmful algal blooms in the presence of toxic substances,” Journal of Mathematics Applied in Medicine and Biology, vol. 19, no. 2, pp. 137–161, 2002.

S. Gakkhar and A. Singa, “Effects of the delay and seasonality on toxin producing phytoplankton-zooplankton system,” International Journal of Biomathematics, vol. 5, no. 5, Article ID 1250047, 21 pages, 2012.

R. R. Sarkar, S. Pal, and J. Chattopadhyay, “Role of two toxin-producing plankton and their effect on phytoplankton-zooplankton system—a mathematical study supported by experimental findings,” BioSystems, vol. 80, no. 1, pp. 11–23, 2005.

R. R. Sarkar, B. Mukhopadhyay, R. Bhattacharyya, and S. Banerjee, “Time lags can control algal bloom in two harmful phytoplankton-zooplankton system,” Applied Mathematics and Computation, vol. 186, pp. 445–459, 2007.

G. Fu, W. Ma, and S. Ruan, “Qualitative analysis of a chemostat model with inhibitory exponential substrate uptake,” Chaos, Solitons & Fractals, vol. 23, no. 3, pp. 873–886, 2005.

J. Hale and S. V. Lunel, Introduction to Functional Differential Equations, Springer, New York, NY, USA, 1993.

Y. Kuang, Delay Differential Equations with Applications in Population Dynamics, Academic Press, New York, NY, USA, 1993.

G. Butler, H. I. Freedman, and P. Waltman, “Uniformly persistent systems,” Proceedings of the American Mathematical Society, vol. 96, no. 3, pp. 425–430, 1986.

J. K. Hale and P. Waltman, “Persistence in infinite-dimensional systems,” SIAM Journal on Mathematical Analysis, vol. 20, no. 2, pp. 388–395, 1989.