Abstract and Applied Analysis

Volume 2014 (2014), Article ID 451757, 10 pages

http://dx.doi.org/10.1155/2014/451757

Research Article

Nonlinear Dynamics of a Nutrient-Plankton Model

Yapei Wang,1,2 Min Zhao,2,3 Chuanjun Dai,2,3 and Xinhong Pan1,2

¹School of Mathematics and Information Science, Wenzhou University, Wenzhou, Zhejiang 325035, China

²Zhejiang Provincial Key Laboratory for Water Environment and Marine Biological Resources Protection, Wenzhou University, Wenzhou, Zhejiang 325035, China

³School of Life and Environmental Science, Wenzhou University, Wenzhou, Zhejiang 325035, China

Received 11 October 2013; Revised 3 December 2013; Accepted 4 December 2013; Published 16 January 2014

Academic Editor: Malay Banerjee

Copyright © 2014 Yapei Wang 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.

We investigated a nonlinear model of the interaction between nutrients and plankton, which was addressed using a pair of reaction-advection-diffusion equations. Based on numerical analysis, we studied a model without diffusion and sinking terms, and we found that the phytoplankton density (a stable state) increased with the increase of nutrient density. We analyzed the model using a linear analysis technique and found that the sinking of phytoplankton could affect the system. If the sinking velocity exceeded a certain critical value, the stable state became unstable and the wavelength of phytoplankton increased with the increase of sinking velocity. Furthermore, band patterns were also produced by our model, which was affected by the diffusion and sinking of phytoplankton. Thus, the change in the diffusion and sinking of phytoplankton led to different spatial distributions of phytoplankton. All of these results are expected to be useful in the study of plankton dynamics in aquatic ecosystems.

1. Introduction

Plankton play an important role in the ecology of the ocean and climate because of their participation in the global carbon cycle at the base of the food chain [1]. In certain environmental conditions, lakes, reservoir, and marine waters may experience plankton or algal blooms [2, 3]. However, the local and global impacts of plankton blooms on water quality, carbon cycling, and climate may be damaging. If nutrient source is abundant, and some conditions are satisfied, blooms may become long-term events that affect ecosystems. Plankton blooms can change the types of species present at the base of the aquatic food web and affect human health. Thus, the study of plankton dynamics is currently of major interest.

In the past years, there were many researches on the model between nutrient and phytoplankton and zooplankton [4–6]. A larger number of researchers have attempted to model the relationship between nutrient and phytoplankton and zooplankton, to investigate the dynamics in plankton model. Truscott and Brindley [7] presented a model for the evolution of phytoplankton and zooplankton populations which resembles models for the behavior of excitable media. Luo [8] investigated phytoplankton-zooplankton dynamics in periodic environments, where eutrophication was considered. El Saadi and Bah [9] modeled phytoplankton aggregation using numerical treatment and explored the asymptotic behavior of the model. Banerjee and Venturino [10] studied a phytoplankton-toxic phytoplankton-zooplankton model and found that the toxic phytoplankton does not drive the zooplankton population towards extinction under a certain mechanism. The result is very important for study on plankton. These works make contributions for the study on plankton.

In recent years, many ecologists have paid increasing attention to spatial processes in a wide variety of practical contexts [11]. For example, theoretical community ecologists have explored ecosystems, including vegetation systems [12] and phytoplankton systems [13]. In particular, the modeling of plankton systems is becoming increasingly important because of their roles in carbon cycling and temperature control, particularly their major impacts on global climate change [14]. These modeling strategies are focused on two areas: (i) studies of large and complex systems, which are eventually used to fit field data or to forecast future changes; and (ii) studies of skeleton models for various mechanisms, which can provide insights or stimulate new experiments [14].

The present study belongs to the latter area. We propose a model on phytoplankton using a pair of reaction-advection-diffusion equations, which allow spatial phenomena, such as sinking, and turbulence to be described directly, thereby enabling spatial structures to be studied. It is known that the sinking and mixing of phytoplankton have pronounced effects on the tendency of different phytoplankton to increase. Experimental studies indicated that most fresh water diatoms and other phytoplankton sink in undisturbed water [15]. Theoretical results also demonstrated the importance of sinking, mixing, and diffusion [16, 17]. Theoretical models predicted that a process with reduced vertical mixing may induce oscillations and chaos in the phytoplankton of the deep chlorophyll maxima, which leads to differences in the timescale between the sinking flux of phytoplankton and the upward flux of nutrients [18]. A remarkable finding was the survival of a sinking phytoplankton population even when the diffusivity in the deep layers could not prevent population washout [19–21]. Mellard et al. demonstrated how externally imposed heterogeneity in the form of resource gradients and mixing interacted with internally generated heterogeneity in the form of competition, population dynamics, and movement to determine the spatial distribution of phytoplankton [22]. Ryabov et al. showed that the upper mixed layer was an important factor that had the potential to shape the spatial distribution and species composition of phytoplankton, but it also provided insights with general ecological importance [23]. van de Koppel et al. studied self-organized spatial patterning in an algae-mussel model, where regular spatial patterns were formed in young mussel beds on soft sediments in the Wadden Sea [24]. Self-organized spatial patterns are of considerable interest to theoretical biology [25–30] because of the basic paper by Turing [31] on the role of nonequilibrium reaction-diffusion prepatterns in biomorphogenesis. Furthermore, recent modeling studies of plankton support the self-organized spatial patterns, such as patchiness [32, 33] and bands [25, 34, 35].

The rest of this paper is organized as follows. In the next section, we present a model based on the theoretical ecology and partial differential equations, which is addressed using a pair of reaction-advection-diffusion equations. In Section 3, we analyze stable behavior of the nonspatial system firstly. What is more, the stable behavior of the spatial system is analyzed. And we obtain the condition under which the steady state becomes unstable. Finally, a series of simulations are given. Using simulation, we investigate the effect of critical factor on the system. In Section 4, discussion and conclusion are presented.

2. The Model

Natural ecosystems of plankton exhibit great variability in space and time. The growth of phytoplankton is dependent mainly on nutrients and light. After the mortality of phytoplankton, nutrients are returned to the system over short time scales with minimal losses [39] through microbial decomposition. In addition, biological factors such as higher predation and physical factors such as the sinking of phytoplankton into the water column also affect the ecosystem, which has been examined previously [40]. Turbulence also affects these systems [3, 19, 41]. Vertical mixing brings nutrients from the lower layers of the ocean into the mixed layer. Based on the previous analysis, the following general structure is obtained:

where is the nutrient density and is the phytoplankton population density.

Dugdale proposed the use of Michaelis-Menten enzyme kinetics to describe nutrient-phytoplankton interactions [42]. The Michaelis-Menten equations have the same form as the well-known Monod equations [43], which are used in the Droop equations, and they have formed the basis of a number of modeling studies that aimed to simulate phytoplankton blooms [44]. Thus, we employed Michaelis-Menten kinetics in terms of “uptake.” Furthermore, a Holling type II functional response has been used widely to describe zooplankton predation in various theoretical studies [45, 46]. It has also been reported that the Holling type II functional response shows good concordance with experimental data [33, 47, 48]. Hence, in the present paper, Holling II functional response is adopted to describe zooplankton grazing on phytoplankton. Therefore, a pair of specific models is defined as follows:

where a vertical water column is considered. Let indicate the depth in the water column; is the width in the water column. For vertical mixing, we assume that is a constant concentration, which includes the nutrient input flowing into the system and the nutrient from below the mixed layer, is the rate of exchange between the lower and upper layers, is the nutrient content of phytoplankton, denotes the maximum growth rate of phytoplankton, is the half-saturation constant for nutrients, is the half-saturation constant for phytoplankton, denotes the maximum predation rate of zooplankton on phytoplankton, is the mortality of phytoplankton, is the proportion of nutrients in dead phytoplankton that is recycled, is the sinking velocity of phytoplankton, and and are the diffusion rates of nutrients and phytoplankton, respectively, which are caused by mixing and turbulence. In addition, is the Laplacian operator. Table 1 provides the parameter values used and their units, which is obtained from published studies [19, 36–38].

Table 1: Parameter values used.

Symbol
Value
Unit

0.5
g·m⁻³

0.08
day⁻¹

0.02
dimensionless

0.5
day⁻¹

0.01
dimensionless

0.24
day⁻¹

2
g·m⁻²·day⁻¹

0.005
g·m⁻²

4
g·m⁻²

1.008
m·day⁻¹

1.038
m²·day⁻¹

1.038
m²·day⁻¹

Note: parameter value was estimated, parameter value was estimated based on de Angelis et al. [36], and the other parameter values were derived from previous studies [19, 36–38].

3. Results

3.1. Stable Behavior of the Nonspatial System

In the nonspatial system (i.e., system (2a), (2b) without spatial derivatives), according to and , vertical isocline and horizontal isocline can be obtained, respectively, as follows. Vertical isocline : . Horizontal isocline :. Obviously, the line, , is asymptote of vertical isocline , and the line, , is asymptote of horizontal isocline . In the following discussion, it is assumed that the condition always holds; otherwise phytoplankton become extinct eventually. For vertical isocline , , which is derivative. There are two roots in when the condition holds, root is . Then, it is obvious that the asymptote, , is on the left of line . Therefore, the vertical isocline is continuous when . And holds when ; holds when . From the horizontal isocline , if , then there is a positive equilibrium in the nonspatial system at least; if , then there is a positive equilibrium in the nonspatial system at least when the condition holds.

When the condition holds, there is no root in if . From , vertical isocline is monotone decreasing when , and holds when ; holds when . According to horizontal isocline , if , then there is a positive equilibrium in the nonspatial system at least when the condition holds; if , then there is no positive equilibrium in the nonspatial system when the condition holds.

When the condition holds, there are two roots in if . From , vertical isocline is monotone increasing, when . And holds, when ; holds, when . According to horizontal isocline , if , then there is a positive equilibrium state in the nonspatial system at least when the condition holds; if , then there is no positive equilibrium state in the nonspatial system when the condition holds.

It is noted that these conditions only confirm that there exists positive equilibrium state in the nonspatial system when these conditions are satisfied. It does not mean that there must be no positive equilibrium state in the nonspatial system when these conditions are not satisfied. In addition, it is not difficult to find that there is always a trivial steady state, , consisting of bare nutrients without phytoplankton in the nonspatial model. The Jacobian matrix of nonspatial system at the equilibrium is

It is obvious that the index of equilibrium is , when the condition holds, which is stable. In particular, when the conditions and hold, if or and , then there is no positive equilibrium in nonspatial system. Under these conditions, equilibrium is locally asymptotically stable. Furthermore, equilibrium is globally asymptotically stable in . Because first quadrant is a positive invariant set according to and , there is no limit cycle because there is no equilibrium in first quadrant.

Based on previous discussion, there exists positive equilibrium in the nonspatial model under some conditions, which is defined by . The Jacobian matrix of nonspatial model at the equilibrium is

From , it is easy to find that the equilibrium is unstable when , which is saddle. When , the index of equilibrium is when the condition holds, and it is locally asymptotically stable using Roth-Hurwitz criterion when the condition holds.

Although the expression of equilibrium can hardly be obtained, the stable behavior of the nonspatial system is determined when some parameters are given in Table 1. In the present paper, our interest is how some factors, such as the nutrient concentration and nutrient cycling effect , affect the system. Hence, the stable behavior of the nonspatial system is analyzed using the graph (see Figure 1). In Figure 1, when the nutrient concentration or nutrient cycling effect increases, there is always a trivial equilibrium consisting of bare nutrients without phytoplankton. In Figures 1(a) and 1(b), when the nutrient concentration is , there is only a trivial steady state in the nonspatial system. When the nutrient concentration , there were two other steady states: one is always unstable (green dashed, saddle), while the other is stable (red solid line, focus). When the nutrient concentration is , the focus disappears, and a node emerges (blue line, node). In Figures 1(c) and 1(d), a similar analysis is obtained, and the difference among Figures 1(a), 1(b), 1(c), and 1(d) is indicated by the grey zone. The trivial steady state and saddle coexist in the grey zone. In the following discussion, the nontrivial homogeneous steady state (focus or node) is defined by .

(a)

(b)

(c)

(d)

Figure 1: The existence and stability of an equilibrium versus the parameter (a, b) and parameter (c, d). The other parameters are given in Table 1.

3.2. Stable Behavior of the Spatial System

In this section, we consider the sensitivity of the system (2a), (2b) to change in the parameter values. A linear analysis technique is employed to focus on the parameters essential for the system behavior [49]. Our interest is how the nutrient concentration , nutrient cycling effect , sinking velocity , and the diffusion rate of phytoplankton affect the system. Symmetry breaking occurred when the stationary homogeneous solution, , is linearly unstable to small spatial perturbations in the presence of diffusion and advection, but that is linearly stable to perturbations in the absence of the diffusion and advection terms. To analyze the spatial system and determine how a small heterogeneous perturbation of the homogeneous steady state developed within a large time period, the following perturbation is considered [41]:

where is the perturbation growth rate, is the wave-number, and is an imaginary unit (). Substituting expression (5) into (2a), (2b) and neglecting all nonlinear terms in and , the following characteristic equation is obtained for the eigenvalues :

where the elements of the Jacobian determinant of the nonspatial system are taken at the stationary homogeneous solution , as follows:

By , it is not difficult to find that , , and , when is positive equilibrium. When , . The characteristic equation (8) can be described as

where and . In the previous analysis, the parameters , , and are allowed to vary, but the other parameters are fixed in Table 1. In Figure 2(a), the value of is given when the parameters and are changed. In addition, in Figure 2(b), in zone III, , so for ; in zone II, , but for , so ; in zone I, , and for . In zone I, to determine the sign of for different values of and , we need to analyze further because of the for .

(a)

(b)

Figure 2: (a) An illustration of versus and . The other parameters are given in Table 1. (b) An illustration of the sign of versus and ,  (m·day⁻¹).

From expression (8), we can obtain

where and .

To analyze the spatial system, the real and imaginary parts of the eigenvalues must be obtained, which are described as follows:

where . The solution is stable when the real parts of all eigenvalues are less than zero; that is, . The solution is unstable when one of the real parts with a finite wave number is greater than zero at least. The critical point is got when . However, the analytical expression for the critical point is difficult to be obtained. Indeed, we only need to consider the maximum value of . Thus, the critical condition can be obtained using , as follows:

In expression (11), the solution is unstable if the right-hand side of the equal sign is always less than zero. Otherwise, a necessary condition for expression (11) to hold is that . We consider the following case: the nutrient concentration is allowed to vary, but the values of other parameters are in Table 1. Then, the sinking velocity is a function related to the nutrient concentration and the wave number . By Figure 2(b), if the diffusion rate of phytoplankton, , is larger than 0.2 cm²·s⁻¹, then . Thus, the right-hand side of expression (11) is positive within , and the sinking velocity has a minimum at the point . Figure 3(a) confirms expression (11). In Figure 3(a), the neutral curve is convex with a unique minimum in the range .

(a)

(b)

(c)

Figure 3: (a) A typical neutral curve , defined using expression (11) for different values of . (b) Numerical calculation of the stability on space. (c) An illustration of the dispersion relation ( versus the wave number ), where  (m²·day⁻¹) and the other parameters are given in Table 1. Blue line:  (m·day⁻¹); red line:  (m·day⁻¹); green line:  (m·day⁻¹).

The effects of sinking velocity and nutrient concentration on the behavior of system (2a), (2b) are shown in Figure 3(b), which shows the transition from a no phytoplankton state through a banded phytoplankton state to a homogeneous phytoplankton state when the nutrient concentration increases and the sinking velocity is fixed. However, when the nutrient concentration is fixed, a homogeneous phytoplankton state becomes a banded phytoplankton state with the increase of sinking velocity , that is, a banded self-organized spatial pattern emerges because of the sinking velocity .

In zone I of Figure 3(b), there is only a trivial steady state, that is, no phytoplankton. In zone II, the steady state is stable, while the steady state becomes unstable in zone III because of the effect of sinking velocity of phytoplankton, which is further confirmed by Figure 3(c), where  (m·day⁻¹) and . It is obvious that the maximal real part of is larger than zero when the sinking velocity of phytoplankton is larger than the critical value of the sinking velocity of phytoplankton; that is, the instability of the steady state will occur. The imaginary value of is not equal to zero. In Figure 3(b), the red line represents the critical value of the sinking velocity of phytoplankton. The critical value of the sinking velocity of phytoplankton increases with the increase of nutrient concentration ; that is, when the nutrient concentration increases, the sinking velocity cannot affect the stability of the stable state.

3.3. Effects of the Parameters on the Wavelength of the Banded Pattern

Section 3.2 has described how a banded phytoplankton state emerges when the sinking velocity has reached a certain critical value. The determination of the wavelength of the banded pattern is a key issue. In particular, how do the parameters, such as the nutrient concentration , and the sinking velocity , affect the change of wavelength? The relationships among the wavelength, the nutrient concentration , and the sinking velocity are shown in Figure 4, which shows that the wavelength increases when the sinking velocity exceeds the critical value , but the wavelength decreases when the nutrient concentration increases.

Figure 4: An illustration of the variation in the pattern wavelength with and , where the symbol (■) represents the critical value of the sinking velocity ,  (m²·day⁻¹), and the other parameters are given in Table 1.

3.4. The Simulation

In the previous sections, we discussed the effects of parameters, including the nutrient concentration and the sinking velocity , on the system (2a), (2b). In this section, we discuss the numerical solution of the system (2a), (2b) in one-dimensional and two-dimensional spaces. In a one-dimensional space, a periodic boundary condition is employed, and system (2a), (2b) is solved on a rectangular spatial grid of 1 × 200 points. In two-dimensional space, system (2a), (2b) is studied in a horizontal (x, z)-plane with zero-flux boundary conditions (left and right) and periodic boundary conditions (top and bottom), which is solved on a rectangular spatial grid of 100 × 300 points. The initial conditions comprise a homogeneous state which is randomly perturbed. Furthermore, we assume that the diffusion rate of phytoplankton is larger or smaller than that of nutrients because of the viscosity and living of phytoplankton. Of course, it is also feasible that the diffusion rate of phytoplankton is equal to that of nutrients.

Firstly, the one-dimensional solution of system (2a), (2b) is shown in Figure 5. In Figure 5, we consider a vertical water column, where the depth of the water column is 120 m and the time is 600 days. We found that oscillations did occur; that is, the stable state became unstable because of spatial effects. Figure 6 shows the analysis of relationship between nutrients and phytoplankton further. In Figure 6(a), the relationship between the spatial distributions of nutrients and phytoplankton is given at the 600th day, which shows that nutrient concentration reaches the minimal value when the density of phytoplankton reaches the maximal value. The nutrient concentration affects the density of phytoplankton, and the density of phytoplankton increases with the increases of nutrient concentration. Thus, eutrophication may explain phytoplankton blooms. Furthermore, the effects of phytoplankton sinking on the relative maxima for nutrients and phytoplankton are shown in Figure 6(b). The relative maxima of phytoplankton increase with the increase of sinking velocity, whereas the relative maxima of nutrients decrease with the increase of sinking velocity. Therefore, the sinking flux has an important role in the increase of the density of phytoplankton.

(a)

(b)

Figure 5: Model simulations based on depth and time: (a) nutrient density (unit g·m⁻³); (b) phytoplankton density (unit g·m⁻³), where (m²·day⁻¹), (m·day⁻¹), and the other parameters are given in Table 1.

(a)

(b)

Figure 6: (a) Spatial distributions of the nutrient density and phytoplankton density on the 600th day, where the left vertical axis denotes the nutrient density, the right vertical axis denotes the phytoplankton density,  (m²·day⁻¹),  (m·day⁻¹), and the other parameters are given in Table 1. (b) Relationship between nutrient density, phytoplankton density, and the sinking velocity , where the right vertical axis denotes the maximal nutrient density, and the left vertical axis denotes the maximal phytoplankton density. The purple field shows that the equilibrium state is stable, where  (g·m⁻³) and  (m²·day⁻¹).

To further analyze the dynamic behavior of system (2a), (2b), we consider the solution of system (2a), (2b) in two-dimensional space. The band patterns are observed in the field, as shown in Figure 7. Figures 7(a), 7(b), and 7(c) show the patterns of phytoplankton at the 1000th day in the two-dimensional space. As discussed in Section 3.2, our numerical results confirm the predictions of the linear analysis that a band pattern of phytoplankton occurs if the nutrient concentration and the sinking velocity of phytoplankton satisfy some conditions. Figure 7(a) shows the emergence of parallel and crossed patterns, which indicate that band patterns with different speeds coexist in the system (2a), (2b) where the wavelength of patterns is different. By contrast, the patterns in Figure 7(b) are much more regular and almost parallel. In the real world, the sinking velocity of phytoplankton varies at different spatial points. Thus, to add more realism to the system, we forced the model to undergo periodic changes in the sinking velocity of phytoplankton, that is, , where denotes the width of water column. The results are shown in Figure 7(c), where the values of the parameters are the same as those used in Figures 7(a) and 7(b), except the sinking velocity . However, it is obvious that the patterns in Figure 7(c) are very different from those in Figures 7(a) and 7(b). Thus, the sinking flux of phytoplankton has an important role in the system.

(a)

(b)

(c)

Figure 7: Simulation in two-dimensional space. The figures show the density levels of the phytoplankton on the 1000th day, where the width is 100 m and the depth is 300 m. (a)  (m²·day⁻¹) and  (m·day⁻¹). (b)  (m²·day⁻¹) and  (m·day⁻¹). (c)  (m²·day⁻¹) and . The other parameters are given in Table 1.

4. Discussion and Conclusion

Banded patterns have been described in several resource-limited ecosystems around the world. In the real world, numerous population patterns have been observed, including banded vegetation, patches, and spiral waves, which can be regular or irregular. Physical factors may cause these types of pattern, such as wind, water flow, and turbulence. Internal factors in populations also force these patterns to occur.

In the present study, we used a nutrient-plankton model with both diffusion and advection to investigate the interaction between nutrient and plankton. Our model was simple because it was only an abstraction of real-world phenomena but the model reproduced many features of real-world phenomena. Our explanation focuses on a predator-prey interaction between phytoplankton and their nutrient source. In particular, how do the sinking of phytoplankton and the input of nutrients affect the interaction? Our analytical results showed that the homogeneous steady state became unstable because of the sinking of phytoplankton. The critical value of the sinking phytoplankton led to an instability in the homogeneous steady state, which depended on the input of nutrients. Our numerical results showed that the homogeneous steady state was unstable against small spatially heterogeneous perturbations.

Figure 1(b) shows that when the nutrient concentration increased beyond a critical value, the increase in the concentration of phytoplankton was stable; that is, the concentration of phytoplankton tended toward a certain stable state. Spatial effects did not influence the stable state when the sinking flux was below a critical value, as shown in Figure 3(b). Thus, an abundance of nutrient inputs flowed into the system, which led to the high-level reproduction of phytoplankton, which may trigger phytoplankton blooms.

Figures 5 and 7 show that oscillation could occur because of the sinking flux. In particular, Figure 5 shows that both spatial and temporal oscillations arose in the nutrients and the phytoplankton. The sinking of phytoplankton can also lead to the increase in the phytoplankton density and wavelength when the input of nutrients is fixed. It is possible that the phytoplankton sinks from the surface of water until it reaches a depth where the nutrient conditions are suitable for growth. Figure 6(a) shows that the relationship between phytoplankton and nutrients is mutually constrained. Thus, abundant nutrition leads to the mass propagation of phytoplankton, which consumes large amounts of nutrient, thereby depleting the nutrient levels. Thus, the sinking of phytoplankton and the input of nutrients can change the spatial distribution of phytoplankton under these conditions, which may promote the increase of phytoplankton density. In particular, eutrophication may promote phytoplankton blooms.

Conflict of Interests

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

Acknowledgments

This work was supported by the Key Program of Zhejiang Provincial Natural Science Foundation of China (Grant no. LZ12C03001), the National Natural Science Foundation of China (Grant no. 31170338), and the National Key Basic Research Program of China (973 Program, Grant no. 2012CB426510).

References

R. Reigada, R. M. Hillary, M. A. Bees, J. M. Sancho, and F. Sagués, “Plankton blooms induced by turbulent flows,” Proceedings of the Royal Society B, vol. 270, no. 1517, pp. 875–880, 2003.

A. Huppert, B. Blasius, R. Olinky, and L. Stone, “A model for seasonal phytoplankton blooms,” Journal of Theoretical Biology, vol. 236, no. 3, pp. 276–290, 2005.

M. Sandulescu, C. López, and U. Feudel, “Plankton blooms in vortices: the role of biological and hydrodynamic timescales,” Nonlinear Processes in Geophysics, vol. 14, no. 4, pp. 443–454, 2007.

C. J. Dai, M. Zhao, and L. S. Chen, “Bifurcations and periodic solutions for an algae-fish semicontinuous system,” Abstract and Applied Analysis, vol. 2013, Article ID 619721, 11 pages, 2013.

S. Abbas, M. Banerjee, and N. Hungerbühler, “Existence, uniqueness and stability analysis of allelopathic stimulatory phytoplankton model,” Journal of Mathematical Analysis and Applications, vol. 367, no. 1, pp. 249–259, 2010.

J. F. Zhang and D. Zhang, “Hopf-Pitchfork bifurcation in a phytoplankton-zooplankton model with delays,” Abstract and Applied Analysis, vol. 2013, Article ID 340174, 5 pages, 2013.

J. E. Truscott and J. Brindley, “Ocean plankton populations as excitable media,” Bulletin of Mathematical Biology, vol. 56, no. 5, pp. 981–998, 1994.

J. H. Luo, “Phytoplankton-zooplankton dynamics in periodic environments taking into account eutrophication,” Mathematical Biosciences, vol. 245, pp. 126–136, 2013.

N. El Saadi and A. Bah, “Numerical treatment of a nonlocal model for phytoplankton aggregation,” Applied Mathematics and Computation, vol. 218, no. 17, pp. 8279–8287, 2012.

M. Banerjee and E. Venturino, “A phytoplankton-toxic phytoplankton-zooplankton model,” Ecological Complexity, vol. 8, no. 3, pp. 239–248, 2011.

E. E. Holmes, M. A. Lewis, J. E. Banks, and R. R. Veit, “Partial differential equations in ecology: spatial interactions and population dynamics,” Ecology, vol. 75, no. 1, pp. 17–29, 1994.

C. A. Klausmeier, “Regular and irregular patterns in semiarid vegetation,” Science, vol. 284, no. 5421, pp. 1826–1828, 1999.

A. D. Barton, S. Dutkiewicz, G. Flierl, J. Bragg, and M. J. Follows, “Patterns of diversity in marine phytoplankton,” Science, vol. 327, no. 5972, pp. 1509–1511, 2010.

H. Malchow, “Spatio-temporal pattern formation in nonlinear non-equilibrium plankton dynamics,” Proceedings of the Royal Society B, vol. 251, no. 1331, pp. 103–109, 1993.

D. A. Bella, “Simulating the effect of sinking and vertical mixing on algal population dynamics,” Journal of the Water Pollution Control Federation, vol. 42, no. 5, pp. 140–152, 1970.

O. Kerimoglu, D. Straile, and F. Peeters, “Role of phytoplankton cell size on the competition for nutrients and light in incompletely mixed systems,” Journal of Theoretical Biology, vol. 300, pp. 330–343, 2012.

R. K. Upadhyay, N. Kumari, and V. Rai, “Wave of chaos in a diffusive system: generating realistic patterns of patchiness in plankton-fish dynamics,” Chaos, Solitons and Fractals, vol. 40, no. 1, pp. 262–276, 2009.

J. Huisman, N. N. Pham Thi, D. M. Karl, and B. Sommeijer, “Reduced mixing generates oscillations and chaos in the oceanic deep chlorophyll maximum,” Nature, vol. 439, no. 7074, pp. 322–325, 2006.

N. Shigesada and A. Okubo, “Effects of competition and shading in planktonic communities,” Journal of Mathematical Biology, vol. 12, no. 3, pp. 311–326, 1981.

D. C. Speirs and W. S. C. Gurney, “Population persistence in rivers and estuaries,” Ecology, vol. 82, no. 5, pp. 1219–1237, 2001.

A. V. Straube and A. Pikovsky, “Mixing-induced global modes in open active flow,” Physical Review Letters, vol. 99, no. 18, Article ID 184503, 2007.

J. P. Mellard, K. Yoshiyama, E. Litchman, and C. A. Klausmeier, “The vertical distribution of phytoplankton in stratified water columns,” Journal of Theoretical Biology, vol. 269, no. 1, pp. 16–30, 2011.

A. B. Ryabov, L. Rudolf, and B. Blasius, “Vertical distribution and composition of phytoplankton under the influence of an upper mixed layer,” Journal of Theoretical Biology, vol. 263, no. 1, pp. 120–133, 2010.

J. van de Koppel, M. Rietkerk, N. Dankers, and P. M. J. Herman, “Scale-dependent feedback and regular spatial patterns in young mussel beds,” The American Naturalist, vol. 165, no. 3, pp. 66–77, 2005.

W. M. Wang, Z. G. Guo, R. K. Upadhyay, and Y. Z. Lin, “Pattern formation in a cross-diffusive Holling-Tanner model,” Discrete Dynamics in Nature and Society, vol. 2012, Article ID 828219, 12 pages, 2012.

M. P. Hassell, H. N. Comins, and R. M. May, “Spatial structure and chaos in insect population dynamics,” Nature, vol. 353, no. 6341, pp. 255–258, 1991.

R. K. Upadhyay, N. K. Thakur, and V. Rai, “Diffusion-driven instabilities and spatio-temporal patterns in an aquatic predator-prey system with Beddington-Deangelis type functional response,” International Journal of Bifurcation and Chaos, vol. 21, no. 3, pp. 663–684, 2011.

J. von Hardenberg, E. Meron, M. Shachak, and Y. Zarmi, “Diversity of vegetation patterns and desertification,” Physical Review Letters, vol. 87, no. 19, Article ID 198101, 4 pages, 2001.

M. Rietkerk, S. C. Dekker, P. C. De Ruiter, and J. Van De Koppel, “Self-organized patchiness and catastrophic shifts in ecosystems,” Science, vol. 305, no. 5692, pp. 1926–1929, 2004.

A. I. Borthagaray, M. A. Fuentes, and P. A. Marquet, “Vegetation pattern formation in a fog-dependent ecosystem,” Journal of Theoretical Biology, vol. 265, no. 1, pp. 18–26, 2010.

A. M. Turing, “The chemical basis of morphogenesis,” Philosophical Transactions of the Royal Society B, vol. 237, pp. 37–72, 1952.

R. Reigada, R. M. Hillary, M. A. Bees, J. M. Sancho, and F. Sagués, “Plankton blooms induced by turbulent flows,” Proceedings of the Royal Society B, vol. 270, no. 1517, pp. 875–880, 2003.

H. Serizawa, T. Amemiya, and K. Itoh, “Patchiness in a minimal nutrient—phytoplankton model,” Journal of Biosciences, vol. 33, no. 3, pp. 391–403, 2008.

R.-H. Wang, Q.-X. Liu, G.-Q. Sun, Z. Jin, and J. Van De Koppel, “Nonlinear dynamic and pattern bifurcations in a model for spatial patterns in young mussel beds,” Journal of the Royal Society Interface, vol. 6, no. 37, pp. 705–718, 2009.

Q. X. Liu, E. J. Weerman, P. M. J. Herman, et al., “Alternative mechanisms alter the emergent properties of self-organization in mussel beds,” Proceedings of the Royal Society B, vol. 279, no. 1739, pp. 2744–2753, 2012.

D. L. de Angelis, S. M. Bartell, and A. L. Brenkert, “Effects of nutrient recycling and food-chain length on resilience,” American Naturalist, vol. 134, pp. 778–805, 1989.

A. M. Edwards, “Adding detritus to a nutrient-phytoplankton-zooplankton model: a dynamical-systems approach,” Journal of Plankton Research, vol. 23, no. 4, pp. 389–413, 2001.

G. L. Bowie and W. B. Mills, Rates, Constants, and Kinetics Formulations in Surface Water Quality Modeling, US Envitonmental Protection Agency, Athens, Greece, 2nd edition, 1985.

M. Ramin, G. Perhar, Y. Shimoda, and G. B. Arhonditsis, “Examination of the effects of nutrient regeneration mechanisms on plankton dynamics using aquatic biogeochemical modeling,” Ecological Modelling, vol. 240, pp. 139–155, 2012.

A. Huppert, B. Blasius, and L. Stone, “A model of phytoplankton blooms,” American Naturalist, vol. 159, no. 2, pp. 156–171, 2002.

M. Lévy, D. Iovino, L. Resplandy et al., “Large-scale impacts of submesoscale dynamics on phytoplankton: local and remote effects,” Ocean Modelling, vol. 43-44, pp. 77–93, 2012.

R. C. Dugdale, “Nutrient limitation in the sea: dynamics, identification, and significance,” Limnology and Oceanography, vol. 12, pp. 685–695, 1967.

D. L. de Angelis, Dynamics of Nutrient Cycling and Food Webs, vol. 9, Springer, 1992.

W. J. O’Brien, “The dynamics of nutrient limitation of phytoplankton algae: a model reconsidered,” Ecology, vol. 55, no. 1, pp. 135–141, 1974.

C. S. Holling, “Resilience and stability of ecological systems,” Annual Review of Ecology, Evolution, and Systematics, vol. 4, pp. 1–23, 1973.

M. Scheffer, “Fish and nutrients interplay determines algal biomass: a minimal model,” Oikos, vol. 62, no. 3, pp. 271–282, 1991.

W. Gentleman, A. Leising, B. Frost, S. Strom, and J. Murray, “Functional responses for zooplankton feeding on multiple resources: a review of assumptions and biological dynamics,” Deep-Sea Research Part II, vol. 50, no. 22–26, pp. 2847–2875, 2003.

M. M. Mullin, E. F. Stewart, and F. J. Fuglister, “Ingestion by planktonic grazers as a function of food concentration,” Limnology and Oceanography, vol. 20, no. 2, pp. 259–262, 1975.

J. D. Murray, “Mathematical biology,” in Biomathematics, vol. 19, Springer, Berlin, Germany, 2nd edition, 1993.