Dynamics Analysis of an Aquatic Ecological Model with Temperature Effect ()
1. Introduction
As we all know, the temperature is an important environmental factor that determines the growth of the algal population, is one of the important factors that affect the growth of algal cells, the composition and content of biological macromolecules in cells, and is also a key ecological factor that affects the growth, reproduction and population succession of aquatic plants [1] [2]. Furthermore, suitable temperature is a necessary condition for algal bloom outbreaks and also an important environmental factor for the replacement of the dominant algal population [3] [4]. Thus, it is important to study the effect of temperature on algal growth and provide the theoretical basis for the prevention and treatment of water eutrophication.
Algae is a ubiquitous photosynthetic organism; the algal growth rate can increase with increasing temperature up to a certain limit; this is because that temperature strongly influences the cellular chemical composition and uptake of nutrients and also plays a significant role in algal growth [5]. The paper [6] pointed out that the ability to model algal productivity under transient conditions of temperature was critical for assessing the profitability and sustainability of full-scale algae cultivation outdoors. The paper [7] investigated the effect of temperature and irradiance on the growth and reproduction of the green macroalga and gave that a suitable temperature range over 21˚C - 29˚C was more favorable for growth and reproduction. The paper [8] proposed that Microcystis aeruginosa had a wide range of adaptation to temperature, and the optimal growth temperature was 25˚C - 30˚C. The paper [9] showed that the temperature condition that was most conducive to algal growth was 25˚C, and the optimal condition for algal toxins release was 28˚C. The paper [10] inquired into the effects of temperature on the germination of micro-propagules via laboratory experiments and indicated that sea temperature played a significant role in the germination of green algae. The paper [11] sought to elucidate the effects of temperature on algal growth rates, biomass accumulation, fatty acid production and composition and pointed out that temperature significantly impacted the overall productivity of algal biofuel systems by influencing species growth rates and fatty acid production. The paper [12] studied the effects of different temperatures and illumination time on algal growth and obtained that the value of the algal growth rate constant was reduced to 0.812d−1 by lowering the water temperature to 16˚C. In conclusion, the algal population has strong adaptability to the environment, and temperature is one of the important ecological factors affecting algal growth.
In an aquatic ecosystem, the ecological effects between algae and biological manipulation predators are mutual, mainly including the harm of algal blooms to fish and the feeding and regulation of algae by biological manipulation predators [13]. The use of biological manipulation predators to control the excessive algal growth in eutrophic lakes was proposed and gradually attracted attention after the eutrophication of lakes became more and more common, and certain results in the control of cyanobacteria blooms have been achieved in some lakes at home and abroad [14] [15]. The paper [16] pointed out that silver carp completely eliminated cyanobacteria Microcystis by size and biovolume reduction. The paper [17] deemed that silver and bighead carps were just suitable for controlling cyanobacteria bloom not total algae biomass, and the application of the fish for algal control in water supply reservoirs should be cautious. The paper [18] believed that it was feasible to use silver carp and bighead carp to control Microcystis in eutrophic water. Furthermore, using ecological models to explore the dynamic relationship between algae and biological manipulation predators has also developed rapidly. The paper [19] showed that Microcystis aeruginosa aggregation could effectively control the dynamic feeding behavior of filter-feeding fish and provide shelter from predators. The paper [20] obtained that the filter-feeding fish population could be a crucial factor in controlling the proliferation of the algae population based on an algae-fish model. Therefore, fully understanding the ecological function relationship between algae and biological manipulation predator is a premise for better implementing algae control strategy with biological manipulation predators and ensuring algae control effect.
In this paper, firstly, the temperature variable is introduced to build an aquatic ecological model which can describe the dynamic relationship between algae and biological manipulation predators. Secondly, mathematical theory works are implemented to obtain some critical threshold conditions by investigating the evolution process of some specific dynamic properties. Finally, the numerical simulation works are carried out to show the evolution process of dynamic properties and the influence mechanism of temperature variables. Generally speaking, the main purpose of this paper is to use an aquatic ecological model with temperature effect to explore the coexistence mode between algae and biological manipulation predators and reveal how temperature variable affect their dynamic relationships.
2. Aquatic Ecological Model
At present, the temperature has been regarded as a key factor for monitoring and predicting algal blooms. Furthermore, the reproduction rate of algae has an important relationship with temperature, and the appropriate temperature is a necessary condition for the rapid propagation and growth of the algae population. In order to better investigate the effect of temperature on the dynamic relationship between algae and biological manipulation predator, we will propose an aquatic ecological model with temperature effect, which can be described as follow:
(2.1)
where
and
are density of algae population and biological manipulation predator (silver carp, bighead carp and anodonta woodiana elliptica) respectively.
is maximum growth rate of algae population at reference temperature
(it is considered to be 25˚C).
is the maximum environmental capacity for algae population. R represents actual temperature during the test run, which can be either a constant or a function of time t.
is the lowest temperature value when the algal growth is 0, and it is generally considered to be 15˚C.
is a respiratory rate,
is a nonpredation-induced mortality of algae population, s is a sedimentation rate,
is intrinsic growth rate of biological manipulation predator,
is mortality rate of biological manipulation predator, e is energy conversion rate, a is capture rate of biological manipulation predator, b is a semi saturation coefficient, and
is assimilated food of catabolic loss during predation period.
For simplicity, we will replace the model (2.1) with the following variables:
Then the model (2.1) can be rewritten as
(2.2)
For model (2.2), we first discuss the existence and stability of all possible equilibrium points and explore the existence of a limit cycle with some key constraints. Then we give some critical conditions to demonstrate the occurrence of transcritical bifurcation and Hopf bifurcation. Finally, some numerical simulation results are implemented to verify the validity and feasibility of theoretical results. At the same time, through some numerical simulation results, the influence of temperature on dynamic relationship between algae and biological manipulation predators will be explored, and then ecological evolution significance represented by bifurcation dynamic evolution behavior is also given.
3. Equilibrium Points and Their Stability
In this section, some preliminary results shall be presented, including the existence and stability of all possible equilibrium points of the model (2.2).
We will consider the following equation to explore all possible equilibrium points of the model (2.2),
(3.1)
It is easy to acquire that the model (2.2) always has trivial equilibrium point
, and a biological manipulation predator extinction equilibrium point
if
. In view of the biological significance and the characteristics of the model (2.2), the interior equilibrium points are conditional.
In order to explore the existence and stability of all possible equilibrium points of the model (2.2), we will gradually give the following Theorems 1 - 6.
Theorem 1. The model (2.2) has a unique interior equilibrium point
if and only if the following conditions hold,
Proof: From the previous analysis, we know that
has two distinct real roots, that is
or
, thus there are
and
. If we take
into
, there are
or
. Similarly, we take
into
, we can obtain
and
. According to the biological significance of the internal equilibrium point,
and
must be established, hence the existence equivalence condition is
. This ends the proof.
Theorem 2. The model (2.2) has a boundary equilibrium point
and it is always unstable.
1) If
,
is an unstable node; 2) if
, then
is a saddle.
Proof: The Jacobian matrix of the model (2.2) evaluated at
is
The eigenvalues of
are
,
. Based on the biological significance of the model (2.2), the value m must be less than 1. so
. Thereby, if
, then
, and
is an unstable node; if
, then
, and
is a saddle. This completes the proof.
Theorem 3. The model (2.2) always has a boundary equilibrium point
1) if
,
is a stable node; 2) if
, then
is a saddle.
Proof: The Jacobian matrix of the model (2.2) evaluated at
is
The eigenvalues of
are
,
. Since
, so when
, there is
, then
is a saddle. When
, there is
, then
is a stable node. This completes the proof.
Theorem 4. Under the condition
, the boundary equilibrium point
is globally asymptotically stable, where
.
Proof: This theorem can be derived by flow analysis. When
, the model (2.2) has no internal equilibrium point, and has only two boundary equilibrium points
and
. It is easy to know that
is always unstable and is a saddle, its unstable manifold is x-axis. Thus, we can divide the positive quadrant into the following three regions,
Here
is derived from the first equation of the Equations (3.1). Considering the biological significance of the model (2.2), we need only discuss the region
and
. In fact, solutions that start in region
will eventually enter region
by crossing the algae population nullcline vertically (downwards) in finite time, for the reason that there is no internal equilibrium point in
and
and
. Once in region
, it is clear to see that solutions are trapped, because
and
, so the solutions in
cannot cross the algae population nullcline, hence we must have
as
. This ends the proof.
Theorem 5. Under the condition that the theorem 1 exists, the equilibrium point
is locally asymptotically stable if
, and is unstable if
.
Proof: The Jacobian matrix of the model (2.2) evaluated at
is given by
The determinant and the trace of matrix
are given by
It is easy to check that
is always positive. Hence, if
, then the equilibrium point
is locally asymptotically stable; if
, then the equilibrium point
is unstable. If we can substitute the expression of
into the trace of the Jacobian matrix, we cannot directly derive the sign of
because of algebraic complexity of the expression; therefore, we will calculate the values of
by using numerical simulation in Section 5.
Theorem 6. If
is satisfied, then the internal equilibrium point
is globally asymptotically stable.
Proof: To prove the global asymptotic stability of the internal equilibrium point
, the following Lyapunov function was constructed:
Based on the basic properties of the function, it can be concluded that
is continuous for all
and
, and is computationally available:
Therefore, the internal equilibrium point
is only extreme value of the
function in the positive quadrant,
For all the
and
, we can get
Then, we can obtain
Furthermore, we know that the internal equilibrium point
satisfies the following equation:
Combined with the above two formulas for the simplification, we can obtain
By choosing
if
, we can obtain
According to the above formula, it is obvious that the internal equilibrium point
makes
and for all other
has
, so
satisfies the global stability theorem of the Lyapunov function. Thus the internal equilibrium point
is globally asymptotically stable when
is feasible and the implicit condition
is satisfied. This ends the proof.
Theorem 7. When the following two conditions are true, there exists at least one limit cycle of the model (2.2).
1)
; 2)
.
Proof: From the conclusion of the Theorem 5, we know that the internal equilibrium point
may be an unstable focus. Now we will prove Theorem 7 by constructing an invariant region
, which consists of the following line
, and
axis,
Through the derivation calculation, we can get
It is easy to verify that
for
and
. Furthermore, if
when M is a sufficiently large positive number, then
for
, which implies that as long as
holds, there must be a large enough positive number M to make
hold. Thus, the model (2.2) has at least one limit cycle by Poincare-Bendixson Theorem [21] [22]. This ends the proof.
4. Bifurcation Analysis
It is well known that the evolution process of bifurcation dynamics has important biological significance in the process of population dynamics. Therefore, we will explore some bifurcation dynamics behaviors of the model (2.2) and give some threshold conditions for specific bifurcation dynamics of the model (2.2).
4.1. Transcritical Bifurcation
Here we will prove that the model (2.2) undergoes a transcritical bifurcation at
. We recall from Theorem 3, the boundary equilibrium
loses its stability at
and one of the eigenvalues of
is zero; it implies that the equilibrium point
becomes non-hyperbolic, so it is possible to undergo a transcritical bifurcation around
.
Theorem 8. The model (2.2) undergoes a transcritical bifurcation when
, where
.
Proof: We use Sotomayor’s theorem [23] [24] to prove that the model (2.2) undergoes a transcritical bifurcation with a bifurcation parameter p. If
, then
, which means that the Jacobian matrix
has a zero eigenvalue. Now, let V and W be the eigenvectors corresponding to the zero eigenvalue of
and
respectively. After a simple calculation they can be given by
Therefore,
Due to
Furthermore, we can obtain
Thus, based on Sotomayor’s theorem we can deduce that the model (2.2) undergoes a transcritical bifurcation as the parameter p passes through a critical threshold
. This completes the proof.
4.2. Hopf Bifurcation
From the analysis of Theorem 5, we know that if
,
is locally asymptotically stable; if
, then
is unstable. Hence, we can easily conclude that the interior equilibrium point
may lose its stability by Hopf bifurcation under sufficient conditions. Considering p as a Hopf bifurcation control parameter, the Hopf bifurcation threshold
can be detected by using
. When the value of the parameter p passes from one side of
to the other side, the stability property of
changes and the periodic orbits can be generated. Thus, we can yield the following theorem.
Theorem 9. Based on Theorem 1, the model (2.2) undergoes a Hopf bifurcation around
when
.
Proof: In order to ensure the existence of Hopf bifurcation, we need to verify the transversality condition; then we have
and
This can guarantee the existence of Hopf bifurcation around
.
Next we will compute the first Lyapunov number to discuss the stability of limit cycle. Firstly, we will transform
to the origin by the transformation
, then we can get
where
and
and
are power series with terms
.
Therefore, the first Lyapunov number
is given by
where
If
, the limit cycle is stable; if
, the limit cycle is unstable. However, the expression for Lyapunov number
is rather cumbersome; we cannot directly judge the sign it, so we will give some numerical simulation results in Section 5.
Based on the mathematical theory, the existence and stability threshold conditions of all possible equilibrium points are given, and the critical conditions for inducing specific bifurcation dynamics of the model (2.2) are analyzed, which can provide a theoretical basis for subsequent numerical simulation work. Furthermore, it should also be emphasized that the key parameter p can seriously affect dynamic evolution characteristics of the model (2.2).
5. Numerical Simulations and Results
Now, we will investigate dynamic properties of the internal equilibrium point
and explore the effect of some parameters on dynamic relationship between
and
. From Theorem 1, we can see that the relation expression of x and y is
and
, thus some numerical simulations are given in Figure 1 with
,
,
,
,
Figure 1. (a) Dynamic relationship between population density x, y and parameter p value; (b) Dynamic relationship between population density x and parameter p value.
. It is easy to find from Figure 1(b) that the density of algae
increases with the increase of parameter p-value, and the growth is slow in the early stage and extremely fast in the later stage, which implies that when the value of parameter p exceeds a certain critical threshold, the density of algae
will become larger and larger. Furthermore, it is obvious to know from Figure 1(a) that when the value of parameter p exceeds a certain critical threshold and is determined, the density of biological manipulation predator
is a concave function of the density of algae
, and
can get a maximum value when the value of
is
, which must satisfy
, this to say
. Therefore, it must be emphasized that the dynamic properties of the internal equilibrium point
and the dynamic relationship between population
and
mainly depends on key parameters m and p.
In order to verify the validity and feasibility of Theorem 5 and 7, the stability of the internal equilibrium point
and the existence of the limit cycle are numerically simulated with
and
. It is easy to find from Figure 2(a) that the internal equilibrium point
is stable, which means that algae and biological manipulation predators can form a constant steady-state coexistence mode. Furthermore, it is obvious to know from Figure 2(b) that model (2.2) has a limit cycle, which indicates that algae and biological manipulation
Figure 2. (a) The stability of the internal equilibrium point
with
,
,
,
,
,
; (b) The existence of limit cycle with
,
,
,
,
,
.
predators can form a periodic oscillation coexistence mode. In a word, the value of parameter p has an important influence on coexistence mode of algae and biological manipulation predators.
In order to better understand how the value of parameter p affects the dynamic behavior evolution of the model (2.2), we give a bifurcation diagram of the model (2.2) in Figure 3. It can be seen clearly from Figure 3 and Figure 4 that if the value of parameter p is larger than
, the model (2.2) has only two boundary equilibrium points
and
,
is stable and
is unstable, which implies that biological manipulation predator will eventually approach extinction and algae will eventually approach the maximum biomass state, that is to say, biological manipulation predator and algae cannot form a final coexistence mode. If the value of parameter p is between
, the model (2.2) goes through a transcritical bifurcation, which will induce the model (2.2) to have an internal equilibrium point, and this internal equilibrium point is asymptotically stable, detailed numerical simulation results are shown in Figure 4(b) and Figure 4(a). Therefore, it can be said that the transcritical bifurcation induces the formation of a constant steady-state coexistence mode between biological manipulation predator and algae. If the value of parameter p gradually decreases and is less than the critical threshold
Figure 3. Bifurcation diagram of the model (2.2) by taking p as bifurcation parameter with
,
,
,
,
, Bifurcation diagram with respect to parameter
,
, where the red circle is a hollow dot, that means the boundary equilibrium point
is not exist, the pink and green lines implicit the critical values for the Hopf bifurcation and transcritical bifurcation respectively. In addition, the blue circle and blue asterisk separately denote the Hopf (PH) and transcritical (TC) bifurcation point.
value
, the internal equilibrium point
loses stability and a stable limit cycle appears. In other words, the model (2.2) undergoes a Hopf bifurcation; the numerical dynamic evolution process is shown in Figure 5. Therefore, it is worth pointing out that the Hopf bifurcation can produce a periodic oscillation coexistence mode between biological manipulation predator and algae. Thus, the numerical simulation results not only prove the validity and feasibility of the theoretical derivation, but also directly show that the value of key parameter p seriously affects bifurcation dynamics evolution characteristics of the model (2.2).
Figure 4. The phase portraits of the model (2.2),
has different dynamics with the value of parameter p varying. (a) the interior equilibrium point
exists and it is a stable node, boundary equilibrium
is a saddle, and
is an unstable node when
. (b) the interior equilibrium point
coincides with boundary equilibrium
, which is a saddle-node and the parabolic sector is on the upper half plane, and
is a saddle when
. (c)
is a stable node,
is a saddle, and there is no interior equilibrium point when
. The green curves represent the stable or unstable orbits, and the red points are some equilibrium points.
To explore in detail how temperature affects the coexistence mode of biological manipulation predators and algae, and analyze the advantages of the model (2.1), some numerical simulation results are shown in Figures 6-8. As we all know, in the laboratory algal-predator culture test, in order to better maintain the growth of algal population, we usually conduct the culture test at a constant temperature degrees. However, in subtropical reservoirs, the temperature has a periodic change in behavior with time, which is not a constant. Therefore, we will conduct numerical simulation with the state of laboratory culture temperature
,
and the state of field natural temperature
with
,
,
,
,
,
,
,
,
,
,
,
and
. It is obvious to find from Figure 6 and Figure 7 that algae and biological manipulation predators can
Figure 6. Time series diagram of the model (2.1) with
.
Figure 7. Time series diagram of the model (2.1) with
.
Figure 8. Time series diagram of the model (2.1) with
.
form a constant steady-state coexistence mode under the constant temperature assumption, however, algae and biological manipulation predator can form a periodic oscillation coexistence mode under the periodic change temperature assumption, In addition, for algae and biological manipulation predator, the center of the periodic oscillation amplitude is approximately a constant steady state value, which shows that the temperature expression does not affect the growth average biomass of algae and biological manipulation predator, but it will seriously affect their growth dynamics. Moreover, it is easy to see from Figure 8 that algae and biological manipulation predators have a periodic oscillation coexistence mode as the value of temperature parameter R increases to 32, and the biomass of biological manipulation predators has been greatly increased. These results not only show that the model (2.1) has experienced a Hopf bifurcation as the value of parameter R increases, but also indicate that the increase in temperature is conducive to rapid growth of biological manipulation predators.
Based on the numerical simulation results, can clearly indicate that the results of theoretical derivation are effective and feasible. Furthermore, it should also be emphasized that temperature not only affects the dynamic evolution characteristics of the model (2.1), but also affects the biomass level of biological manipulation predator. Moreover, the model (2.2) has specific bifurcation dynamic behaviors (transcritical bifurcation and Hopf bifurcation) under the influence of the value of key parameter p; these two bifurcation dynamics behaviors lead to a constant steady-state coexistence mode and a periodic oscillation coexistence mode of algae and biological manipulation predator respectively.
6. Conclusions
Under the conceptual framework of biological control of cyanobacteria in eutrophic lakes and reservoirs, based on the fact that temperature is an extremely important factor in determining ecology, which has an important relationship with algae proliferation rate, an aquatic ecological model with temperature effect is proposed to explore the coexistence modes of algae and biological manipulation predator and investigate how temperature affects their dynamic evolution. Suppose temperature parameter R is a constant variable, which can approximately describe ecological culture system of algae and biological manipulation predator under laboratory conditions if temperature parameter R is a periodic function variable, which mainly represents the natural ecosystem of algae and biological manipulation predators in a naturally eutrophic lake.
Based on dynamic population theory, some threshold conditions are given to guarantee the existence and stability of all possible equilibrium points, and some critical conditions for the occurrence of transcritical bifurcation and Hopf bifurcation are also deduced. Furthermore, some key parameters affecting the dynamic evolution characteristics of the model (2.2) are found through theoretical derivation and numerical simulation. All in all, these results are the theoretical basis for subsequent numerical simulation work and abstractly display the influence of some parameters on the dynamic evolution of the model (2.2).
Through the numerical simulation test on dynamic behaviors of the model (2.1), the influence mechanism of temperature on the stable succession of aquatic ecosystem is discovered in Figure 6 and Figure 8, the coexistence mode of algae and biological manipulation predator can change from a constant steady-state mode to a periodic oscillation mode with the temperature increasing gradually, which also indirectly indicates that the appropriate temperature range is one of the key factors for algae and biological manipulation predator to form a stable coexistence mode, and the periodic oscillation coexistence mode is more favorable to control the growth rate of algae population by biological manipulation. Based on the bifurcation dynamics evolution analysis of the model (2.2), it is worth pointing out that transcritical bifurcation can induce the appearance of the internal equilibrium point
, which represents the coexistence of algae and biological manipulation predator in a periodic oscillation mode, and completely changes the dynamic coexistence nature of algae and biological manipulation predator. Furthermore, when the value of control parameter p decreases and falls below a critical threshold, the coexistence mode of algae and biological manipulation predator has changed fundamentally again, periodic oscillation coexistence mode will replace the constant steady-state coexistence mode through a Hopf bifurcation. These results directly show that the value of key parameter p plays an important role in the bifurcation dynamic behavior evolution of the model (2.2). In general, some theoretical and numerical results obtained in this paper can provide a certain theoretical basis for the formation of a healthy and stable aquatic ecosystem and also provide certain numerical support for the feasibility of biological manipulation technology.
In the follow-up research works, firstly, we will introduce Arrhenius exponential temperature function and partial normal distribution temperature function into ecological modeling and investigate the impact of different temperature function manifestations on the dynamic behavior of the model (2.1). Secondly, we will continue to deepen the environmental impact factors of such aquatic ecological models and then further explore the influence of various environmental factors on the dynamic relationship between algae and biological manipulation predators. Finally, we will further explore the dynamic pattern behavior of the model (2.1) with the help of these papers [25] [26] [27]. In a word, all these results are expected to be useful in studying the dynamic behavior of aquatic ecosystems.
Acknowledgements
This work was supported by the key International Cooperation Projects of China (Grant No: 2018YFE0103700), the National Natural Science Foundation of China (Grant No: 31570364 and 61871293), key R&D Program Projects in Zhejiang Province of China (Grant No: 2021C03166), by the Science and Technology Program of Wenzhou, China (Grant No.X20210097).