A dynamic model describing heterotrophic culture of chorella and its stability analysis

. Chlorella is an important species of microorganism, which includes about 10 species. Chlorella USTB01 is a strain of microalga which is isolated from Qinghe River in Beijing and has strong ability in the utilization of organic compounds and was identiﬁed as Chlorella sp. (H. Yan etal, Isolation and heterotrophic culture of Chlorella sp., J. Univ. Sci. Tech. Beijing , 2005, 27 :408-412). In this paper, based on the standard Chemostat models and the experimental data on the heterotrophic culture of Chlorella USTB01, a dynamic model governed by diﬀerential equations with three variables ( Chlorella , carbon source and nitrogen source) is proposed. For the model, there always exists a boundary equilibrium, i.e. Chlorella -free equilibrium. Furthermore, under additional conditions, the model also has the positive equilibria, i.e., the equilibira for which Chlorella , carbon source and nitrogen source are coexis- tent. Then, local and global asymptotic stability of the equilibria of the model have been discussed. Finally, the parameters in the model are determined ac- cording to the experimental data, and numerical simulations are given. The numerical simulations show that the trajectories of the model ﬁt the trends of the experimental data well.


1.
Introduction. Chlorella is an important species of microorganism, which includes about 10 species (see, for example, [3]). A single Chlorella cell is of sphere with the diameter from 2µm to 12µm. Chlorella has the following characteristics. (a) It grows quickly and can efficiently photosynthesize. (b) It includes rich protein, highly unsaturated fatty acids, carotenoid, lutein, astaxanthin and a variety of vitamins (see, for example, [3], [36] - [38] and [40]). In particular, lutein in Chlorella has the functions of restraining the decreasing of sight and the blindness caused by age and of inhibiting on the growth of tumor etc. (see, for example, [12], [13], [15], [19] and [23]). In addition, Chlorella also contains Chlorella growth factor which has the effect of increasing immune response of T cells and B cells. Chlorella also has the effect of decomposing and eliminating harmful substances in environment (see, for example, [10] and [33]). Furthermore, Chlorella has the effect of absorbing the heavy metal copper etc.. Hence, in environmental science field, Chlorella is usually used to remove organic pollutants and heavy metals (see, for example, [34] and [35]).
1.1. Culture of Chlorella. The culture of Chlorella is divided into autotrophic culture and heterotrophic culture. Heterotrophic culture is divided into batch culture, feed-batch culture and continuous culture. Batch culture is extensively used in applications. In batch culture, microorganisms are put into the culture vessel which contains certain nutrients, under appropriate temperature, salinity and pH. Then, after a period of time of culture, the microorganisms in the culture vessel are harvested one-time. Batch culture is simple, but it is usually difficult to achieve higher biomass of microorganisms. Furthermore, the increasing of the initial concentration of nutrients may result in the inhibition on the growth of microorganisms. In feed-batch culture, nutrients are added into the culture vessel for every fixed time. Hence, feed-batch culture can reduce the inhibition on the growth of microorganisms and achieve higher biomass of microorganisms. Continuous culture is divided into semi-continuous culture and continuous culture. Semi-continuous culture is that nutrients are added into the culture vessel and, at the same time, microorganism is harvested for every fixed time. In continuous culture, nutrients are put into the culture vessel continuously and, at the same time, microorganisms flow out continuously. Hence, a balance between the inflow rate and the outflow rate can be maintained. Furthermore, continuous culture can also decrease the inhibition on the growth of microorganism caused by higher concentration of nutrients and accumulation of harmful substances in the culture vessel.
1.2. Heterotrophic culture of Chlorella sp. USTB01. It has been found that some species of Chlorella can not only be cultured in suitable organisms even in the absence of light (see, for example, [8] and [30]), but also contain richer lutein (see, for example, [8], [21] and [27]). Usually, carbon source and nitrogen source are main nutrients in culture of Chlorella. In [25], [26] and [31], the impact of concentration of glucose on the growth of Chlorella has been studied. However, till now, there are few species of Chlorella which can grow quickly (see, for example, [18], [22], [36] and [38]).
Recently, a species of Chlorella, named Chlorella sp. USTB01, which can be cultured quickly by the method of heterotrophic culture has been successfully sieved out from Qinghe river in Beijing (see, for example, [36]). Then, the impact of carbon-to-nitrogen mass ratio on the growth of Chlorella sp. USTB01 has also been investigated in details in [31], [36] and [37]. In the experiments, nitrogen sources are urea, ammonium chloride and potassium nitrate, and the culture time is limited to 36 hours. The experimental data show that, while carbon-to-nitrogen mass ratio is 20:1 and nitrogen source is ammonium chloride, the biomass of Chlorella sp. USTB01 in the culture vessel is higher (see Fig.1), but the protein contained is lower. If carbon-to-nitrogen mass ratio is still 20:1, but nitrogen source is replaced by potassium nitrate, then, the biomass of Chlorella sp. USTB01 in the culture vessel is not higher (see Fig.1), but the protein contained is higher.
Furthermore, if carbon-to-nitrogen mass ratio is 10:1 and nitrogen source is ammonium chloride, the biomass of Chlorella USTB01 in the culture vessel is higher (see Fig.2), and the protein contained is also higher.   The purpose of the paper is to model the growth of Chlorella USTB01 by constructing a dynamic model described by differential equations and to predict the growth status of Chlorella sp. USTB01 for longer culture time. The organization of the paper is as follows.
In Section 2, firstly, standard Chemostat models are introduced. Secondly, a dynamic model governed by differential equations is proposed based on standard Chemostat models and the experimental data obtained for heterotrophic culture of Chlorella sp. USTB01. Thirdly, for the model, global existence and nonnegativity of the solutions with appropriate initial conditions and classification of the equilibria are discussed. In Section 3, local and global asymptotic stability of two classes of the equilibria (i.e., Chlorella-free equilibrium and the equilibria for which Chlorella, carbon sources and nitrogen sources are coexistent) of the model in Section 2 are discussed by stability theory of differential equations. In Section 4, the parameters in the model in Section 2 are determined according to the experimental data and then numerical simulations are carried out. The numerical simulations suggest the following interesting facts. (i) The trajectories of the model in Section 2 fit the experimental data very well. (ii) Chlorella sp. UTSB01 is in the exponential growth phase from 36 hours to 40 hours. (iii) If the culture time is extended to 60 hours, the biomass of Chlorella sp. USTB01 in the culture vessel can reach at least 95% of the theoretical biomass. The proofs of the main results in Section 3 are finally given in Appendices A and B.

2.1.
Chemostat models. Dynamical models which describe continuous culture of microorganisms are also called Chemostat models (see, for example, [28]). Chemostat models can be used to the researches on the growth of algae and plankton in lakes and oceans, as well as the growth of endangered species (see, for example, [4] and [28]). The standard Chemostat model with a single limiting resource and a single microorganism is usually described by the following differential equations (see, for example, [11] and [28]), where S(t) and X(t) are the concentrations of the limiting resource (substrate) and microorganism at time t, respectively. S 0 > 0 denotes the concentration of input limiting resource and is assumed to be constant. The constant D > 0 is the dilution rate. The death rate of microorganism is assumed to be insignificant compared to the dilution rate and is ignored. The function µ(S) denotes how X consumes S and is called the uptake function. In most cases, µ(S) is chosen as Monod (1950) function, that is µ(S) = µ m S/(K m + S), where µ m is the maximum growth rate of microorganism, K m > 0 is the Michaelis-Menten (or half-saturation) constant. In more general case, µ(S) is assumed to be a monotone function which satisfies the conditions: is called the yield coefficient with respect to the nutrient S. Usually, δ(S) is chosen as constant or linear function of the nutrient S: δ(S) = AS +B, where A ≥ 0 and B > 0 are constants (see, for example, [4], [24], [28], [29] and the references therein). Based on different biological meanings, the model (1) has been modified and extended in the following aspects in recent years (see, for example, [1], [2], [5] - [7], [11], [16], [17], [20], [24], [28], [29], [32] and the references therein): (i) The models in which multiple microorganisms compete for single limiting resource, and the dilution rates of microorganisms may be different; (ii) The models in which single or multiple microorganisms compete for multiple resources, and the resources are perfectly substitutable or perfectly complementary; (iii) The models in which uptake functions are nonmonotone and the inhibiting effects of higher concentrations of substrates to the growth of microorganism are considered; (iv) The models in which there are time delays in the nutrients conversion processes; (v) The models in which the effects of inhibitors (internal or external) on microbial competition (such as plasmid-bearing, plasmid-free competition etc) are considered. (vi) The models in which the input rates of limiting resources are not constants (such as periodic functions of time t etc).

2.2.
Dynamical model in the heterotrophic culture of Chlorella. Chlorella sp. USTB01 ( Chlorella in short hereafter) has a strong ability in the utilization of organic compounds and includes rich protein (see, for example, [36]). In this subsection, based on the basic Chemostat model (1) and the experimental data in [36] in the heterotrophic culture of Chlorella, a dynamical model is proposed. Let us first look back the experimental procedures in [36]. (i) First, an appropriate amount of Chlorella is added into the culture vessel which has definite volume and is fed with an appropriate amount of carbon source and nitrogen source. (ii) Then, after the beginning of the culture, an appropriate amount of carbon source and nitrogen source are added into the culture vessel in batches (every 6 hours), and the quantity of the carbon source and nitrogen source added are dependent on the consumption of the carbon source and nitrogen source and the biomass of Chlorella in the culture vessel. Furthermore, the quantity of the carbon source and nitrogen source added every 6 hours is not less than the quantity added at the beginning of the experiment, and proportional to the biomass of Chlorella in culture vessel (see Fig 3). (iii) The culture time is limited 36 hours and the outflow of carbon source, nitrogen source and Chlorella in the culture vessel are not considered. (iv) The experimental data show that the mortality of Chlorella in the culture vessel increases as the biomass of Chlorella increases (i.e., product inhibition). Hence, the growth of Chlorella in the culture vessel is density-dependent.  Let X(t), C(t) and N (t) denote the mass of Chlorella, carbon source and nitrogen source, respectively, in the culture vessel at time t. We would point out here that, for simplicity, it is assumed that input of carbon source and nitrogen source is continuous. Hence, the dynamic model describing the growth of Chlorella in the heterotrophic culture can be expressed as following, where the constant d 1 > 0 is the mortality rate of Chlorella. The constants d 2 > 0 and d 3 > 0 are the attrition rates of carbon source and nitrogen source, respectively. The term f X 2 in which f > 0 is constant indicates that the growth of Chlorella in the culture vessel is density-dependent. The terms αX and βX in which α > 0 and β > 0 are constants indicate that the quantity of carbon source and nitrogen source added every 6 hours are proportional to the biomass of Chlorella in culture vessel, respectively. The constants a > 0 and b > 0 represent the constant input of carbon source and nitrogen source, respectively. Since carbon source and nitrogen source are perfectly complementary in the culture of Chlorella, the term rµ 1 (C)µ 2 (N ) represents the growth rate of Chlorella, and the terms r 1 µ 1 (C)µ 2 (N )X and r 2 µ 1 (C)µ 2 (N )X represent the quantity of decreasing of carbon source and nitrogen source, respectively, where r > 0, r 1 > 0 and r 2 > 0 are constants, and the functions µ 1 (C) and µ 2 (N ) are nonnegative and continuous for C ≥ 0 and N ≥ 0. δ i (i = 1, 2) are yield coefficients, which are defined as δ i = mass of organism f ormed mass of substrate consumed (i = 1, 2).
The constants θ 1 > 0, θ 2 > 0 and θ 3 > 0 are relative growth rates of Chlorella, carbon source and nitrogen source, respectively. For simplicity of theoretical analysis, in this paper, the functions µ 1 (C) and µ 2 (N ) are chosen as Monod type functions, i.e., where k 1 > 0 and k 2 > are the half-saturation constants with respect to carbon source and nitrogen source, respectively. The yield coefficients δ 1 and δ 2 are assumed to be constants. Hence, without loss of generality, r 1 /δ 1 and r 2 /δ 2 are still denoted by r 1 and r 2 , respectively. Therefore, the dynamic model (2) can be rewritten in the following simpler form, 2.3. Existence of equilibria. In this subsection, let us consider nonnegativity of the solutions and existence of the equilibria of (3).
First, according to biological meanings, the initial condition of (3) is given as where the constants X 0 , C 0 and N 0 represent the initial mass of Chlorella, carbon source and nitrogen source, respectively. From Lemma 1 in Appendix B and standard theory on existence of solutions of ordinary differential equations (see, for example, [9]), it can be easily shown that the solution (X(t), C(t), N (t)) of (3) with the initial condition (4) is existent, unique and nonnegative for all t ≥ 0. Let (X, C, N ) be any equilibrium of (3). Then, (X, C, N ) satisfies the following nonlinear algebraic equations, There are two cases to be discussed.
(i) (3) always has the boundary equilibrium E 0 = (X, C, N ) = (0, a/d 2 , b/d 3 ). The existence of E 0 indicates that, if there is no Chlorella to be added into the culture vessel at the beginning of the culture, the concentrations of carbon source and nitrogen source in the culture vessel always maintain the constants values a/d 2 and b/d 3 , respectively. The equilibrium E 0 is also called Chlorella-free equilibrium.
(ii) Let (X, C, N ) be any equilibrium with X > 0, C > 0 and N > 0. From (5), Clearly, X should satisfy the conditions and br − [(r 2 d 1 − βr)X + r 2 f X 2 ] > 0. (8) Substituting the second and the third equations of (6) into the first equation of (6) gives a fifth order algebraic equation. Hence, under suitable conditions, there may be at most five different positive roots for the fifth order algebraic equation. Let X = X * be any such positive root which also satisfies the conditions (7) and (8). Then, from (6), C = C * > 0 and N = N * > 0 can be obtained. Therefore, the dynamical model (3) at most has five equilibria of the type of E * = (X, C, N ) = (X * , C * , N * ) with X * > 0, C * > 0 and N * > 0. The equilibrium E * = (X * , C * , N * ) is also called positive equilibrium which indicates that Chlorella, carbon source and nitrogen source are coexistent for any time t ≥ 0.
From the first equation in (6), we obtain the function which implies that d 1 < r is a necessary condition for a positive equilibrium to exist. In fact, E 0 is globally asymptotically stable when d 1 ≥ r (see Theorem 3.3 in Section 3 below). General speaking, it is not easy to give a complete theoretical analysis for the existence of all the positive equilibria. But, according to the experiment data provided in [36], the approximate estimations of all the positive equilibria can be obtained by Matlab programming.
3. Stability analysis of equilibria. Stability analysis of the equilibria of the dynamic model (3) is very important for understanding the growth status of Chlorella. In this section, some sufficient conditions are given to ensure global asymptotic stability of boundary equilibrium E 0 = (0, a/d 2 , b/d 3 ) and local asymptotic stability of the positive equilibrium E * = (X * , C * , N * ) by standard stability theory of ordinary differential equations (see, for example, [9], [28]).
First of all, for local asymptotic stability of E 0 = (0, a/d 2 , b/d 3 ), we have the following then the boundary equilibrium E 0 = (0, a/d 2 , b/d 3 ) is locally asymptotically stable.

The proof of Theorem 3.2 is also given in Appendix A.
Remark 2. Note that sufficient conditions for (10) to hold are a ≥ d 2 C * and b ≥ d 3 N * . Hence, for fixed constants a > 0 and b > 0, if the attrition rates d 2 and d 3 of carbon source and nitrogen source, respectively, are small enough and the growth rate r of Chlorella is large enough, the conditions (7), (8) and (10) can be satisfied. Therefore, it follows from Theorem 3.2 that the positive equilibrium E * = (X * , C * , N * ) is locally asymptotically stable. This implies that Chlorella, carbon source and nitrogen source in the culture vessel are coexistent and their masses tend to constant values as time t increases.
Remark 3. There may exist multiple positive equilibria for the dynamical model (3). Since complexity of the expression of the positive equilibrium E * = (X * , C * , N * ), detailed theoretical and numerical analysis on multi stability are omitted. In Section 4 below, according to the experimental data in [36], the parameters in the dynamical model (3) are chosen, and then the value of the positive equilibrium E * = (X * , C * , N * ) is computed numerically. Theorem 3.2 and numerical simulations show that the positive equilibrium E * = (X * , C * , N * ) is asymptotically stable and that the trajectories of the dynamical model (3) fit the trends of the experimental data well.
From the point of view in both mathematics and biology, theoretical analysis on global asymptotic stability of the equilibria of the dynamical model (3) has more important significance. Global asymptotical stability of the equilibria E 0 or E * implies that the asymptotic properties of Chlorella, carbon source and nitrogen source in culture vessel are not dependent on the initial values X 0 , C 0 and N 0 of Chlorella, carbon source and nitrogen source.
For global asymptotically stability of the boundary equilibrium E 0 , we have the following Theorem 3.3. If d 1 ≥ r, or d 1 < r and then, the boundary equilibrium E 0 is globally asymptotically stable.
The proof of Theorem 3.3 is given in Appendix B.
Remark 4. The condition d 1 ≥θ in Theorem 3.3 implies the condition d 1 > θ in Theorem 3.1. Furthermore, numerical simulations strongly suggest that the boundary equilibrium E 0 should also be globally asymptotically stable even if the condition d 1 ≥ θ holds. Therefore, the condition (11) in Theorem 3.3 may be further improved.
As for global asymptotic stability of the positive equilibrium E * , detailed discussions shall be given in other paper.

Numerical simulations and discussions.
4.1. Numerical simulations. In this subsection, let us discuss how the trajectories of the dynamical model (3) fit the experimental data in [36] based on Theorem 3.2 and numerical simulations. To observe the growth status of Chlorella in the culture vessel, two sorts of nitrogen sources (i.e., ammonium chloride and potassium nitrate) are used. Carbon source and nitrogen source are added into the culture vessel for every 6 hours with two kinds of carbon-to-nitrogen mass ratios 20:1 and 10:1.
Case (IV) Nitrogen source is potassium nitrate and carbon-to-nitrogen mass ratio is 10:1. The parameters in (3) Fig.7 shows that the numerical simulation trajectory of (3) with the initial value (X 0 , C 0 , N 0 ) = (1.8, 0, 0) fits the experimental data in [36] well. the growth of Chlorella is in delay phase from 0 hour to 12 hours, in the accelerated phase from 12 hours to 24 hours, in exponential growth phase from 24 hours to 40 hours, and in deceleration phase from 40 hours to 60 hours. Furthermore, we have from Fig.s 4-7 that, for Cases (I), (II), (III) and (IV) in Subsection 4.1, the saturated biomass of Chlorella are about X * ≈ 62.5998, X * ≈ 47.3691 , X * ≈ 58.1376 and X * ≈ 58.4398, respectively. When the culture time is limited 36 hours (see, for example, [36]), the experimental data shows that the biomass of Chlorella in the culture vessel are approximately 53.51% (33.5 Remark 5. The purpose of the paper is to model the growth of Chlorella USTB01 which have some special characteristics in applications and proposes an ordinary differential equation model (3) according to the experimental procedures in [36]. The model (3) can be revised as the following more general form with time delays, In (12), the constants ρ 1 ≥ 0 and ρ 2 ≥ 0 are the rate constants at which carbon source and nitrogen source are recycled because of the death of Chlorella. The constant τ ≥ 0 denotes the time delay involved in the conversion of nutrients to Chlorella. The constant σ ≥ 0 is a fixed time during which carbon source and nitrogen source are released completely from dead Chlorella. The factor e −d1τ represents the approximate proportion of Chlorella that remain in the culture vessel during the conversion process. Theoretical analysis on stability of the equilibria of (12) and its applications in the research of the control strategies of algal toxins of blooms of the local lakes and rivers in Beijing city shall be given in the other papers.
Proof of Theorem 3.2. The corresponding characteristic equation at the positive equilibrium E * = (X * , C * , N * ) can be expressed as the following form, where a 1 = −(a 11 +a 22 +a 33 ), a 2 = a 11 a 22 +a 11 a 33 +a 22 a 33 −a 12 a 21 −a 13 a 31 −a 23 a 32 , a 3 = a 11 a 23 a 32 + a 12 a 21 a 33 + a 13 a 31 a 22 − a 11 a 22 a 33 − a 12 a 23 a 31 − a 13 a 21 a 32 . Here, the definition of a ij is the same as that in the proof of Theorem 3.1, but E(X, C, N ) is replaced by E * (X * , C * , N * ). Since we have that Hence, Here, the assumptions in Theorem 3.2 and the following equalities been used. Hence, a 2 > 0. Furthermore, we have that we have that X(t) ≡ 0 for any t ∈ [0, ρ) if X 0 = 0, and that X(t) > 0 for any t ∈ [0, ρ) if X 0 > 0. If there is a t 2 ∈ [0, ρ) such that C(t 2 ) = 0, we have from the second equation of (3) thatĊ(t 2 ) = αX(t 2 ) + a > 0, which clearly implies that C(t) ≥ 0 for any t ∈ [0, ρ). Similarly, we have that N (t) ≥ 0 for any t ∈ [0, ρ).
For t ∈ [0, ρ), we have thaṫ Consider the comparison system, Since the solution of the comparison system is existent and nonnegative for any t ≥ 0, we have from comparison principle that X(t) ≤ Z(t) for any t ∈ [0, ρ). Thus, X(t) ≤ M 0 for some positive constant M 0 and all t ∈ [0, ρ). Therefore, for t ∈ [0, ρ), we have from (3) thatĊ Again, we have from (13), (14) and comparison principle that C(t) and N (t) are bounded for t ∈ [0, ρ). Hence, the solution (X(t), C(t), N (t)) can be continued to +∞ and is unique and nonnegative. Note that lim t→+∞ Z(t) = 0 for r ≤ d 1 , and lim t→+∞ Z(t) = M for r > d 1 . Thus, we have from (13), (14) and comparison principle that, for r ≤ d 1 , lim t→+∞ X(t) = 0, lim t→+∞ C(t) = a/d 2 and lim t→+∞ N (t) = b/d 2 ; and that, for r > d 1 , lim t→+∞ X(t) ≤ (r −  Proof. According to Lemma 4.1, it only needs to show that X(t) ≤ M , C(t) ≤ M 1 and N (t) ≤ M 2 for any t ≥ 0 and (X(0), C(0), N (0)) ∈ G. In fact, if there exits some t 3 > 0 such that X(t 3 ) > M , then we have from the Lagrange mean-value theorem that there is some t 4 ∈ (0, t 3 ) such that X(t 4 ) > M andẊ(t 4 ) > 0. On the other hand, we have from (3) thaṫ which is a contradiction. Thus, X(t) ≤ M for any t ≥ 0. Therefore, we have from (3) that for any t ≥ 0, from which we easily have that C(t) ≤ M 1 and N (t) ≤ M 2 for any t ≥ 0. This completes the proof of Lemma 4.2.
Proof of Theorem 3.3. If r ≤ d 1 , we have from Lemma 4.1 and Theorem 3.1 that the non-Chlorella equilibrium E 0 is globally asymptotically stable. Let us consider the case of r > d 1 . Since the condition d 1 ≥θ in Theorem 3.3 implies the condition d 1 > θ in Theorem 3.1, hence, the Chlorella-free equilibrium E 0 is locally asymptotically stable. Consider the nonnegative function V = x on G. It is clear that V = x is continuous onḠ and that the derivative along the solutions of (3) satisfieṡ Hence, V = x is a Liapunov function of (3) on G.
Define the subset of G as E = {(X, C, N )|(X, C, N ) ∈ G,V = 0}, and let H be the largest invariant set of (3) in E. Since E 0 ∈ H, we have that H is nonempty and that E ⊆ {(X, C, N ) ∈ G|X = 0}. From the invariance of H and (3), it is easy to show that H = {E 0 }. Therefore, it follows from well known Lyapunov-LaSalle invariance principle that E 0 is globally attractive. This completes the proof of Theorem 3.3.