Stability and persistence analysis of a microorganism flocculation model with infinite delay

1.   Introduction

In the research of microbial growth dynamics, a rich variety of differential equation models and theoretical research of chemostats have been greatly developed, and fruitful research results have been achieved ^[1,2,3,4]. Flocculants are widely used in the sedimentation/coagulation/collection of microorganisms in petri dishes/reaction tanks, such as the removal of harmful microorganisms in polluted water and the development and utilization of beneficial microorganisms (such as microbial algae health food and algae-based microbial clean energy). Therefore, the research and development of flocculants (organic, inorganic, microbial, etc.) with the advantages of high efficiency, no pollution and low cost remains one of the important research topics in the field of microbial application technology. In literature ^[5,6,7,8,9], a dynamics model of microbial continuous culture with flocculation, time delay and different functional reaction functions was established based on the classical chemostat model, and the existence, local and global stability of the forward and backward bifurcation of the equilibrium, as well as the consistent persistence of the model were studied by constructing an appropriate Lyapunov functional and other methods.

Based on the microorganism flocculation models in ^[6,7,8] with discrete time delays, in this paper, we intend to further consider the following microorganism flocculation model with infinite delay:

Here, $n(t)$, $x(t)$ and $s(t)$ represent the concentrations of nutrients, microorganisms and flocculant at time $t$, respectively. The positive constants $n^{0}$ and $s^{0}$ indicate the input concentrations of nutrients and flocculant, respectively. The positive constant ${D}$ indicates the outflow and inflow rate of the substances in the chemostat. The non-negative constants $h_{1}$ and $h$ indicate the consumption ratio of nutrients and the growth ratio of microorganisms, respectively. The non-negative constants $h_{2}$ and $h_{3}$ indicate the flocculating ratio of microorganisms and the consumption ratio of flocculant, respectively. The distribution function $f(\theta)$ is non-negative

and satisfies $\int_{0}^{\infty} f(\theta)d\theta = 1$. $\mu$ is a non-negative constant, and $e^{-\mu\theta}$ is the survival rate of microorganisms.

As we know, infinite delay was introduced early into dynamics models with a single nutrient and single-species microorganism growth problems. It can describe the cumulative cycle of microbial decomposition into nutrients (from an early time in the past) and the accumulated time it takes for nutrients to be absorbed/stored by microorganisms until they are converted into actual biomass ^{[10,11,12,13,14,15,16]}. Further, its use can be extended to the models with multiple microorganisms, substitutable nutrients, non-substitutable nutrients (fully complementary nutrients), etc. ^[17,18,19]. As pointed out in ^[15], the differential equations with distributed infinite delay have been successfully applied in biological modeling since the work of Volterra, and they are regarded to be more realistic than the differential equations with discrete and finite delay.

Similarly, infinite delay is also widely used to characterize the dynamic modeling of practical issues in many fields, such as population growth ^{[20,21,22,23,24]}, infectious diseases in population ^[25,26], virus infection and immunity ^[27,28,29].

The main purpose of this paper is to analyze the local and global stability of the boundary equilibrium (microorganism-free equilibrium) and the positive equilibrium (microorganism co-existent equilibrium) of model (1.1) and provide the corresponding criteria by using the stability theory of infinite delay differential equations and other related methods. At the same time, for the persistence of model (1.1), an explicit estimate of the eventual lower bound of its arbitrary positive solutions is given.

2.   Dissipativeness and classification of the equilibria

For simplicity, we define the following parameter transformation:

The dimensionless system for (1.1) is given by

Here, for the convenience of writing, the corresponding symbols are still used after the "overline" is removed. We assume that $\int_{0}^{\infty} e^{-\mu \theta} f(\theta)d\theta = q\leq1$ and $\int_{0}^{\infty } \theta e^{-\mu\theta}f(\theta)d\theta < \infty$.

According to ^{[23,30,31,32,33,34]}, the phase space of model (2.1) can be taken as

where $\alpha$ is a positive constant; the norm $\|\phi\| = \sup_{\theta \leq 0}|\phi(\theta)| e^{\alpha \theta}$. The nonnegative cone of $\mathcal{C}$ is $\mathcal{C}_{+} = \left\{\phi \in \mathcal{C} : \phi(\theta)\geq0, \; \theta \in(-\infty, 0]\right\}$.

For $\phi\in\mathcal{C_+}$, we define $\phi_{t}\in\mathcal{C_+}$, $\phi_{t}(\theta) = \phi(t+\theta), \; \theta \in(-\infty, 0]$.

Considering the biological significance of model (2.1), for any initial function $\phi = (\phi_1, \phi_2, \phi_3)^T$, we restrict $\phi = (\phi_1, \phi_2, \phi_3)^{T}\in BC_+ \subset \mathcal{C_+}$, where

Thus, we have the following theorem.

Theorem 2.1. For $\phi = (\phi_1, \phi_2, \phi_3)^{T}\in BC_+$, the solution of model (2.1) with the initial condition $(N(\theta), X(\theta), S(\theta))^{T} = (\phi_1(\theta), \phi_2(\theta), \phi_3(\theta))^{T}$ $(\theta \in(-\infty, 0])$ is nonnegative and unique, and it satisfies $N(t) > 0(t > 0)$, $X(t)\geq0(t\geq0)$, $S(t) > 0(t > 0)$, $\underset{t\rightarrow \infty }{\lim \sup } N(t)\leq 1$, $\underset{t\rightarrow \infty }{\lim \sup} X(t)\leq \xi,$ $\underset{t\rightarrow \infty }{\lim \sup } S(t)\leq 1$,

where $\xi$ is a positive constant not depending on the initial function.

Proof. From the first and third equations in model (2.1), we can easily get $\underset{t\rightarrow \infty }{\lim \sup } N(t)\leq 1$ and $\underset{t\rightarrow \infty }{\lim \sup } S(t)\leq 1$. Define

Here, the convergence of $\int_{0}^{\infty } \theta e^{-\mu \theta} f(\theta)d\theta$ should be noticed. Then, for $t\geq0$, the derivative of $V$ along any solution of model (2.1) is

From Theorem A in ^[35], we obtain that the solution of model (2.1) is ultimately uniformly bounded. Furthermore, there is a positive constant $\xi$ independent of the initial function such that $\underset{t\rightarrow \infty }{\lim \sup } X(t)\leq \xi$.

According to the discussion in ^[5,6], the equilibria of model (2.1) can be classified as follows.

The model (2.1) always has a boundary equilibrium (microorganism-free equilibrium) $E_{0} = (1, 0, 1)^{T}$.

To consider the existence of the positive equilibria (microorganism co-existent equilibrium) for $X > 0$, we consider the equation of the variable $X$:

where $R_{0} = \frac{rq-1}{h_{2}}, \; a = h_{2}h_{3}, \; b = h_{1}h_{2}+h_{1}.$ If $rq\leq 1(R_{0}\leq0)$, (2.1) has no positive equilibria. If $rq > 1$, the following conclusion holds:

1) $a\leq b$. If $R_{0} > 1$, there exists a unique positive equilibrium $E_{1}^{\ast } = (N_{1}^{\ast}, X_{1}^{\ast}, S_{1}^{\ast})^{T}$. If $R_{0}\leq1$, there are no positive equilibria.

2) $a > b$. If $R_{0} > 1$, there exists a unique positive equilibrium $E_{1}^{\ast } = ({N_{1}}^{\ast }, X_{1}^{\ast }, S_{1}^{\ast })^{T}$.

If $R_{0} < 1$, let $F(z) = a^{2}\left[(z-\frac{b}{a})^{2}+\frac{4h_{1}}{a}(z-1)\right]$. By $F(z) = 0$, we have

and we know that $z_{1}$ satisfies $\frac{b}{a} < z_{1}\equiv \omega < 1$.

If $R_{0} = \omega$, there exists a unique positive equilibrium $E_{\omega}^{\ast} = (N_{\omega}^{\ast}, X_{\omega}^{\ast}, S_{\omega}^{\ast})^{T}$.

If $\omega < R_{0} < 1$, there exist two positive equilibria $E_{1}^{\ast} = (N_{1}^{\ast}, X_{1}^{\ast}, S_{1}^{\ast})^{T}$, $E_{2}^{\ast} = (N_{2}^{\ast}, X_{2}^{\ast}, S_{2}^{\ast})^{T}$.

If $R_{0} < \omega$, there are no positive equilibria.

3.   Stability of the boundary equilibrium

This section mainly considers the local and global stability of the boundary equilibrium $E_0$.

3.1.   Local stability of the boundary equilibrium

The linearized system of model (2.1) at any equilibrium is

Then, we have the following theorem.

Theorem 3.1. If $R_{0} < 1$, then the boundary equilibrium $E_{0}$ is locally asymptotically stable. If $R_{0} > 1$, then the boundary equilibrium $E_{0}$ is unstable.

Proof. The characteristic equation at the boundary equilibrium $E_{0}$ is

There are always eigenvalues $\lambda _{1, 2} = -1 < 0$. Consider the following transcendental equation:

Let us assume that Eq (3.2) has a root $\lambda$ with a non-negative real part, i.e., $\text{Re}\lambda \geq0$; this then implies that

Obviously, the inequality Eq (3.3) contradicts $R_0 < 1$. Therefore, all roots of Eq (3.2) have negative real parts when $R_0 < 1$.

If $R_0 > 1$, we see that $\psi(0) = h_{2}(1-R_{0}) < 0$. For $\lambda\geq0$, we have

Therefore, $\psi (\lambda)\rightarrow +\infty$ as $\lambda\rightarrow +\infty$. According to the intermediate value theorem of continuous functions, the equation $\psi (\lambda) = 0$ has at least one positive eigenvalue $\lambda$.

3.2.   Global stability of the boundary equilibrium

Theorem 3.2. If $R_0 < \frac{1}{1+h_{3}\xi}$ $(<1)$, the boundary equilibrium $E_{0}$ is globally asymptotically stable.

Proof. Noting Theorem 3.1, we only need to prove that the boundary equilibrium $E_{0}$ is globally attractive.

For the arbitrary bounded solution of model (2.1), we consider the following differentiable function

Since $\underset{t\rightarrow \infty }{\lim \sup} X(t)\leq \xi$, then the third equation from model (2.1) yields $\underset{t\rightarrow \infty }{\lim \inf} S(t)\geq \frac{1}{1+h_{3}\xi}.$ There exists a sufficiently large $t_1 > 0$; we have $S(t) > \frac{1}{1+{h_{3}}\xi}-\varepsilon > 0\;$for $t > t_1.$ Let $R_{0} < \frac{1}{1+{h_{3}}\xi}-\varepsilon > 0$ for sufficiently small $\varepsilon > 0$. Therefore, for $t > t_1$,

Integration at both sides of the above inequality yields

We can obtain $X(t)$ as integrable on $[t_1, +\infty)$. In addition, considering the initial function $\phi = (\phi_1, \phi_2, \phi_3)^{T}\in BC_+$ and the boundedness of the solution, we have that the derivative function $X'(t)$ is bounded on $[t_1, +\infty)$ from the second equation of model (2.1). Therefore, $X(t)$ is uniformly continuous on $[t_1, +\infty)$. Using the familiar Barbalet's lemma, we have $\underset{t\rightarrow \infty }{\lim} X(t) = 0$. Further, from the first and third equations of model (2.1), we can get $\underset{t\rightarrow \infty }{\lim} N(t) = \underset{t\rightarrow \infty }{\lim} S(t) = 1$.

4.   Stability of the positive equilibria

For the local asymptotic stability of the positive equilibrium $E_{1}^{\ast }$, there is the following theorem.

Theorem 4.1. The following conclusions hold$:$

$1)$ Let $a\leq b$ $($the forward bifurcation$)$. If $R_{0} > 1$, then the positive equilibrium $E_{1}^{\ast }$ is locally asymptotically stable.

$2)$ Let $a > b$ $($the backward bifurcation$)$.

$(i)$ If $R_{0} > 1$, then the positive equilibrium $E_{1}^{\ast }$ is locally asymptotically stable.

$(ii)$ If $\omega < R_{0} < 1$, then the positive equilibrium $E_{1}^{\ast }$ is locally asymptotically stable and the positive equilibrium $E_{2}^{\ast}$ is unstable.

$(iii)$ If $R_{0} = \omega$, then the positive equilibrium $E_{\omega}^{\ast}$ is linearly stable.

Proof. The characteristic equation of the linearized system at any positive equilibrium $E^{\ast }$ is

where

First, we consider the local asymptotic stability of $E_{1}^{*}$. Here, we notice the equations $N^{\ast} = \frac{1}{1+h_{1}X^{\ast}}$, $S^{\ast} = \frac{1}{1+h_{3}X^{\ast}}$ and $1+h_{2}S^{\ast } = rqN^{\ast }$. Similar to the discussion of Theorem 3.2 in ^[5], it is easy to obtain $J(0) = \tilde{J}(0) = N^{\ast }X^{\ast }H(X^{\ast })$, where $H(X^{\ast }) = h_{1}^{2}h_{3}X^{\ast 2}+2h_{1}h_{3}X^{\ast}+h_{3} +rqh_{1}-rqh_{3} > 0$. Therefore, $J(\lambda) = 0$ does not have zero root $\lambda = 0$.

Rewrite the equation $\tilde{{J}}(\lambda) = 0$ as

Let us assume that Eq (4.1) has roots $\lambda = x\pm iy$ $(x\geq0, \; y\geq0\; and\; |x|+|y| > 0)$. Taking the modulus on both sides of (4.1), it follows that

On the other hand, we have

where $G(x) = (x+1+h_{1}X^{\ast })(x+rqN^{\ast }+h_{3}X^{\ast}) -(rqN^{\ast})(x+1+h_{3}X^{\ast }).$ For $x > 0,$ it has from $G'(x) = 2x+1+h_{1}X^{\ast }+h_{3}X^{\ast } > 0$ that $G(x)\geq G(0) = N^{\ast}X^{\ast}H(X^{\ast}) > 0.$ Therefore, $F(x, y) > 0$ for $x\geq0$, $y\geq0$ and $|x|+|y| > 0$. Clearly, this is a contradiction. Thus, it is proved that all roots of Eq (4.1) have negative real parts.

Next, we consider the instability of $E_{2}^{\ast }$. Completely similar to the discussion of Theorem 3.3 in ^[5], at the positive equilibrium $E_{2}^{\ast }$, $\tilde{J}(\lambda)$ satisfies $\tilde{J}(0) < 0$.

Then for $\lambda\geq0$,

Therefore, $\tilde{J}(\lambda)\rightarrow +\infty$ as $\lambda\rightarrow +\infty$. Furthermore, it can be easily obtained that the equation $J(\lambda) = 0$ has at least one positive eigenvalue, implying that the positive equilibrium $E_{2}^{\ast}$ is unstable.

Finally, for the positive equilibrium $E_{\omega}^{\ast}$, completely similar to the discussion of Theorem 3.3 in ^[5], one obtains that the equation $\tilde{J}(\lambda) = 0$ has a single root $\lambda = 0$, and all other roots have negative real parts. Therefore, the positive equilibrium $E_{\omega}^{\ast}$ is linearly stable.

5.   Permanence of the model

For the initial function $\phi = (\phi_1, \phi_2, \phi_3)^{T}\in BC_+$ and $\phi_1\leq1$ and $\phi_3\leq1$, it is not difficult to show that the solution ${(N(t), X(t), S(t))^{T}}$ of model (2.1) satisfies ${N(t)}\leq1$ and ${S(t)}\leq1$ for $t\geq0$. Therefore, in the discussion of persistence, the initial function $\phi = (\phi_1, \phi_2, \phi_3)^{T}\in BC_+$ can be defined to satisfy $\phi_1\leq1$ and $\phi_3\leq1$.

Theorem 5.1. If $R_{0} > 1$, for $\phi = (\phi_1, \phi_2, \phi_3)^{T}\in BC_+$, $\phi_1\leq1$, $\phi_3\leq1$ and $\phi_2(0) > 0$, model (2.1) is uniformly persistent and the solution $(N(t), X(t), S(t))^{T}$ of model (2.1) satisfies

where $\tau_{1}$ is a positive constant, and

Proof. Let ${(N(t), X(t), S(t))^{T}}$ $(t\geq0)$ be a solution of model (2.1). By the conclusion of Theorem 2.1, we have $\underset{t\rightarrow \infty }{\lim \sup} X(t)\leq \xi;$ then, it follows from the first and third equations of model (2.1) that $\underset{t\rightarrow \infty }{\lim \inf} N(t)\geq \upsilon_1,$ $\underset{t\rightarrow \infty }{\lim \inf} S(t)\geq \upsilon_3$.

For $t\geq0,$ we define auxiliary functions $V(t)$:

Noticing that the boundedness of the solution ${(N(t), X(t), S(t))^{T}}$ $(t\geq0)$ and $\int_{0}^{\infty } \theta e^{-\mu\theta}f(\theta)d\theta < \infty$, we obtain that the function $V(t)$ is bounded and continuously differentiable. For $t\geq0$, we have

$S(t)\leq1$ $(t\geq0)$ is used here.

There exists a sufficiently large $T > 0$ such that $N(t)\geq \frac{1}{2(1+h_1\xi)}$ for $t\geq T$.

For every $t_{0}\geq T > 0$, it is shown below that the inequality $X(t)\leq X_{1}$ ($t\geq t_{0}$) cannot hold.

If not, then there exists some $t_{0}\geq T$ such that $X(t)\leq X_{1}$ for $t\geq t_{0}$. From the first equation of model (2.1), we have $N'(t)\geq 1-(1+h_{1}X_{1})N(t) = 1-\frac{1}{A}N(t)$ for $t\geq t_{0}$. Noticing that $\frac{1}{2(1+h_1\xi)}-A < 0$, therefore, for $t\geq t_{0}+T_{0}$, we obtain

Then, for $t\geq t_{0}+T_{0}$, we have $V'(t)\geq ({q}r\rho -1-h_{2})X(t)$.

It follows from the conditions in Theorem 5.1 that $r\int_{0}^{\tau _{1}}e^{-\mu \theta} f(\theta)d\theta > r({q}-\varepsilon) > r\rho ({q}-\varepsilon) > 1+h_{2}$.

Let $X_{m} = \underset{-\tau_{1}\leq \theta \leq 0}{\min } X(t_{0}+T_{0}+\tau _{1}+\theta);$ the following proof shows that $X(t)\geq X_{m}$ for all $t\geq t_{0}+T_{0}$.

If not, then there exists a $T_{1}\geq 0$ such that $X(t)\geq X_{m}$ for $t_{0}+T_{0}\leq t\leq t_{0}+T_{0}+\tau_{1}+T_{1}\equiv \bar{t}$, $X(\bar{t}) = X_{m}$ and $X'(\bar{t})\leq 0.$ It follows immediately from the second equation of model (2.1) that

$S(\bar{t})\leq1$ is used here. Clearly, this is a contradiction. It indicates that $X(t)\geq X_{m}$ for all $t\geq t_{0}+T_{0}$. Thus, for $t\geq t_{0}+T_{0}$, we have

Then, $V(t)\rightarrow +\infty$ as $t\rightarrow +\infty$, this is a contradiction with the boundedness of $V(t)$.

The following two cases are discussed.

Case (Ⅰ) The inequality $X(t)\geq X_{1}$ holds for all large $t$, implying that $\underset{t\rightarrow \infty }{\lim \inf }X(t)\geq X_1 \geq \upsilon_2$.

Case (Ⅱ) The function $X(t)$ oscillates infinitely above and below $X = X_{1}$ for all large enough t. At this time, we only consider interval $[t_{1}, t_{2}]$, where $t_{2} > t_{1}\geq T$, $X(t_{1}) = X_{1},$ $X(t_{2}) = X_{1},$ $X(t) < X_{1}$ $(t_{1} < t < t_{2})$.

If $t_{2}-t_{1}\leq T_{0}+\tau _{1}$, it follows from the second equation of model (2.1) that $X'(t)\geq -(1+h_{2})X(t)$ for $t\in [t_{1}, t_{2}]$. Therefore, by integrating, we have

If $t_{2}-t_{1} > T_{0}+\tau _{1}$, this case is exactly similar to the derivation of Case (Ⅰ), which make it easy to obtain $X(t)\geq \upsilon_2$ for $t\in \lbrack t_{1}, t_{1}+T_{0}+\tau_{1}]$.

It is further proved that there is still $X(t)\geq \upsilon_2$ for $t\in \lbrack t_{1}+T_{0}+\tau _{1}, t_{2}]$.

If this were not true, then there exists a $T_{2}\geq 0$ such that $X(t)\geq \upsilon_2$ for $t_{1}\leq t\leq t_{1}+T_{0}+\tau _{1}+T_{2}\equiv \hat{t}$, $X(\hat{t}) = \upsilon_2$ and $X'(\hat{t})\leq 0$. Similarly, it also follows from the second equation of model (2.1) that

$S(\hat{t})\leq1$ is also used here. Therefore, a contradiction is obtained. This shows that for $t\in \lbrack t_{1}, t_{2}]$, $X(t)\geq \upsilon_2$ always holds. Since the interval $[t_{1}, t_{2}]$ is arbitrarily chosen, this proves that for all sufficiently large $t$, there must be $X(t)\geq \upsilon_2$. Thus, it follows that $\underset{t\rightarrow \infty }{\lim \inf }X(t)\geq \upsilon_2.$

6.   Conclusions

In this paper, we have obtained the dissipation of the solution and the local asymptotic stability of the equilibrium as well as the uniform persistence of model (2.1) with infinite delay through fine analysis on the distribution of the roots of the characteristic equations in the complex plane and the asymptotic properties of the solutions, combined with the construction of appropriate Lyapunov functionals. Furthermore, a sufficient condition for the global stability of the boundary equilibrium is given. In particular, the conclusion of uniform persistence suggests that microorganism flocculation and collection are sustainable as long as the threshold parameter $R_0 > 1$. Also, the explicit estimation expressions for the eventual lower bound of the evolution of biomass (nutrient-$N (t)$, microorganism-$X (t)$, flocculant-$S (t)$) with time are given. These results further extend the corresponding conclusions in ^[5,6].

However, it should be pointed out that (2.1) is an infinite time-delay system, and the information of the initial function is always stored in the time-delay term, which makes it more difficult to obtain an explicit expression of $\xi$ in Theorem 2.1 than in the case of bounded time-delay systems, and we will leave this issue for further discussion. Second, the sufficient condition $R_0 < \frac{1}{1+h_{3}\xi}$ $(<1)$ given in Theorem 3.2 for the global stability of the boundary equilibrium (microorganism-free equilibrium) $E_0$ is also conservative, at least for the case of the forward bifurcation. In addition, we have not considered the global stability of the positive equilibrium $E^{*}$ (microorganism co-existent equilibrium) of the model (2.1), since new Lyapunov functionals with infinite delay may need to be constructed (see, for example, ^{[10,11,12,13,14,15]}). Furthermore, we see from Theorem 5.1 that the parameter $\xi$ in Theorem 2.1 plays an important role in the values of the parameters $\upsilon_1$, $\upsilon_2$ and $\upsilon_3$.

Finally, similar to ^[16], let us choose the following weak kernel and strong kernel:

The corresponding parameters $q$ can be easily calculated as

Therefore, the condition "$\int_{0}^{\tau_{1}} e^{-\mu \theta} f(\theta)d\theta > {q}-\varepsilon$" in Theorem 5.1 becomes the following simpler form:

Especially, consider when (2.1) degenerates to the following one with bounded time delay:

where $\tau\geq0$ is a constant. Then similar to ^[5,6], for $t\geq0$, let us consider the following differentiable function:

It follows that for $t\geq \tau$,

Thus, it has

Acknowledgments

This paper has been supported by the National Natural Science Foundation of China (No. 11971055) and the Fundamental Research Funds for the Central Universities (No. FRF-BY-17-017).

Conflict of interest

The authors declare there is no conflict of interest.