Proof: The existence and uniqueness of the solution of the system (12) are given by the Cauchy-Lipschitz theorem walter , since all functions involved in this system are globally Lipschitzian.

In order to show the positivity of the solution of the system (12), let $X_{0}\geq 0$, $,B_{0}\geq 0$, $s_{0}\geq 0$ and $P_{0}\geq 0$ be given. We proceed by step, and show, first, that $B(t)\geq 0$, then we deduce that $X(t)\geq 0$, and then finally the positivity result for $s(t)$ and $P(t)$.

1. ⁢1/

   If $B_{0}=0$, then the second equation of the system (12) gives $B(t)=0$ for all $t>0.$ If $B_{0}>0$ and $X_{0}$ and $s_{0}$ in $\hbox{${\mathbb{R}}$}_{+}$, by the uniqueness theorem of the solutions of the Cauchy problem walter , we show that $B(\cdot)$ cannot take negative values. Indeed, we suppose there exists $t^{{}^{\prime}}>0$ such that $B(t^{{}^{\prime}})=0$. We know (see walter ) that there exists one and only one value of the initial condition $B_{0}$ for which the solution satisfies $B(t^{{}^{\prime}})=0$. But we have already shown that for the initial condition $B_{0}=0$, the equality $B(t^{{}^{\prime}})=0$ is verified. The solution being unique, this is a contradiction with $B_{0}>0$. Therefore, for any vector $(X_{0},B_{0},s_{0},P_{0})$, with non-negative components, we have $B(t)\geq 0$ for all $t>0.$

2. 2/

   Using the positivity of $B(\cdot)$, any solution of the first equation, with $(X_{0},B_{0},s_{0},P_{0})$ in $\hbox{${\mathbb{R}}$}_{+}^{3}$ satisfies $\displaystyle\frac{dX(t)}{dt}\geq-K_{H}X(t)$ for all $t>0$. We deduce, by integration, that $X(t)\geq X_{0}e^{-K_{H}t}$, for all $t>0$. Therefore, for any vector $(X_{0},B_{0},s_{0},P_{0})$, with non-negative components, we have $X(t)\geq 0$, for all $t>0$.

3. 3/

   From the third equation of (12) and by the positivity of $X(t)$ we have

   $\displaystyle\displaystyle\frac{ds(t)}{dt}\geq-\frac{1}{Y_{B/s}}\>\>(\mu(s(t))+m_{s})\>\>B(t).\>\>\>\>\>\>\forall t>0.$

   (16)

   Let us now consider $S$ the solution of the differential equation

   $\displaystyle\left\{\begin{array}[]{rcll}\displaystyle\frac{dS}{dt}&=&-\displaystyle\frac{1}{Y_{B/s}}\>\>(\mu(s)+m_{s})\>\>B.\\ S(0)&=&s_{0}.\end{array}\right.$

   (19)

   We prove, by the same argument as before, that the solution $S(t)$ of (19) cannot cross the plane $S=0$, nor reach 0 in finite time. By the comparison theorem LAKSH 1969 , since $s(t)\geq S(t)$, for all $t>0$, $s(t)$ is also non-negative for all $t>0$, and cannot reach 0 in finite time, whatever $(X_{0},B_{0},s_{0},P_{0})$ in $\hbox{${\mathbb{R}}$}_{+}^{3}$.

4. 4/

   By Hypothesis 1, Hypothesis 2 and the positivity of $B(\cdot)$, all the coefficients of the last equation of (12) are non-negative. We deduce that $\displaystyle\frac{dP}{dt}$ is non-negative and so $P(t)$ is increasing. Therefore according to the fact that $P(0)$ is in $\hbox{${\mathbb{R}}$}_{+}$, $P(t)$ is also non-negative for all $t>0$.

An equilibrium point of the non-linear system (12), $E=(X^{*},B^{*},s^{*},P^{*})$, is a solution of the system

$$A(V)V=0,$$

(20)

where $A$ and $V$ are defined in (14).

From the second equation of the system (20), we derive that either $B=0$ or $\mu(s)=k_{d}$.

• ✓

  If $B=0$ then by the first equation, we have $X=0$ whatever $s\in\hbox{${\mathbb{R}}$}$, and the last two equations are automatically verified for any $s,P\in\hbox{${\mathbb{R}}$}$. The vector $E=(0,0,s,P)$, for $s,P\in\hbox{${\mathbb{R}}$}$ is then an equilibrium point.

• ✓

  If $B\neq 0$, the the fourth equation of (20) gives $\displaystyle\frac{k_{d}}{Y_{p/s}}=-m_{p},$ which is not possible because $m_{p}>0$ and $\displaystyle\frac{k_{d}}{Y_{p/s}}\geq 0$.

So, we have a continuum of equilibrium points given by $E=(0,0,s,P)$, for $s$ and $P$ in ${\mathbb{R}}$, whatever the function $\mu$.

Usually, using the linearisation method, the stability of equilibria of ODEs is determined, by the sign of real part of eigenvalues of the Jacobian matrix. For a point $E=(0,0,s^{*},P^{*})$ these eigenvalues are given by

$\displaystyle\lambda_{1}=-K_{H},\qquad\lambda_{2}=\mu(s^{*})-k_{d},\qquad\lambda_{3}=0,\qquad\lambda_{4}=0.$

(21)

We cannot conclude on the stability or instability, because these equilibrium points are not hyperbolic (the eigenvalues $\lambda_{3}$ and $\lambda_{4}$ are zero). In such situations, one could use the Lyapunov function (see walter , page 319), but it is generally difficult to find such a function, especially for a complex non-linear system. We show in the next section, using Barbalat’s lemma barbalat , that the solutions of the system tend to an equilibrium point when $t$ tends to $+\infty$, which makes it a globally and asymptotically stable equilibrium.

Note that $\lambda_{2}\leq 0$ if and only if $\mu(s^{*})\leq K_{d}$. In the study of basins of attraction, we will introduce the set of the values of $s$ such $\mu(s)\leq K_{d}$ defined by

$\displaystyle\mathcal{S}:=\{s\in\hbox{${\mathbb{R}}$}_{+}:\mu(s)\leq k_{d}\}.$

(22)

Under Hypothesis 1 and Hypothesis 2, this set has a non-empty interior. For the Monod law, we have

$\displaystyle\mathcal{S}=[0,\lambda]\>\>\>\>\>\>\mbox{with}\>\>\>\>\>\>\lambda=\frac{k_{d}k_{s}}{\mu_{max}-k_{d}}.$

(23)

In this section, we describe the asymptotic behaviour of the solution of the system (12) and deduce the stability of the equilibrium points. We begin by showing that when time goes to infinity, the solution tends to a point of the form $(0,0,s^{*},P^{*})$ for some values $s^{*}$ and $P^{*}$ which depend on the parameters of the problem and the initial conditions.

Proof: Consider the function $Z(t)$ defined by

$\displaystyle Z(t)=\left(1+\frac{m_{s}Y_{B/s}}{k_{d}}\right)X(t)+\alpha B(t)+\alpha Y_{B/s}s(t).$

(26)

This function satisfies some properties that are essential for the rest of the proof. We start by mentioning them before evaluating the limits at infinity of the components of the solution of the system (12).

• •

  $Z(t)\geq 0$, for all $t\geq 0$.

• •

  The derivative of $Z(t)$ is given by

  $$\frac{dZ}{dt}=K_{H}\left[\left(\alpha Y_{B/s}-1\right)-\frac{m_{s}Y_{B/s}}{k_{d}}\right]X,$$

  (27)

  which makes it possible to write $Z(t)$ also in the form

  $$Z(t)=Z(0)+K_{H}\left[\left(\alpha Y_{B/s}-1\right)-\frac{m_{s}Y_{B/s}}{k_{d}}\right]\int^{t}_{0}X(\tau)d\tau$$

  (28)

• •

  There exists $0\leq Z\infty<+\infty$ such that

  $$\lim_{t\to+\infty}Z(t):=Z\infty,$$

  (29)

  indeed $Z(t)$ is decreasing (by Hypothesis 2, $\left(\alpha Y_{B/s}-1\right)-\dfrac{m_{s}Y_{B/s}}{k_{d}})<0$) and minorized by zezo.

• •

  The second derivative of $Z(t)$ is given by

  $$\frac{d^{2}Z}{dt^{2}}=K_{H}\left((\alpha Y_{B/s}-1)-\frac{m_{s}Y_{B/s}}{k_{d}}\right)\left(-K_{H}X+\alpha k_{d}B\right).$$

  (30)

1. 1.

   We begin by showing that the solution of (12) converge asymptotically to an equilibrium point.
   

   Limits of $X(t)$ and $B(t)$
   With respect to the definition of $Z(t)$, the limit (29) implies that the variables $X(t)$, $B(t)$ and $s(t)$ are bounded (since they are non-negative). We deduce that the second derivative (30) is bounded, therefore $\displaystyle\frac{dZ}{dt}$ is uniformly continuous on $\hbox{${\mathbb{R}}$}^{+}$. Barbalat’s Lemma barbalat allows us to assert that $\displaystyle\underset{t\rightarrow+\infty}{\lim}\frac{dZ}{dt}=0$. By (27), we obtain $\displaystyle\lim_{t\to+\infty}X(t)=0$.
   Similarly, $\displaystyle\frac{d^{2}X}{dt^{2}}$, which is expressed as a function of $X(t)$ and $B(t)$, is also bounded, and then $\displaystyle\frac{dX}{dt}$ is uniformly continuous and by Barbalat’s lemma, we deduce $\displaystyle\underset{t\rightarrow+\infty}{\lim}\frac{dX}{dt}=0$, and therefore, by the first equation of the system and knowing that $\displaystyle\lim_{t\to+\infty}X(t)=0$, we have $\displaystyle\underset{t\rightarrow+\infty}{\lim}B(t)=0$.
   

   Limit of $s(t)$
   By (26) and (29) and as we have already shown that $\displaystyle\lim_{t\to+\infty}X(t)=\displaystyle\lim_{t\to+\infty}B(t)=0$, we can affirm that there exists a positive real $s^{\star}$ such that

   $$\lim_{t\to+\infty}s(t)=s^{\star}=\frac{Z_{\infty}}{\alpha Y_{B/s}}.$$

   (31)

   The function $Z$ being decreasing, we have $Z_{\infty}\leq Z(0)$, and then $s^{\star}\leq\displaystyle\frac{Z(0)}{\alpha Y_{B/s}}$.
   As $Z(0)=\left(1+\displaystyle\frac{m_{s}Y_{B/s}}{k_{d}}\right)X_{0}+\alpha B_{0}+\alpha Y_{B/s}s_{0}$, we get (24).
   

   Limit of $P(t)$
   Consider $Z(t)$ written in the form (28). The function $Z(.)$ being bounded, we deduce that it is the same for $\displaystyle\int^{t}_{0}X(\tau)d\tau$, then by using the integration of the system (12) between $0$ and $t$ for $t>0$ gives

   	$\displaystyle X(t)$	$\displaystyle=$	$\displaystyle X(0)-K_{H}\int_{0}^{t}X(\tau)d\tau+\alpha\,k_{d}\int_{0}^{t}B(\tau)d\tau$
   	$\displaystyle B(t)$	$\displaystyle=$	$\displaystyle B(0)+\int_{0}^{t}\mu(s(\tau))B(\tau)d\tau-k_{d}\int_{0}^{t}B(\tau)d\tau$
   	$\displaystyle s(t)$	$\displaystyle=$	$\displaystyle s(0)-\displaystyle\frac{1}{Y_{B/s}}\int_{0}^{t}\mu(s(\tau))B(\tau)d\tau-m_{s}\int_{0}^{t}B(\tau)d\tau+K_{H}\int_{0}^{t}X(\tau)d\tau$
   	$\displaystyle P(t)$	$\displaystyle=$	$\displaystyle P(0)+\int_{0}^{t}\frac{1}{Y_{P/s}}\mu(s(\tau))B(\tau)d\tau+m_{P}\int_{0}^{t}B(\tau)d\tau.$

   We deduce by cascade that $\displaystyle\int^{t}_{0}B(\tau)d\tau<+\infty$, $\displaystyle\int^{t}_{0}\mu(s(\tau))B(\tau)d\tau<+\infty$ and that $\displaystyle\lim_{t\to+\infty}P(t)$ is finite and positive as a limit of a positive and bounded increasing function.

   Taking the limit in the last expressions of $X(t)$, $B(t)$ and $P(t)$, when $t$ tends to $+\infty$, and having $\displaystyle\lim_{t\to+\infty}X(t)=\displaystyle\lim_{t\to+\infty}B(t)=0$ and (24) gives

   	$\displaystyle 0$	$\displaystyle=$	$\displaystyle X_{0}-K_{H}\int^{+\infty}_{0}X(\tau)d\tau+\alpha k_{d}\int^{+\infty}_{0}B(\tau)d\tau$
   	$\displaystyle 0$	$\displaystyle=$	$\displaystyle B_{0}-k_{d}\int^{+\infty}_{0}B(\tau)d\tau+\int^{+\infty}_{0}\mu(s(\tau))B(\tau)d\tau$
   	$\displaystyle s^{*}$	$\displaystyle=$	$\displaystyle s_{0}-\displaystyle\frac{1}{Y_{B/s}}\int_{0}^{+\infty}\mu(s(\tau))B(\tau)d\tau$
   			$\displaystyle\qquad\qquad\qquad-m_{s}\int_{0}^{+\infty}B(\tau)d\tau+K_{H}\int_{0}^{+\infty}X(\tau)d\tau$
   	$\displaystyle\lim_{t\to+\infty}P(t)$	$\displaystyle=$	$\displaystyle P_{0}+\int_{0}^{+\infty}\frac{1}{Y_{P/s}}\mu(s(\tau))B(\tau)d\tau+m_{P}\int_{0}^{+\infty}B(\tau)d.\tau$		(32)

   Linear combinations of the first three equations give

   $$\int^{t}_{0}B(\tau)d\tau=\frac{B_{0}+Y_{B/s}\left[(X_{0}+s_{0}-s^{*})\right]}{k_{d}\left[1-Y_{B/s}(\alpha-\frac{m_{s}}{k_{d}})\right]}$$

   (33)

   $$\int^{t}_{0}\mu(s(t))B(\tau)d\tau=\frac{Y_{B/s}\left[B_{0}(\alpha-\frac{m_{s}}{k_{d}})+(X_{0}+s_{0}-s^{*})\right]}{\left[1-Y_{B/s}(\alpha-\frac{m_{s}}{k_{d}})\right]}$$

   (34)

   By replacing the expressions (33) and (34) in (32), we obtain the expression (25). Note that by Hypothesis 2, $\left[1-Y_{B/s}(\alpha-\dfrac{m_{s}}{k_{d}})\right]$ is strictly positive.
   

2. 2.

   Let us now show that $s^{\star}$ cannot be equal to 0 when $X_{0}$ or $s_{0}$ are not zero. Suppose that $X_{0}\neq 0$ or $s_{0}\neq 0$ and that $s(t)$ tends to 0 when t tends to infinity. Then by continuity, the function $\mu(s(t))$ will tend towards zero too, and there would exist $T$ such that for all $t\geq T$ we have

   $$\mu(s(t))\leq Y_{B/s}\alpha k_{d}-Y_{B/s}m_{s}.$$

   (35)

   The sum of the first and third equations of the system (12) allows us to write

   $$\frac{d}{dt}(X(t)+s(t))=\frac{1}{Y_{B/s}}(Y_{B/s}\alpha k_{d}-Y_{B/s}m_{s}-\mu(s))B.$$

   (36)

   By (35) we conclude that $\dfrac{d}{dt}(X(t)+s(t))\geq 0$ for all$t>T.$ By the same argument as the one used in the demonstration of the positivity of the solution, we conclude that if $X(0)+s(0)>0$, the variable $X(\cdot)$ or $s(\cdot)$ cannot reach 0 in finite time. We then have $X(t)+s(t)\geq X(T)+s(T)>0$ for all $t>T$, which is a contradiction with the fact that $X(t)+s(t)$ tends to $0$ when $t$ tends towards $+\infty$.
   

3. 3.

   We prove now, that $s^{*}$ belongs to $\mathcal{S}$. Let $(X_{0},s_{0},B_{0},P_{0})$ be in $\hbox{${\mathbb{R}}$}_{+}^{4}$ with $X_{0}>0$ and $B_{0}>0$. Suppose that $s^{*}$, defined by (24), does not belong to $\mathcal{S}$. Then there would exist $\tilde{T}>0$ such that for all $t>\tilde{T}$ we would have : $\mu(s(t))-k_{d}>\eta:=\frac{\mu(s^{*})-k_{d}}{2}.$ By (12), we would then have $\frac{dB(t)}{dt}>\eta B(t).$ Therefore we would have $B(t)>B(s)\exp^{\eta(t-s)}$ for all $t,s>\tilde{T},\>\text{such that}\>\>t\geq s,$ and then $B(.)$ cannot converge asymptotically to $0$, which is incompatible with the previous results .

As we have already pointed out, the equilibria of the system (12) are not hyperbolic. We cannot conclude on their stability by using directly the linearization technique and the central variety theorem PERCO . However, using a suitable change of variable we prove the existence of an invariant stable variety in the following theorem.

Proof: The proof of this theorem is based on two changes of variables, the first one is for the variable $s$ and the second for the variable $P$. We obtain a new system associated with the system (12) and whose equilibria are now hyperbolic. To define the change of variable for $s$, we fix $s^{*}>0$ and $P^{*}>0$ such that $\mu(s^{*})<k_{d}$. Consider the solutions of (12) with $B(t)>0$, for $t\in[0,+\infty)$, and let

$$z(t)=\dfrac{X(t)+(s(t)-s^{*})}{B(t)}+\varphi,\quad\text{with}\quad\varphi:=\dfrac{\dfrac{\mu(s^{*})}{Y_{B/s}}-(\alpha k_{d}-m_{s})}{\mu(s^{*})-k_{d}},$$

(37)

which is equivalent to

$$s(t)=s^{*}-X(t)+(z(t)-\varphi)B(t).$$

(38)

The function $z(t)$ satisfies the equation (see Appendix, page A, for more details)

$$\frac{dz}{dt}=-\gamma(\mu(s)-\mu(s^{*}))-(\mu(s)-k_{d})z$$

(39)

with $\displaystyle\gamma:=\left[\frac{k_{d}(1-\alpha Y_{B/s})+Y_{B/s}m_{s}}{Y_{B/s}(k_{d}-\mu(s^{*}))}\right]>0$, and using (38) we can write $\mu(s)$ in terms of the new variable and denote by $F(z,X,B)$ the following expression

$$F(z,X,B):=\mu\biggl{(}s^{*}-X+(z-\varphi)B\biggr{)}$$

(40)

The second change of variable is defined by setting $P^{*}>0$ and $s^{*}>0$ such that $\mu(s^{*})<k_{d}$. For the solutions of (12) with $B(t)>0$ for $t\in[0,+\infty)$, we define

$$W(t)=\frac{P(t)-P^{*}}{B(t)}+\omega,\quad\text{with}\quad\omega:=\dfrac{\dfrac{-\mu(s^{*})}{Y_{P/s}}-m_{P}}{\mu(s^{*})-k_{d}},$$

(41)

which is equivalent to

$$P(t)=P^{*}+(W(t)-\omega)B(t).$$

(42)

The variable $W(t)$ satisfies the equation (see Appendix, page A, for more details)

$$\frac{dW}{dt}=\phi(\mu(s)-\mu(s^{*}))-(\mu(s)-k_{d})W,$$

(43)

with $\displaystyle\phi:=\left[\frac{k_{d}+Y_{P/s}m_{P}}{Y_{P/s}(k_{d}-\mu(s^{*}))}\right]>0$ .
The system (12) becomes with the variables $z$ and $W$, by using the notation (40)

$$\left\{\begin{array}[]{rcll}\displaystyle\frac{dz}{dt}&=&-\gamma\biggl{(}F(z,X,B)-\mu(s^{*})\biggr{)}-\biggl{(}F(z,X,B)-k_{d}\biggr{)}z,\vspace{0.1cm}\\ \displaystyle\frac{dX}{dt}&=&-K_{H}X+\alpha k_{d}B,\vspace{0.1cm}\\ \displaystyle\frac{dB}{dt}&=&\biggl{(}F(z,X,B)-k_{d}\biggr{)}B,\vspace{0.1cm}\\ \displaystyle\frac{dW}{dt}&=&\phi\biggl{(}F(z,X,B)-\mu(s^{*})\biggr{)}-\biggl{(}F(z,X,B)-k_{d}\biggr{)}w,\\ \end{array}\right.$$

(44)

defined in the domain given by

	$\displaystyle\mathcal{D}$	$\displaystyle:=$	$\displaystyle\{(z,X,B,W)\in\hbox{${\mathbb{R}}$}\times\hbox{${\mathbb{R}}$}_{+}\times\hbox{${\mathbb{R}}$}_{+}^{*}\times\hbox{${\mathbb{R}}$}\quad\text{such that}$
			$\displaystyle\qquad s^{*}-X+(z-\varphi)B\geq 0\quad\text{and}\quad P^{*}+(W-\omega)B\geq 0\}.$

By the condition $s^{*}\in int\>\mathcal{S}$, we have $\mu(s^{*})\neq k_{d}$, and we can check that $(0,0,0,0)$ is the only equilibrium of (44) in $\mathcal{D}$. Indeed, if $E=(z^{*},X^{*},B^{*},W^{*})$ is an equilibrium point, from the $3^{th}$ equation, we derive that either $B=0$ or $F(z^{*},X^{*},B^{*})=k_{d}$.

• •

  If $F(z^{*},X^{*},B^{*})=k_{d}$, from the $1^{th}$ and the $4^{th}$ equations of the system (44) we deduce that $\mu(s^{*})=k_{d}$ which contradicts the condition $s^{*}\in int\>\mathcal{S}$.

• •

  If $B=0$, by the $2^{th}$ equation, we have $X=0$ whatever $s^{*}\in\hbox{${\mathbb{R}}$}_{+}$. In this case, from the first equation we deduce that we necessarily have $z=0$ since we cannot have $\mu(s^{*})=k_{d}$.

• •

  The same reasoning done for the $4^{th}$ equation gives $W=0$.

Thus the vector $E=(0,0,0,0)$ is the only equilibrium point of (44) in $\mathcal{D}$ and the Jacobian matrix associated is given by

$\displaystyle\mathfrak{J}_{(0,0,0,0)}=\begin{pmatrix}-(\mu(s^{*})-k_{d})&\gamma\mu^{{}^{\prime}}(s^{*})&\gamma\varphi\mu^{{}^{\prime}}(s^{*})&0\\ 0&-K_{H}&\alpha k_{d}&0\\ 0&0&\mu(s^{*})-k_{d}&0\\ 0&-\phi\mu^{{}^{\prime}}(s^{*})&-\phi\varphi\mu^{{}^{\prime}}(s^{*})&-(\mu(s^{*})-k_{d})\end{pmatrix}$

(45)

It admits a double eigenvalue $r_{1}$ and two single eigenvalues $r_{2}$ and $r_{3}$:

$$r_{1}=-(\mu(s^{*})-k_{d}),\qquad r_{2}=-K_{H},\qquad r_{3}=\mu(s^{*})-k_{d},$$

(46)

all of them are real and non-zero. According to the Central Variety Theorem Wiggins page 35) and under the condition $\mu(s^{*})<K_{d}$, the point $(0,0,0,0)$ is thus an hyperbolic equilibrium point with a stable two-dimensional variety $\mathcal{L}$ and an unstable two-dimensional variety $\mathcal{U}$ which are respectively positive and negative invariants. So, any trajectory $(z(.),X(.),B(.),W(.))$ in $\mathcal{D}$ of (44) in the stable two-dimensional invariant variety $\mathcal{L}\cap\mathcal{D}$ converges asymptotically to $(0,0,0,0)$.