Dynamics Analysis of a Nutrient-Plankton Model with a Time Delay

We analyze a nutrient-plankton system with a time delay. We choose the time delay as a bifurcation parameter and investigate the stability of a positive equilibrium and the existence of Hopf bifurcations. By using the centermanifold theorem and the normal form theory, the direction of the Hopf bifurcation and the stability of the bifurcating periodic solutions are researched. The theoretical results indicate that the time delay can induce a positive equilibrium to switch from a stable to an unstable to a stable state and so on. Numerical simulations show that the theoretical results are correct and feasible, and the system exhibits rich complex dynamics.


Introduction
Plankton is the basis of all aquatic food chains and its importance for marine ecosystems is widely recognized [1].The dynamics of zooplankton-phytoplankton systems have been discussed by many authors [2][3][4][5][6][7][8].A wealth of studies have shown that a time delay may have a complex effect on the dynamics of such systems, with effects that include stability switches for equilibria, the existence of Hopf bifurcation periodic solutions, and the direction and stability of Hopf bifurcation [9][10][11][12][13][14][15][16][17][18].
Algal blooms are a common feature of marine ecosystems, and a high nutrient concentration has an important influence on algal blooms [19].Plankton-nutrient interaction models can provide a better understanding of plankton dynamics.The model formulated by Huppert et al. [20] consists of two variables, nutrient level  and phytoplankton biomass , in the following system: Ṅ =  −  ()  − , Ṗ =  ()  − . ( Considering the biological significance, all the constants, , , , , and , are assumed to be positive.The phytoplankton growth rate  can be represented as a periodic function () = ( + ), where  is the forcing period affected by seasonal conditions such as light, salinity, and water temperature.
Harmful algal blooms have become a serious environmental problem worldwide and have been widely studied [8,21].On the basis of field studies and mathematical modeling, Chattopadhayay et al. [22] proposed the following phytoplankton-zooplankton model: ( The authors studied the existence and local stability of steady states and the existence of Hopf bifurcation for system (2) by taking various combinations of () and () [22].The results showed that toxin-producing plankton is helpful in terminating planktonic blooms by decreasing the grazing pressure of zooplankton.
To reflect maturation time, capture time, and other factors, a time delay is often included in mathematical models of population dynamics.Incorporation of a time delay can provide a better understanding of the dynamics of biological models.In recent years, many researchers have studied the fact that time delay plays important roles in biological dynamical systems [18,[23][24][25][26][27][28].Following Huppert et al. [20] and Chattopadhayay et al. [22], we consider the following plausible nutrient-phytoplankton-zooplankton system with a time delay: Ṅ =  −  −  ≜  1 (, , ) , Ṗ =  ( − )  −  ℎ +  −  −  2 ≜  2 (, , ) , where  is the nutrient concentration,  and  are the density of the phytoplankton and zooplanktons population, respectively,  is the external source of nutrients flowing into the system,  is the small loss rate to reflect sinking of nutrients from the epilimnion down to the hypolimnion,  and  are the capture rates for phytoplankton on nutrient and zooplankton on phytoplankton, and  and  denote the rates of biomass conversion.ℎ is the half-saturation constant,  and  are the death rates for phytoplankton and zooplankton, and  is the interspecies competition coefficient for phytoplankton, which reflects phytoplankton self-limitation.The parameter  ≥ 0 denotes the time delay, which arises because of the time required by phytoplankton to absorb nutrients.The remainder of the paper is organized as follows.In Section 2, we determine the stability of the positive equilibrium and the existence of Hopf bifurcation.Section 3 illustrates the direction of Hopf bifurcation and the stability of bifurcating periodic solutions.To verify the theoretical analysis, we present numerical simulations in Section 4. Finally, Section 5 provides some conclusions.

Stability and Existence of Hopf Bifurcation
In this section, we study the local stability of the positive equilibrium of system (3).The positive equilibrium  * ( * ,  * ,  * ) corresponds to the coexistence of the three species, where It is obvious that  * exists if and only if From ( 5), we know that the positive equilibrium  * of system (3) exists if external nutrients flowing into the system exceed a certain critical value, and the ratio of biomass consumed by zooplankton is greater than the sum of the zooplankton mortality rate and the rate of toxic substance production by phytoplankton.Hereafter we assume that the conditions stated in (5) always hold.
Using the translations we translate the positive equilibrium  * to the origin.Then the linearized system for system (3) near ( * ,  * ,  * ) is where The associated characteristic equation of system (7) is where It is well known that the positive equilibrium  * of system (3) is locally asymptotically stable if all roots of (9) have negative real parts and is unstable if (9) has a root with positive real parts.Now, we shall discuss the distribution of the roots of (9).First, when  = 0, the characteristic equation becomes According to the Routh-Hurwitz criterion,  * is locally asymptotically stable if and only if By Corollary 2.4 in the paper of Ruan and Wei [29], we know that if instability occurs for a particular value of time delay, a characteristic root of (9) must intersect the imaginary axis.Hence, suppose that  ( > 0) is a purely imaginary root of (9).Inserting this into (9) and separating the real and imaginary parts, we obtain Squaring and adding these two equations, we have where  1 =  2 − 2 and  2 =  2 − 2 −  2 .Let  =  2 ; then (14) becomes Form (15), we define a function When Δ 1 > 0, the equation   () = 0 has two real roots: Noticing that lim  → +∞ () = +∞ and (0) =  2 > 0, then we introduce the following results, which have been proved by many authors [30,31].

Lemma 2. If (17) holds, then
Proof.Differentiating both sides of ( 9) with respect to , we have Then where As (15) has two positive roots, denoted by  1 ,  2 , and ( 0 ) is the local minimum value.Assuming  1  <  2  , we obtain   ( 1 ) > 0 and   ( 2 ) < 0. Hence, On the basis of the above analysis, we have the following theorem.
Theorem 3. Suppose that ( 5) and ( 12) hold; then we obtain the following: 3) is locally asymptotically stable for all values of  ≥ 0.

Numerical Simulations
In this section, we verify the theoretical results proved in previous sections using numerical simulations for the following parameter values:

(63)
According to the analysis above, the existence and stability of the interior equilibrium for system (3) both change with the value of the external source of nutrients , as illustrated in Figure 1.In Figure 1(a),  (the conversion frequency of stability switches for  * ) changes with .From (5), the condition for the existence of  * is  > 0.1416.Region I in Figure 1(a) indicates that system (3) has no positive equilibrium if the external source of nutrients flowing into the system is less than the critical value, 0.1416.In region IV,  * is unstable when condition (12) is not satisfied.This shows that when  is greater than a critical value, 0.755,  * is always unstable and a time delay will no longer affect system (3).When the value of  is in region II, condition (17) is not satisfied and ( 14) has no positive roots.Hence,  * is stable only when  = 0.In region III, Theorem 3 (ii) holds and there exists a nonnegative integer  such that  * is locally asymptotically stable when and it is unstable when To simulate clear stability switches, Figure 1   On the basis of the above analysis, we present examples in Figures 2-5.The equilibrium point  * is locally asymptotically stable in Figure 2. When  passes through critical values,  * loses its stability, and bifurcating periodic solutions occur and are stable, as shown in Figures 3-4.In Figure 5, the bifurcating solutions disappear and chaos occurs when  passes through critical values.
The numerical results in Figures 2-5 show that oscillating solutions are strongly affected by .To generalize the dynamical behavior of system (3) influenced by parameter , we make numerical simulation using a wide range of values of parameter  to show the bifurcation diagram (see Figure 6(a)).Results for the nutrient density and phytoplankton density as a function of  are similar and are not shown here.From Figure 6(a), the order-2 periodic solution appears in the system (3) by increasing parameter .Keeping increasing , chaos occurs.Subsequently, an order-3 periodic solution bifurcates from chaos.Finally, chaos occurs again with increase of .To investigate the influence of parameter  on amplitudes of the oscillations observed when the positive equilibrium undergoes Hopf bifurcations, we simulate the stability switches of positive equilibrium in the component  (see Figure 6  1 2 = 50.4243,respectively.The black solid lines denote that the positive equilibrium component  * = 16.92598 is stable when  ∈ [0,  1 0 ) ∪ ( 2 0 ,  1 1 ) ∪ ( 2 1 ,  1 2 ).The blue and red curves denote the maximum and minimum  values, respectively.The black dashed lines indicate that when  ∈ ( 1 0 ,  2 0 ) ∪ ( 1  1 ,  2 1 ) ∪ ( 1 2 , +∞), the equilibrium loses its stability and oscillations occur.Therefore, the numerical simulations agree with the theoretical predictions.

Conclusion
In this paper, we have studied the dynamics of a nutrientplankton system with a time delay.We have investigated the stability of the positive equilibrium and the existence of Hopf bifurcation.By using center manifold theory and the normal form method, we determined the direction of Hopf bifurcation and the stability of bifurcating periodic solutions.
In detail, we have found that if some conditions are satisfied, the phenomenon of stability switches arises from system (3).When the time delay passes through some critical values, there exists a nonnegative integer  such that the positive equilibrium switches  times from stability to instability to stability and so on, and Hopf bifurcations occur at the positive equilibrium.The numerical simulations in Figure 1 indicate that the theoretical results are correct, and the number of stability switches changes with the value of the external source of nutrients .This result indicates that the external source of nutrients flowing into the system has important influence on the complex dynamics of system (3).
Numerical simulations also showed that as the time delay further increases, the periodic solutions disappear and chaos appears.All these results not only will help in further investigating the dynamics of pelagic ecosystem in theory but also are very useful to understand the complex phenomena really happening in marine ecosystem.

Figure 1 :
Figure 1: (a) Relationship between the conversion frequency  of stability switches for positive equilibrium  * and the external source of nutrients , where the blue spots represent the value of , the black line represents the boundary of the existence of  * , the green line represents the boundary of the appearance of , and the red line represents the disappearance of .(b) Interval values for stability switches and Hopf bifurcation about  * , where the green diamonds and the blue dots represent  1  (0 ≤  ≤ ) and  2  (0 ≤  ≤  − 1), respectively, and the honeydew and deep pink rectangles, respectively, represent the stable and unstable interval of  * .

Figure 1 (
Figure 1(b), it is easy to see that the stability switches of positive equilibrium  * occur with the change of time delay .On the basis of the above analysis, we present examples in Figures2-5.The equilibrium point  * is locally asymptotically stable in Figure2.When  passes through critical values,  * loses its stability, and bifurcating periodic solutions occur and are stable, as shown in Figures3-4.In Figure5, the bifurcating solutions disappear and chaos occurs when  passes through critical values.The numerical results in Figures2-5show that oscillating solutions are strongly affected by .To generalize the dynamical behavior of system (3) influenced by parameter , we make numerical simulation using a wide range of

Figure 6 :
Figure 6: (a) Bifurcation diagram of system (3) with  = 0.5 and initial point (0.6, 0.7, 16); (b) the stability switches of positive equilibrium in the component  with respect to time delay  for  = 0.2991 and  * = 16.92598,where the green point denotes the Hopf bifurcation point, the black solid line denotes the stable positive equilibrium, the black dashed line denotes the unstable positive equilibrium, and the blue and red curves denote the maximum and minimum values of zooplankton density , respectively.