Stability analysis of a chemostat model with maintenance energy

Accepted Manuscript Stability analysis of a chemostat model with maintenance energy

Tonghua Zhang, Tongqian Zhang, Xinzhu Meng

PII: DOI: Reference:

S0893-9659(16)30358-5 http://dx.doi.org/10.1016/j.aml.2016.12.007 AML 5144

To appear in:

Applied Mathematics Letters

Received date : 27 October 2016 Revised date : 4 December 2016 Accepted date : 5 December 2016 Please cite this article as: T. Zhang, et al., Stability analysis of a chemostat model with maintenance energy, Appl. Math. Lett. (2016), http://dx.doi.org/10.1016/j.aml.2016.12.007 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

*Manuscript Click here to view linked References

Stability analysis of a chemostat model with maintenance energy Tonghua Zhanga , Tongqian Zhangb,c , Xinzhu Mengb,c,∗ a Department of Mathematics, Swinburne University of Technology, Melbourne, VIC 3122 Australia of Mathematics and Systems Science, Shandong University of Science and Technology, Qingdao 266590, China c State Key Laboratory of Mining Disaster Prevention and Control Co-founded by Shandong Province and the Ministry of Science and Technology, Shandong University of Science and Technology, Qingdao 266590, PR China b College

Abstract In this paper, we dedicate ourselves to the study of a diffusive model for unstirred membrane reactors with maintenance energy and subject to a homogeneous Neumann boundary condition. It shows that the unique constant steady state is globally asymptotically stable when it exists. This result further implies the non-existence of any spatial patterns. Keywords: Diffusion, Global stability, Local stability, Membrane reactor 2010 MSC: 35B35, 92D25

1. Introduction Motivated by works in [1–4], in this letter we consider a spatially generalised version of chemostat models. By removing the condition of “well stirred” and including a term of “maintenance energy” as in [5, 6], we reach a model in the form of ( µ(S) ∂S V ∂t 0 = DS ∆S + F (S0 − S) − α V X − V mS X, (1) ∂X V ∂t0 = DX ∆X + βF (X0 − X) + V Xµ(S) − V kd X, where mS is known as maintenance coefficient representing the energy used for functions other than cell growth and such energy is generated by consumption of the substrate[5], µ(S) is known as the response function in [2, 5, 7] and DS and DX are diffusive coefficients. Depending on the types of nutrient and microorganisms, µ(S) may take different forms, for example authors used Monod growth rate in [5–10]. Following the suggestion of [11], we consider a growth rate in the form of S2 µ(S) = Kµ2m+S 2 . For the convenience of mathematical analysis, we perform linear transformations S

S = k1 u, X = k2 v, t0 = k3 t, τ = with k1 = KS , k2 = αKS , k3 = d1 =

1 µm ,

V F

and let

k3 D S k3 DX ∗ τ αmS ∗ S0 X0 , d2 = , τ = , m∗S = , kd = kd k3 , s0 = , X∗ = . V V k3 µm k1 0 k2

Then after dropping the asterisks for notational simplicity we reach a nondimensional model ( 1 u2 v ∂u ∂t = d1 ∆u + τ (s0 − u) − 1+u2 − mS v, β ∂v u2 v ∂t = d2 ∆v + τ (X0 − v) + 1+u2 − kd v,

(2)

where all parameters are positive and 0 ≤ β ≤ 1, and τ1 is known as the residence time [5, 12]. In what follows, we assume, without loss of generality, the above model is subject to the following initial condition and homogeneous Neumann boundary condition [1] u(x, 0) = u0 (x) > 0, v(x, 0) = v0 (x) > 0, x ∈ Ω, ∂u ∂v = = 0, x ∈ ∂Ω. ∂ν ∂ν

(3) (4)

∗ Corresponding

author Email addresses: [email protected] (Tonghua Zhang), [email protected] (Xinzhu Meng)

Preprint submitted to Applied Mathematics Letters

December 5, 2016

In the rest of this paper, we aim to investigate the local and global stabilities of the constant steady state E ∗ of model (2) subject to initial condition (3) and boundary condition (4). We first discuss the existence of the uniform steady states of the system. Then Section 3 devotes to the study of ideal membrane reactor. More precisely, Section 3.1 dedicates to the study of local stability, and the global stability of E ∗ will be discussed in Section 3.2, where we first show the existence of a positively invariant set and then demonstrate the global stability by employing a Lyapunov function. Finally, at the end of Section 3, we extend our results for ideal membrane reactor model to some non-ideal membrane models with sterile reactors. 2. Uniform steady states of (2) In this case we have β = 0. Then it is easy to verify that system consisting of (2)-(4) with β = 0: (a) always has an equilibria E0 (s0 , 0), which is known as a washout equilibrium; and (b) a positive equilibrium E ∗ (u∗ , v ∗ ) q s2

∗ 0 when kd < 1+s 2 where u = 0 Jacobian at the equilibrium is

kd 1−kd

J≡

and v ∗ =



J11 J21

J12 J22

s0 −u∗ τ (kd +mS ) .



=

When without diffusive effect, namely di = 0, i = 1, 2 the

2uv (1+u2 )2 2uv (1+u2 )2

− τ1 −

2

u − 1+u 2 − mS u2 1+u2

− kd

!

.

Straightforward calculation shows at the washout equilibrium point (s0 , 0) the two eigenvalues are s20 1 − kd λ1 = − , λ2 = τ 1 + s20 which imply that when the positive equilibrium exists the washout equilibrium is a saddle point, otherwise it is a stable node. Since, in practice we are only interested in positive equilibrium in what follows we always assume s20 ∗ is stable. kd < 1+s 2 . And we know in this case the positive equilibrium E 0

3. Stability analysis of the positive equilibrium 3.1. Local stability analysis of E ∗ Let L=



d1 ∆ + J11 J21

 J12 . d2 ∆ + J22

(5)

Then we can prove the following theorem for the initial-boundary problem (2)-(4). s2

0 Theorem 3.1. When 0 < kd < 1+s 2 system governed by equations (2)-(4) with β = 0 has a unique positive 0 equilibrium E ∗ . And when it exists, E ∗ is uniformly asymptotically stable in the sense of [13].

Proof. The existence of the uniform steady state has been discussed in previous section. Next, we prove the local stability by verifying that all eigenvalues of the linear operator associated with (2) have negative real part. For the convenience of discussion, we first revisit some notations used in [14]. Assume λi , i = 0, 1, 2 · · · are eigenvalues P ∂2 of −∆ on its domain Ω with Neumann boundary condition, where ∆ = j ∂x 2 is the Laplacian operator and j

λi+1 > λi > λ0 = 0. Furthermore, we use E(λi ) to denote the eigenspace associated with λi and to represent the orthonormal basis of E(λi ) by Xi . Then the solution space, X = {(u, v)} of (2)-(4) can be decomposed as the direct L∞ sum X = i=0 Xi . By Peng and Wang [14], it is easy to see that Xi is an invariant set under the linear operator L defined in (2) and eigenvalues of L on Xi are equivalent to that of matrix   −d1 λi + J11 J12 Mi = . J21 −d2 λi + J22 Since at the positive equilibrium E ∗ we have J11 < 0, J12 < 0, J21 > 0 and J22 = 0, the determinant of Mi satisfies −d λ + J11 J12 = d1 d2 λ2i − d2 J11 λi − J12 J21 > −J12 J21 > 0 det Mi = 1 i J21 −d2 λi + J22

and the trace tr Mi = −(d1 + d2 )λi + J11 < J11 < 0 for all i = 0, 1, 2, · · · , respectively. Hence, we know that E1 is uniformly asymptotically stable. Moreover, the stability of E ∗ is local in the sense of Henry [13]. 2

Remark 3.1. Theorem 3.1 shows the local stability of the positive equilibrium of system (2)-(4), which implies the system does not have non-constant positive steady state in a neighbourhood of E ∗ . It physically means the local non-existence of spatial patterns. s2

0 Notice that at the washout equilibrium point E0 (s0 , 0), we have J11 = − τ1 , J12 = − 1+s 2 −mS , J21 = 0 and J22 =

s20 1+s20

0

− kd . Therefore, when E ∗ exists, the Jacobian at E0 has two eigenvalues

η1i = −d1 λi + J11 < 0 and η2i = −d2 λi + J22 on each Xi , which implies η20 = J22 > 0 in this case. Then using the similar argument we can prove the following conclusion. Proposition 3.1. For the washout equilibrium E0 of system (2), it is asymptotically stable if kd > kd <

s20

1+s20

s20 ; 1+s20

while

it is unstable.

Remark 3.2. This proposition says that when the washout equilibrium is the only uniform steady state, it is uniformly, asymptotically stable; otherwise, if there is a positively uniform steady state, then the washout one becomes unstable. 3.2. Global stability of the positive equilibrium E ∗ In previous section, we have mainly discussed the local stability of E ∗ . This section devotes to the proof of the global stability. We start with proving the following lemma. Lemma 3.1. For system defined by equations (2)-(4) with β = 0, (i) when mS + kd ≥ τ1 and τ mS < 1, there is a constant D satisfying s0 < D ≤ τ sm0S such that it has a positively invariant set Γ, which is the region enclosed by the positive axes, line u = s0 and line l defined below; (ii) when mS +kd < τ1 , there is a constant D satisfying τ (mSs0+kd ) < D ≤ τ sm0S such that it has a positively invariant set Γ, which is the region enclosed by the axes, line u = s0 and line l defined below. Proof. Rewrite the first equation of (2) as ∂u 1 u2 v − mS v. − d1 ∆u = (s0 − u) − ∂t τ 1 + u2 Then we can easily see that u(·, t) < s0 uniformly for t large enough. Similarly, the second equation of (2) u2 v ∂v − d2 ∆v = − kd v ∂t 1 + u2 shows that v(·, t) remains nonnegative. Define line l := u + v − D = 0, where D > 0 is a constant to be determined later. By selecting D carefully, we next show that all solutions of system (2)-(4) will remain below l if the initial value lies below l too. In fact, along the trajectory defined by system (2) we have dl ∂u ∂v 1 1 = + = (s0 − u) − (mS + kd )v = (s0 − u) − (mS + kd )(D − u) dt ∂t ∂t τ  τ s0 1 = − (mS + kd )D − − (mS + kd ) u < 0 τ τ if s0 − τ (mS + kd )D < (1 − τ (mS + kd ))u. That for all 0 < u(·, t) < s0 we have one of the following: (a) if (1 − τ (mS + kd )) ≤ 0 then D > s0 ; or (b) if (1 − τ (mS + kd )) > 0 then D >

s0 τ (mS +kd )

> s0 .

3

Then (a) and (b) yield the conditions in (i) and (ii), respectively. Furthermore, we can verify that the positive q ∗ kd ∗ ∗ equilibrium E is in Γ by carefully selecting D. In fact, notice that u = 1−kd and v ∗ = τ (ks0d−u +mS ) and for case (i) we have τ (mS ) + kd ) ≥ 1 and τ mS < 1. Then for all Ds satisfying s0 < D < τ sm0S we obtain u∗ + v ∗ − D =

s0 − u∗ + u∗ − D < s0 − u∗ + u∗ − D = s0 − D < 0 τ (kd + mS )

which implies E ∗ is in Γ; Similarly, for case (ii), noticing τ (mS + kd ) < 1 yields u∗ + v ∗ − D = if

s0 τ (mS +kd )

≤D≤

s0 τ mS .

s0 − u∗ s0 s0 + u∗ − D < − u∗ + u∗ − D = −D ≤0 τ (kd + mS ) τ (mS + kd ) τ (mS + kd )

Again we just showed that E ∗ locates in Γ. This completes the proof.

Remark 3.3. Using the invariant sets theory developed by Weinberger [15], in Lemma 3.1 we proved the global existence and uniqueness of solutions of system governed by equations (2)-(4). We next prove the global stability of the positive equilibrium. Theorem 3.2. The positive equilibrium, E ∗ of system (2)-(4) is globally asymptotically stable when kd < β = 0.

s20 1+s20

and

Proof. The existences of the solution and the positive equilibrium of the system have been discussed above. Hence, here we mainly focus on the construction of the Lyapunov function and the proof of the global stability. To this end, denote the solution of (2) with positive initial value (3) and boundary condition (4) by (u(x, t), v(x, t)) . Inspired by the work of Hsu [16] and of Hattaf and Yousfi [17], we construct a Lyapunov function as follows. Define ∗ (u)−k) u2 Q(u) = τ v (f where f (u) = 1+u 2 + mS and k = mS + kd , and let s0 −u W (u, v) = R

Z

u

Q(ξ)dξ +

u∗

Z

v

v∗

η − v∗ dη. η

Then the Lyapunov function is E(t) = Ω W dx. Since for differentiable function h(u) and u satisfying the Neumann boundary condition on ∂Ω we have Z Z Z Z Z Z ∂u h(u)∆udx = h(u)∇2 udx = − ∇h(u) · ∇udx + h(u) =− ∇h(u) · ∇udx = − h0 (u)|∇u|2 dx, ∂n Ω Ω Ω ∂Ω Ω Ω the straightforward calculation along the flow generated by (2) yields Z dE(t) = (Wu ut + Wv vt )dx dt Ω Z    v − v∗   1 u2 v u2 v = Q(u) d1 ∆u + (s0 − u) − − m v + d ∆v + − k v dx 2 S d τ 1 + u2 v 1 + u2 Ω Z   d2 v ∗ = − d1 Q0 (u)|∇u|2 + 2 |∇v|2 dx v Ω Z  1  v − v ∗  u2 v  u2 v + Q(u) (s0 − u) − − mS v + − kd v dx. τ 1 + u2 v 1 + u2 Ω Next, we show that

dE(t) dt

(7)

< 0, which together with Lemma 3.1 implies the globally asymptotical stability of

E ∗ . Since Q(u) can be written as Q(u) = Q0 (u) =

(6)

2

u τ v ∗ ( 1+u 2 −kd )

s0 −u

, the derivative of Q with respective to u is

τ v ∗ (2u(s0 − u) + Q1 (u)) , with Q1 (u) = (1 − kd )u4 + (1 − 2kd )u2 − kd . (1 + u2 )2 (s0 − u)2

Obviously, Q1 is a quadratic polynomial in terms of u2 , with the coefficient of the leading term satisfying 1−kd > 0. 2kd −1 Then Q1 (u), at u = 2(1−k < u∗ , has a minimal value, which implies that Q1 (u) > Q1 (u∗ ) = 0 for all u > u∗ . d) Hence, Q0 (u) > 0 4

which implies the integral over Ω in (6) is strictly less than zero when di 6= 0. Furthermore, we claim that  v − v ∗  u2 v  τ v f (u) − k)  1 1  u2 v ∗ − m v + − k v = f1 = Q(u) (s0 − u) − (s − u) − v f (u) ≤ 0. (8) S d 0 τ 1 + u2 v 1 + u2 (s0 − u) τ If this claim is not true, then we have two subcases ( f (u) − k > 0, 1 ∗ τ (s0 − u) − v f (u) > 0, or

(

(9)

f (u) − k < 0, 1 ∗ τ (s0 − u) − v f (u) < 0,

(10) 2

u since u < s0 and v > 0. In what follows, we prove the case of (9) can not happen. Notice that f (u) = 1+u 2 + mS is q kd ∗ ∗ increasing about u and f (u ) − k = 0. Then f (u) − k > 0 implies that u > u = 1−kd . From the second equation ∗

s0 −u s0 −u ∗ of (9), we have v ∗ < τsf0 −u (u) < τ f (u∗ ) < τ (kd +mS ) = v . This contradiction implies that case (9) can not happen. We then show (10) can not happen either. Otherwise, from the first equation we have f (u) < k = f (u∗ ), which implies s0 −u∗ ∗ > that 0 < u < u∗ . From the second equation, we have v ∗ > τsf0 −u (u) τ (kd +mS ) = v . Again, this is a contradiction implying that (10) is not true. Hence, f1 ≤ 0. Therefore the integral in equations (7) is nonpositive. Then from ∗ ∗ the above analysis, we have proven that dE(t) dt < 0 which implies that (u , v ) is globally asymptotically stable.

Remark 3.4. Theorem 3.2 shows the global non-existence of non-constant positive solution, namely globally system (2)-(4) has no spatial patterns. For the non-ideal membrane reactor or continuous flow reactor model, namely for the case of 0 < β < 1, by using the similar argument, we can easily get the following results when the reactor is sterile, which implies X0 = 0. s2

∗ 0 Theorem 3.3. When 0 < kd + βτ < 1+s of the initial-boundary problem (2)-(4) 2 , the unique positive equilibrium E 0 is uniformly asymptotically stable in the sense of [13].

Lemma 3.2. For initial-boundary problem (2)-(4), and τ mS < 1, there is a constant D satisfying s0 < D ≤ τ sm0S such that it has a positively (i) when mS + kd ≥ 1−β τ invariant set Γ, which is the region enclosed by the positive axes, line u = s0 and line l defined below; s0 s0 (ii) when mS + kd < 1−β τ , there is a constant D satisfying τ (mS +kd )+β < D ≤ τ mS such that it has a positively invariant set Γ, which is the region enclosed by the axes, line u = s0 and line l defined below.

Theorem 3.4. The positive equilibrium, E ∗ of system (2)-(4) is globally asymptotically stable when kd + βτ <

s20 . 1+s20

Remark 3.5. By comparing results for the ideal membrane reactor with the non-ideal, sterile reactor, we can easily to conclude that the ideal reactor can sustain larger washout rate, kd in order to guarantee the stability of the positive equilibrium. Acknowledgement This work is supported by the National Natural Science Foundation of China (11371230, 11501331), the SDUST Research Fund (2014TDJH102), Shandong Provincial Natural Science Foundation, China (ZR2015AQ001, BS2015SF002), Joint Innovative Center for Safe And Effective Mining Technology and Equipment of Coal Resources, Shandong Province, a Project of Shandong Province Higher Educational Science and Technology Program of China (J13LI05).

5

References [1] J. H. Merkin, V. Petrov, S. K. Scott, K. Showalter, Wave-induced chaos in a continuously fed unstirred reactor, Journal of the Chemical Society, Faraday Transactions 92 (16) (1996) 2911–2918. [2] T. Zhang, H. Zang, Delay-induced turing instability in reaction-diffusion equations, Physical Review E 90 (2014) 052908. [3] S.-B. Hsu, J. Jiang, F.-B. Wang, On a system of reaction-diffusion equations arising from competition with internal storage in an unstirred chemostat, Journal of Differential Equations 248 (2010) 2470–2496. [4] S.-B. Hsu, J. Shi, F.-B. Wang, Further studies of a reaction-diffusion system for an unstirred chemostat with internal storage, Discrete and Continuous Dynamical Systems B 19 (10) (2014) 3169–3189. [5] M. I. Nelson, T. B. Kerr, X. Chen, A fundamental analysis of continuous flow bioreactor and membrane reactor models with death and maintenance included, Asia-Pacific Journal of Chemical Engineering 3 (2008) 70–80. [6] T. Zhang, Global analysis of continuous flow bioreactor and membrane reactor models with death and maintenance, Journal of Mathematical Chemistry 50 (2012) 2239–2247. [7] H. L. Smith, P. Waltman, The theory of the chemostat: Dynamics of microbial competition, Cambridge University Press, Cambridge, 2008. [8] Z. Chen, T. Zhang, Long time behaviour of a stochastic model for continuous flow bioreactor, Journal of Mathematical Chemistry 51 (2013) 451–464. [9] T. Zhang, Y. Xing, H. Zang, M. Han, Spatio-temporal dynamics of a reaction-diffusion system for a predator– prey model with hyperbolic mortality, Nonlinear Dynamics 78 (1) (2014) 265–277. [10] T. Zhang, Y. Song, Global stability of the positive equilibrium of a mathematical model for unstirred membrane reactors, Bulletin of the Korean Mathematical Society, Accepted for publication on 4 May 2016. [11] J. D. Murray, Mathematical Biology I: An Introduction, 3rd Edition, Springer-Verlag, 2002. [12] M. L. Shuler, F. Kargi, Bioprocess Engineering: Basic Concepts, 2nd Edition, Prentice hall International Series in the Physical and Chemical Engineering Sciences, Prentice Hall PTR, Upper Saddle River, NJ 07458, 2002. [13] D. Henry, Geometric theory of semilinear parabolic equations, in: Lecture Notes in Mathematics, 3rd Edition, Vol. 840, Springer-Verlag, Berlin, New York, 1993. [14] R. Peng, M. Wang, Global stability of the equilibrium of a diffusive holling-tanner prey-predator model, Applied Mathematics Letters 20 (2007) 664–670. [15] H. F. Weinberger, Invariant sets for weakly coupled parabolic and elliptic systems, Rendiconti di Matematica 8 (1) (1975) 295–310. [16] S.-B. Hsu, A survey of constructing lyapunov functions for mathematical models in population biology, Taiwanese Journal of Mathematics 9 (2) (2005) 151–173. [17] K. Hattaf, N. Yousfi, Global stability for reaction-diffusion equations in biology, Computers and Mathematics with Applications 66 (8) (2013) 1488–1497.

6