This study considers the spatiotemporal dynamics of a reaction-diffusion phytoplankton-zooplankton system with a double Allee effect on prey under a homogeneous boundary condition. The qualitative properties are analyzed, including the local stability of all equilibria and the global asymptotic property of the unique positive equilibrium. We also discuss the Hopf bifurcation and the steady state bifurcation of the system. These results are expected to help understand the complexity of the Allee effect and the interaction between phytoplankton and zooplankton.

1. Introduction

The upper layer of the ocean contains large volumes of drifting plankton, which can be divided into phytoplankton and zooplankton. Phytoplankton is the autotrophic component of the plankton community, which is consumed by zooplankton, most of which are too small to observe individually with the naked eye. Zooplankton, which are heterotrophic organisms in oceans, are also mostly invisible to the naked eye. Therefore, it is difficult and expensive to quantify plankton directly. Plankton not only play an important role in the marine system because they are at the bottom level of the food chain that supports commercial fisheries, but also play important roles in the cycling of many chemical elements, such as carbon, which may affect climate change [1]. Furthermore, when plankton such as blue algae and dinoflagellates are present in large concentrations, the water appears to be discolored or murky, which is known as a red tide, and this can result in the death of marine and coastal species of fish, mammals, and other organisms [1]. Thus, analyzing the dynamics of plankton using mathematical models is beneficial for understanding the features of plankton populations, which have enormous economic and ecological value.

However, the mechanism that leads to the occurrence of red tides is still an unsolved issue. Many models and theories have been proposed by mathematicians and ecologists to explain this phenomenon, but a general and correct explanation still remains a distant goal [25].

The popular mathematical model called a Gause-type predator-prey model is used to consider the phytoplankton-zooplankton interaction in the following form: where and are the population densities or biomass of phytoplankton and zooplankton at time , respectively. describes the intrinsic per capita growth rate of the phytoplankton, which may be a logistic growth function, exponential growth function, or other functions, where is known as the environmental carrying capacity. is the per unit-predator consumption rate of prey, which is commonly called the functional response. Some conventional forms of functional response include Holling types I, II, III, and IV and Ivlev type [68]. is the conversion coefficient and represents the per capita predator mortality, which is assumed to have a linear form, although other forms are possible [9, 10]. The global dynamics of model (1) with a logistic growth rate have been studied during the last three decades based on theoretical analysis and numerical simulations, and many results have been reported [1115].

In recent years, the Allee effect has been the focus of increasing interest and it is recognized to be an important phenomenon in many fields of ecology and conservation biology by more and more people [1621]. The Allee effect is named after W.C. Allee [22] and it describes a positive correlation between the density or number of population and individual fitness of population [16]. Standard population models assume that the fitness of population increases as the population density or size declines [1115, 23, 24], whereas Allee effect states that when a population is below a critical density or size, the population cannot sustain itself and this leads to extinction. Thus, the Allee effect increases the likelihood of extinction [25]. Stephens et al. distinguished between a component Allee effect and a demographic Allee effect [17]. However, conservation biologists are usually more interested in the demographic Allee effect because it ultimately governs the probability of the extinction or recovery of populations with low abundances [16].

Very recent ecological research has shown that two or more Allee effects can act on a single population simultaneously, which is known as the multiple (double) Allee effect [26, 27].

There are many ways of describing the Allee effect [28], including the following differential equation: where describes the growth rate and is an auxiliary parameter where ,  . Indeed, it is considered that represents the size of a fertile population and is the nonfertile population, such as juvenile or oldest individuals [29]. In this case, is called the Allee threshold because when the population density or size is below this threshold, the population is destined for extinction. When , (2) describes a strong Allee effect [8, 30, 31]. In this case, the population growth rate decreases if the population size is below the threshold and the population goes to extinction [29]. In addition, (2) describes a weak Allee effect [29, 3133] for .

It is obvious that (2) is equal to González-Olivares et al. [29] state that (3) describes a double Allee effect, that is, once in the factor and the second time in the term [34, 35].

In a marine environment, the plankton populations tend to move in horizontal and vertical directions due to the strong water current. This movement is usually modeled by a reaction-diffusion equation. In this study, we consider the following reaction-diffusion model with constant diffusion coefficient as well as a strong Allee effect in different spatial locations within a fixed smooth bounded domain . We assume that the response function of the zooplankton follows the law of mass action [15]: where is the Allee threshold, are the diffusion coefficients of phytoplankton and zooplankton, respectively, is the Laplacian operator, and is the spatial habitat of two species, and we assume that the system is a close ecosystem and with a no-flux boundary condition.

This paper is structured as follows. In Section 2, we analyze the basic dynamics of (4) including estimates of the solution and the local and global stability of equilibria. In Section 3, we provide the analysis of the Hopf bifurcation and the steady state bifurcation. A brief discussion and summary are given in Section 4.

2. Main Results

2.1. Basic Dynamics

Suppose that , , , and is a bounded domain; then we obtain the following results.

Theorem 1. The system (4) has a unique solution and the solution is bounded. Furthermore, the solution of (4) satisfies if , .

Proof. Let then , for , which implies that (4) is a mixed quasi-monotone system [36]. We define , , where is the solution of the system where ,  . We find that Therefore, Similarly, This implies that and are the lower solution and upper solution to (4), respectively (Figure 2). Therefore, Theorem in [36] shows that system (4) has a unique solution that satisfies According to the strong maximum principle and the boundary condition, , for and .
From the first equation of (7), we can see that for and for . Thus, . To estimate , let , ; then (12)   (13) and with boundary conditions, we obtain From the proof above, we know that . Thus, for any , , when By integrating (15), , we have This implies that .
If , we can add the two equations in (4) and we have Since we know that the solution of the equation is for . The comparison argument implies that .

Theorem 2. If , then or as .

Proof. From the proof of Theorem 1, if , then and consequently as . This completes the proof.

From a biological viewpoint, this implies that if the initial population density is below the threshold , the phytoplankton become extinct so the zooplankton would become extinct.

2.2. Local and Global Stability of Equilibria

System (4) has four nonnegative steady state solutions: , , , and , where ,  .

The local stability of the steady state solutions can be analyzed as follows.

Theorem 3. (1) is locally asymptotically stable.
(2)  is unstable.
(3)  is locally asymptotically stable when and is unstable for .
(4)   is locally asymptotically stable for and is unstable for , where

Proof. The linearization of (4) at solution can be expressed as where According to Theorems and from [37], we know that if all the eigenvalues of the operator have negative real parts, then the solution is asymptotically stable; if there is at least one eigenvalue with a positive real part, then the solution is unstable; if some eigenvalues have zero real parts, then the stability cannot be determined using this method.
Let be the eigenvalues of on under a homogeneous Neumann boundary condition and , . Thus, it is known that is an eigenvalue of if and only if is the eigenvalue of the matrix for some .
Consider the characteristic equation where
(1) If , then This implies that is locally asymptotically stable.
(2) If , then .
For , one of the eigenvalues is , which implies that is an unstable point.
(3) If , then .
When , This implies that is locally asymptotically stable.
When , for , , which implies that has at least one eigenvalue with positive real part. This implies that is unstable.
(4) If , then , where and implies that . Let be the largest root of .
When , then and This implies that is locally asymptotically stable.
When , then and for , which implies that has at least one eigenvalue with a positive real part. This implies that is unstable.

Theorem 4. If , then the positive constant steady state is globally asymptotically stable.

Proof. Let us consider a Lyapunov function as where and will be determined next.
By taking the time derivative of , we have Due to the Neumann boundary condition, it can easily be derived that Further, If we choose arbitrarily and , then we can obtain Therefore, if , then and if . This completes the proof.

3. Bifurcation Analysis

In this section, we mainly analyze the stability of the steady state and take as the bifurcation parameter (or equivalently take as a parameter). In particular, we assume that all of the eigenvalues of are simple.

We know from the proof of Theorem 2 that the stability of is determined by the trace and determinant of . Let We refer to as the Hopf bifurcation curve and as the steady state bifurcation curve [38] (Figures 4 and 5).

First, is equal to . We can summarize the properties of as follows, which are easy to prove so the proof is omitted.

Lemma 5. Consider ; then exists such that the following hold.(1)If , then in and .(2)If , then such that in , in , and .(3) in , in , and .
Second, is equal to . Let Since according to the continuity of , there exists a , such that . Thus, we can summarize the properties of as follows.

Lemma 6. If , then is decreasing and is increasing in .

Proof. Differentiating with respect to , we have Therefore, . From the definition of , when , and .
However, , , and  imply that , and Since and , then . Therefore, , , and when . If , then for . Thus, . This implies that .

3.1. Hopf Bifurcation Analysis

In this section, we mainly analyze the properties of the Hopf bifurcation for system (4). According to [38], we know that a Hopf bifurcation point must satisfy the following conditions.(A)There exists such that and the unique pair of complex eigenvalues , , exist and are continuously differentiable in , with .

Theorem 7. If holds, then exists such that the system (4) undergoes a Hopf bifurcation at and a smooth curve of positive periodic orbits of (4) bifurcates from . The bifurcating periodic orbits from are spatially homogeneous and the Hopf bifurcation at is supercritical and backward if .

Proof. It can be verified that and for , which implies that the potential Hopf bifurcation point must be in . However, , and ,   for . Therefore is a Hopf bifurcation point. If holds, we know that decreases strictly in . For every , let be the solution of , so we have where is the largest eigenvalue for . Therefore, and for , . Geometrically, from we can determine that is a parabola of with in coordinate system. However, implies that , implies that or , and Lemma 6 implies that is decreasing and is increasing in if . Thus, the curves and have only one intersection point, which is noted as in . Then, if and if , . According to Theorem 2.1 in [13], a smooth curve of positive periodic orbits of (4) bifurcates from . According to Theorem 3.1 in [8], if holds.

3.2. Steady State Bifurcation Analysis

In this section, we consider the steady state bifurcation of system (4). The nonnegative steady state solutions of (4) satisfy the following system: Apparently (43) has spatially homogeneous solutions , , , . First, we discuss the nonnegative steady state solutions of (4). Recall the maximum principle [13].

Lemma 8. Let be a bounded Lipschitz domain in , and let . If is a weak solution of the inequalities and if there is a constant such that for , then a.e. in .

From Lemma 8, it can easily be derived that all nontrivial solutions of equation satisfy .

Theorem 9. The solutions of system (43) are in the form of either or satisfying

Proof. If exists such that , from the strong maximum principle and the boundary condition , we can derive in , so satisfies (45). Similarly, if exists such that , we can derive and consequently . If, for all , and hold, according to the above discussion and the strong maximum principle we have for all . After adding the two equations in (43), we have The maximum principle and Green formula imply that Therefore, we can show the nonexistence of positive steady state solutions when the diffusion coefficients are large.

Theorem 10. Let be a fixed constant. Then, another constant exists such that if and hold, (43) has no nonconstant positive solution.

Proof. Assume that is a positive solution of (43). For convenience, we denote Let then , . By multiplying by the first equation in (43) and then integrating on , we have Similarly, where is defined in Theorem 9.
Therefore, Using the Poincáre inequality, we can obtain Under the assumption , a sufficiently small exists such that . Let ; then we can have , . This completes the proof.

We assume that all of the eigenvalues of are simple and the corresponding eigenfunctions are denoted by . Reference [13] gives an example of with and . Let According to [13], we know that a steady state bifurcation point must satisfy the following conditions.(B) exists such that

Theorem 11. exists such that system (43) undergoes a steady state bifurcation at if holds.

Proof. Apparently (B) is not established for . It is known that is a degree 3 polynomial of and there are at most three for . In particular, if the parameters are selected such that   holds, then, for each , there exists a unique such that . Thus, there is at most one bifurcation point. However, . From the proof of Lemma 6 we know that for . Therefore, holds if .

4. Discussion

Reaction-diffusion phytoplankton-zooplankton models with Allee effects have been studied extensively in recent years. In this study, we rigorously considered a Gause-type predator-prey model with a double Allee effect on prey, which was formulated as (4). It is known that the predator-prey model with the most usual form of Allee effect has a unique limit cycle, but the existence of two limit cycles was proved by González-Olivares et al. [29] with a double Allee effect. Thus, the double Allee effect produces different results with different mathematical expressions.

The paper [15] found that the system without Allee effect was always stable and without fluctuations, but in this paper the results of the stability of the equilibrium and the bifurcation analysis based on a rigorous theoretical analysis show that this system has complex spatiotemporal dynamics: for , the phytoplankton is destined to become extinct and leads to the extinction of zooplankton; after considering the strong Allee effect in phytoplankton, extinction for both species is always a locally stable equilibrium. But for , which is the condition in which the interior equilibrium exists, the interior equilibrium is globally stable in some case and there always exist some other spatiotemporal patterns in other cases (Figure 3).

Overall, our results indicate that the impact of the Allee effect increases the spatiotemporal complexity of the system. The mathematical form which expresses the double Allee effect has a strong impact of the dynamics of system. Thus, we think it is important for ecologists to be aware of the difference of the selection on the forms of Allee effect.

The limitations of our study are that we only consider a simplest phytoplankton-zooplankton interaction and the special formalization to describe the Allee effect. What is more is that, compared with the ODE dynamics, the results shown here are still coarse. Therefore, further research is still needed to elaborate a general theory on the influence of this ecological phenomenon.

Conflict of Interests

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


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