Alaaeddine HAMMOUDI Oana IOSIFESCU

(Dedicated to Philippe G.Ciarlet on the occasion of his 80th birthday)

1 Introduction

Chemotaxis is the ability of some bacteria to direct their movement according to the gradient of chemicals contained in their environment.In soil,some bacteria microorganisms that degrade organic carbon(SOC for short)are motile and chemotactic.This phenomenon is observed in experiments(see[1])and onfield.Nevertheless to our best knowledge no model of terrestrial carbon cycle adresses this issue.Indeed,these models are essentially compartimental corresponding naturally to systems of ordinary differential equations(e.g.Century,RothC,MOMOS)(see[2]).They are used globally to estimate soil CO2emissions in local land management and crop optimization,among other things.

Very few prototypes of spatial soil organic model have been proposed.Some of them use systems of partial differential equations:Balesdent et al.[3]combined vertical directed transport of organic carbon with a degradation phenomenon and diffusion.More recently,Goudjo et al.[4]proposed a three dimensional model for dissolved organic matter using also a system of PDEs.Other authors opted for afinite sequence of interconnected systems of ordinary differential equations each localized in a soil layer(see[5]).

We previously studied the model MOMOS proposed by Pansu[6—7],which is a nonlinear system of ordinary differential equations(see[8])written as

where

and

In these equations the unknown u models the alive microbial biomass,whereas the unknowns v and w are soil organic matters with distinct decomposition rates.

In reality,the nonnegative functions ki,i∈{1,2,3,4},q and f depend not only on time but also on space because of the variability in soil clay content.The phenomena described by MOMOS can also be subjected to the influence of transport and sedimentation through transport and diffusion.

In order to test the effect of soil heterogeneity we studied in[9]the following reactiondiffusion-advection initial problem:

where Ω is a domain in Rnrepresenting the soil,Aiis a diffusion matrix and Bia transport vector,for each i=1,2,3.

In[9]the boundary conditions were either of Dirichlet type(γ =0,βi≡ 1)or of Neumann-Robin type(γ=1,βi(t,x)≥ 0).The right hand side term of(1.1)was

where we replaced the term q(t,x)with q(t,x)|u1|u1for more accuracy,since q(t,x)u1corresponds to a kinetic coefficient that cannot be negative.We assumed there that the diffusion matrices Aiwere bounded,symmetric and coercive:

and the transport vectors Biwere bounded on QT:

Also we assumed that the functions kj,βiand q were nonnegative and bounded,i.e.,for all j=1,2,3,4 and i=1,2,3,

where constant Cmax＞0.Finally we assumed that the initial data and input were nonnegative and bounded:

In[9]we provedfirst that this model had a unique weak solution.We were looking for weak solutions,because initial inputs were not regular enough to give rise to more“regular” solutions.Second,for periodic data,we proved the existence of a maximal and a minimal periodic solution of this system.In some particular cases,the minimal and the maximal periodic solutions coincide and this function becomes a global attractor for any bounded solution of the periodic system.

In the present work a new PDEs model is considered to take account of chemotaxis.The chemotactic movement of bacteria to root exudates is well known to play an important role in rhizosphere colonisation.Field studies with tracers and laboratory experiments using soil columns were both used to demonstrate the effect of chemotaxis on microbial movements.So,the model proposed here can represent the spatial heterogeneity of soil microbial biomass,highlighted by recent observations at submicron scale(see[1]).

The new model derived from a simplified MOMOS ODEs model,which comprised only two differential equations instead of the three originally,where the microbial biomass was u and the organic matter was v.As additional simplifying hypothesis soil temperature,soil moisture,soil texture and organic input were considered to be isotropic and constant with time.Hence,the simplified ODEs model can be expressed as

with the initial conditions(u0,v0),where k1is the microbial mortality rate,k2is the soil carbon degradation rate,q is the metabolic quotient and f is the soil organic carbon input.It can be proved that the unique positive steady state()is stable(see[8]).

The chemotaxis-type model wasfinded following the conventional Keller-Segel approach(see[10]),using an advection-diffusion system.This comprised two parabolic equations in a smooth domain with no-flux boundary conditions.The advection term was controlled by the gradient of the chemo-attractant.Applying the same principles to our problem leads to the following reaction-diffusion-chemotaxis system(Ph):

The parameter β is the chemotaxis sensitivity,a and d are the diffusion coefficients of microbes and soil organic carbon respectively,Ω is a smooth and bounded domain,and h(·)is a continuous function,involved in the modelling of chemotaxis.As bacteria can release exoenzymes to avoid overcrowding,the function h can be selected to limit overcrowding,as required.This new model is,therefore,a new variation of the Keller-Segel approach(see[10])with the reaction part modified to fit the MOMOS model.In the first equation of(Ph),we change again the term qu2by q|u|u(see[9]).

We prove here(see Appendix 1)the existence of Turing patterns that may provide possible explanations for the formation of soil aggregations,for the bacterial(see[21])and microorganisms spatial organizations(hotspots in soil)or justify the formation of the microscopic patterns observed by Vogel et al.[1].Although spatial heterogeneity can be verified visually in a numerical simulation(see Appendix 2),formal mathematical analysis is required to confirm its emergence and to provide a mathematical proof of the necessary conditions.The mathematical criteria are based on matrices derived from equations and analysed using conditions on the determinant,trace and eigenvalues.

Keller-Segel model was the earliest mathematical system involving chemotaxis(see[10]).Many other models emerged specially in biology and ecology.Most authors focused their efforts essentially on existence and on asymptotic behaviour of solutions in one or two dimensional domains in order to avoid blow-up of solutions(see[11—15]and references therein).

Unlike the classical Keller-Segel model,where equations are coupled only by the chemotactic term,the system of partial differential equations(Ph)is also coupled through the reaction term.More specifically,the organic matter will not only attract microorganisms,but part of it will be“transformed”,under a degradation process,to microorganisms.This mechanism introduces a supplementary linear coupling term in thefirst equation of this model.Many authors(see[11—13]and references therein)already considered reaction coupling terms,but under some restrictive conditions,which are not verified here.This feed-back in the chemotactic equation is not compatible with mass conservation of microorganisms,unlike in[12,15].Furthermore,neither the boundedness of the microorganisms total mass nor the positivity and the boundedness(existence of threshold)of the solution remain immediate,unlike in[11,13—14].

Our main concern here is to prove the existence of a unique solution to this minimal MOMOS model improved by adding chemotaxis effect.We consider two chemotactic functions h,a“classical” one,h(u)=u,and a second one which prevents overcrowding of microorganisms,h(u)=u(M−u)if 0≤u≤M and zero otherwise,proposed by Wrzosek[15].

This paper is organized as follows.Section 2 introduces some notations,results and tools used throughout the paper.Section 3 presents sufficient conditions to get global solutions,and to prove the existence of an exponential attractor,in the case where h(u)=u.Section 4 is concerned with the second chemotactic function,where the chemotactic term cancels when u achieves the threshold M,which helps to prove that any local solution is actually global.In Sections 3—4 the domain is two dimensional.In Section 5,still keeping the second form of h and for domains of dimension less than or equal to 3(the dimension 3 is particularly interesting in applications),we prove the existence of a unique solution,with less restrictions on the initial conditions and forcing term than in Section 4.In Appendix 1 we prove that chemotactic term in system(Ph)is mandatory to obtain Turing patterns and in Appendix 2 we give some numerical simulations.

2 Mathematical Preliminary and Notations

Unless it is explicitly indicated,Ω is a bounded region in R2of C3class,the constants a,β,q,d,k1and k2are nonnegative,and f is a nonegative function belonging to an admissible space to befixed later.In all that follows C denotes a positive constant which may vary from line to line.

We recall here some known results(see[16—17]and references therein)that will help afterwards.

Interpolation spaceFor 0≤ s0＜ s＜ s1＜ ∞,Hs(Ω)is the interpolation space[Hs0(Ω),Hs1(Ω)]θwith s=(1−θ)s0+θs1between Hs0(Ω)and Hs1(Ω).Furthermore,we have

Embedding theoremWhen 0＜ s＜ 1,Hs(Ω)⊂ Lp(Ω)forwith the estimate

When s=1,H1(Ω)⊂ Lq(Ω)for any 1≤ q＜ ∞ and

where 1≤p≤q＜∞.

When s＞ 1,Hs(Ω)⊂with continuous embedding.

Fractional power of the Laplace operator(see[12]and[17,Chapter 2.7])Let a0,a1＞0 be constants and L=−a1Δ+a0be the Laplace operator equipped with the Newman boundary conditions,with the domain D(L)=Thus L is a positive definite self-adjoint operator of L2(Ω).For θ＞ 0,the fractional power on L is defined and noted Lθand Lθis also a positive definite self-adjoint operator on L2(Ω).More

with the norm equivalence.

There are some useful inequalities.

Biler’s lemma(see[18])Let 0 ≤ u ∈ H1(Ω)and:= ‖(u+1)log(u+1)‖L1.For any η＞ 0,

where p(·)denotes here some increasing function.

Local existenceWe needfirst to prove the existence of local solution of(Ph).For this purpose,we use the result obtained by Yagi and based on the Galerkin method(see[12,17]).

Let V and H be seperable Hilbert spaces with dense and compact embedding V⊂H.Let V′be the dual space of V and identify H and H′to get

The duality product between V and V′is denoted by 〈·,·〉.It coincides with the scalar product on H denoted by(·,·).

Consider the following Cauchy problem of a semilinear abstract differential equation:

in the space V′.

Here,A is the positive definite self-adjoint operator of H defined by a symmetric sesquilinear

Assumptions on a(·,·)

with constants δ,M ＞ 0.The operator A is also bounded from V to V′.

Assumptions on G(·)G(.)is a continuous function from V to V′,which satisfy

(g.i)For each ζ＞ 0,there exists an increasing continuous function φζ:[0,∞) → [0,∞)such that

(g.ii)For each ζ＞ 0,there exists an increasing continuous function ψζ:[0,∞) → [0,∞)such that

Finally F(·) ∈ L2(0,T;V′)is a given function and U0∈ H is an initial value.Then,we have the following result(see[12]).

Theorem 2.1Under Assumptions(a.i),(a.ii),(g.i)and(g.ii)and for everyF(·) ∈ L2(0,T;V′)andU0∈H,there exists a unique local solution U of(2.7)such that

HereT(U0,F)is determined by the norm‖U0‖Hand‖F‖L2(0,T;V′).

3 First Case:h(u)=u

3.1 Local existence and positivity

Let ε0arbitrarily fixed,ε0∈ (0,1).

Theorem 3.1Letu0∈ L2(Ω),v0∈ H1+ε0(Ω)andf ∈ L2(0,T;Hε0(Ω))be nonnegative functions.Then(Ph)has a unique nonnegative local solution on an interval[0,T0]such that

whereT0depends only on‖f‖L2(0,T;Hε0(Ω)),‖u0‖L2(Ω)and‖v0‖H1+ε0(Ω).

Proof First StepConstruction of a unique local solution.

with U=(u,v)∈V.

Then(Ph)is the following semilinear differential equation:

In order to apply the existence result of Theorem 2.1 to problem(3.1),let us verify the assumptions on a(·,·)and G(·).

The assumptions on a(·,·)are classically satisfied(see for example[12]).

For the conditions on G we have that for an arbitrary U=(u,v)∈ V and δ＞ 0,

and

Finally it is clear that

All these inequalities show that the condition(g.i)is fullfiled.

From the embedding theorem,we have

and using the interpolation theorem and Young inequality we obtain

for an arbitrary ζ＞ 0.On the other hand we have that

All these inequalities permit to show that condition(g.ii)is fullfiled too.

Second StepPositivity of the solution.

Now let us take the following semilinear system:

where A,F and Y0are defined as previously,and the mapping→ V′is defined by

By Theorem 2.1.there exists a local solution U=(u,v)on[0,T0]×Ω with T0depending only on U0and F.Let us define u+=max(u,0)and u−=max(−u,0).

We multiply thefirst equation by−u−and we integrate in space.So

for 0＜t≤T0.

Using Young inequality we get

with ε＞ 0 small enough and Cε＞ 0.

Remark 3.1If initial conditions U0and data f are not positive,this theorem proves anyway the existence of a local solution.However,as this is an ecology model,only nonnegative solutions make sense.

With minor changes due to our different problem(Ph),we prove as in[13]the following theorems.

Theorem 3.2LetU0=(u0,v0)∈ H1(Ω)× H2N(Ω)andf ∈ L2(0,T;H1(Ω)).Then there exists a unique local solutionU=(u,v)of(3.1)on an interval[0,TU0,f]such that

whereTU0,fis determined by‖f‖L2(0,T;H1(Ω)),‖u0‖H1(Ω)and‖v0‖H2(Ω).

3.2 Global existence

This section is devoted to proving the following result.

Theorem 3.4Letε0∈ (0,1)and letu0∈ L2(Ω),v0∈ H1+ε0(Ω)andf ∈ L2(0,T,H1(Ω))∩L∞(0,T;L2(Ω))be nonnegative functions.Then there exists a unique global and nonnegative solution(u,v)for the system(Ph)withh(u)=usuch that

ProofWe proceed in two steps.

First StepWe show that ‖v‖H1(Ω)and

We consider the function log(u+1).Since∇log(u+1)=,it follows that log(u+1)∈L2(0,T0;H1(Ω)).Noting that

we obtain from thefirst equation of(Ph)multiplied by log(u+1)that

So using Stokes theorem we deduce

Since

if we denote Ψ(t)= ‖(u+1)log(u+1)− u‖L1(Ω),we get

with arbitary ε,η ＞ 0.

From the second equation of(Ph)multiplied respectively by v and Δv,we obtain that

with arbitrary A,B＞0 and

with arbitrary C,D＞0.

By addition of thefirst two equation of(Ph)it follows that

which implies that for all t∈[0,T0]we have the inequality

Second StepWe take t1∈ (0,T0)so that v(t1)∈(Ω)and u(t1)∈ H1(Ω)and we set u(t1)=u1and v(t1)=v1.From Theorem 3.1 we already know that such a time t1exists,arbitrary small.In this step t varies in[t1,T0].From thefirst equation of(Ph)we have

From Young inequality and interpolation inequality(2.1)we get

with η ＞ 0 arbitrary.

Therefore(3.7)together with this yields that

In addition

with χ＞0 arbitrary.

Using Biler’s lemma(2.3)we verify from(3.7)that

with p a positive increasing function,depending on ‖f‖L∞(0,T;L2(Ω)),‖u0‖L2(Ω)and ‖v0‖H1+ε0(Ω)as well as the constant C＞0.

Thus we deduce the following inequality:

with p a positive increasing function depending on ‖f‖L∞(0,T;L2(Ω)),‖u0‖L2(Ω)and ‖v0‖H1+ε0(Ω),ξ＞0 an arbitrary constant.

On the other hand,we consider v as a solution of the Cauchy problem

in the space H1(Ω).Since k1u+f ∈ L2(t1,T0;H1(Ω))and v1∈ D(A2)=(Ω)it follows

Therefore

As D(A32)⊂ H3(Ω),we obtain

with some δ＞ 0.Let a1=min(a,k1)＞ 0.We now sum up(3.9)and(3.8)multiplied bywhere C＞0 is the constant appearing in(3.9).Then it follows that

with some constant C1＞ 0 independent of T0.Choosing ξ small enough we conclude that

with some constant C2＞ 0 dependent on ‖f‖L∞(0,T;L2(Ω))and the initial condition U0through‖u0‖L2and ‖v0‖H1+ε0,but independent of T0.The norms ‖u‖L2(t1,T0;H1(Ω))and ‖v‖L2(t1,T0;H3(Ω))do not depend on T0and hence those of‖u‖C([t1,T0];L2(Ω))and ‖v‖C([t1,T0];H2(Ω))do not depend either.

In particular this shows that the solution(u,v)can be extended as a weak solution beyond the T0.

3.3 Exponential attractor

Suppose that f is a positive constant function.Then we have the following result.

Proposition 3.1Letu0∈(Ω)andv0∈(Ω)be nonnegative functions.Letu,vbe the global solution of(Ph).Then,with some continuous increasing functionp(·)the following estimate holds:

for0＜t＜∞.

ProofUsing(3.10)we deduce the existence of two constants σ＞0 and C＞0 such that

Multiplying the first equation of(Ph)by Δu and integrating over Ω gives

where ε,ε′and Cε′ are positive constants derived from Young inequality.Using technical inequalities proved in[13,Proposition 4.1],we obtain

Taking η= ε2,ε=leads to

Take the second equation of(Ph)operated by Δ,choose Δ2v as a test function and integrate the product in Ω.After some calculations as in[13],we have

We sum(3.14)multiplied by γ and(3.13).Thus we obtain

Then for γ and ε′small enough,there exists a positive constant σ′such that

So,we can find χ ＞ 0 such as(3.11)is valid when σ = χ and

We verify also that

Finally,taking thefirst equation of(Ph)operated by∇ and multiplied by∇Δu,as in[13],gives

that is,

The termsRΩ|∇(∇ ·(u∇v)|2dx andRΩ|u∇u|2dx of(3.18)can be estimated(see[13,Proof of Proposition 4.1,Step 6])by

with an arbitrary η ＞ 0.Thus we obtain

Hence we can find a constant χ ＞ 0 such that(3.16)is valid and

To prove the existence of an exponential attractor,we will use the following result.

Proposition 3.2Letu0∈ L2(Ω),v0∈ H1+ε0(Ω)be nonnegative functions.Then there exists a continuous increasing functionp(·)independent ofu0andv0,such that

ProofSince the proof follows exactly the same ideas and technical difficulties as in the proof of[13,Theorem 4.6],we skip it here.

We can now prove the existence of an exponential attractor:Let H=L2(Ω)× H1(Ω)and consider the initial value problem

in H,with A as in Subsection 3.1 and D(A)=

We have proved already the existence of a unique global solution U=(u,v)continuous with respect to the initial condition U0.We define then a continuous semigroup{S(t)t≥0}on K by S(t)U0=U(t).For afixed t＞ 0,S(t)maps K into K ∩D(A).

Denote Br:={(u,v)∈ K;‖u0‖L2+‖v0‖H1+ε0≤ r}a bounded ball of K with radius r ＞ 0.

Proposition 3.3There exists a universal constantC＞0such that the following statement holds:For eachr＞0there exists a timetr＞0such that

ProofFix 0＜r＜∞.By trand Crwe denote some time and positive constant which depend on r but are uniform in U0∈Br,respectively.By Proposition 3.2,there exist a time trand a constant Crsuch that for t≥tr,

The desired estimate will be established step by step.

Let us add thefirst equation of(Ph)and the second one multiplied by 2 and let us integrate in space the result.If Φ(t):= ‖u(t)‖L1+2‖v(t)‖L1,we obtain

Thus

and we deduce

with C,c＞0 universal constants and Cr＞0 a constant depending on r.This shows that there exists a time denoted by trsuch that for all t≥tr,

with C＞0 a universal constant.

From(3.6)and(3.22)it follows that

Then there exists another time trand another universal constant C＞0 such that

From(3.11)and(3.21)we deduce that

and that there exist another time trand another constant C＞0,such that

Finally using(3.16),(3.20)and repeat the argument wefinish the proof.

Theorem 3.5(see[20,Theorem 3.1])LetΓ(t,U0)=S(t)U0be a mapping from[0,T]×HintoH.If G satisfies

andΓis such that

for eachT＞0,then there is an exponential attractorMfor({S(t)},H).

Thus we obtain the following result.

Theorem 3.6There exists an exponential attractorMof the dynamical system({S(t)}t≥0,H)inH

ProofSince the forcing term f is constant and the reaction coupling of thefirst equation of(E)is linear in U:k2v,the proof is the same as provided in[13,Theorem 5.1].

4 Second Case:h(u)=u(M−u)

Let M be a positive constant and consider a continuous function~h of h such as

Then we have the following result.

Proposition 4.1Letε0＞0andfbe a nonnegative function inL2(0,T;H1(Ω))∩′L∞(ΩT).For each nonnegative initial condition(u0,v0)inL2(Ω)× H1+ε0(Ω)there exists a constantT0such that0＜T0≤Tand a unique nonnegative solution(u,v)of(P~h)such that

ProofThe proof is essentially the same as in Subsection 3.2.

Moreover we can prove the following result.

Lemma 4.1Suppose thatMcondition(u0,v0)satisfies almost everywhere inΩthe following inequalities:

then the solution(u,v)of(P~h)satisfies

almost everywhere inΩT.

and

which completes the proof.

Remark 4.1(i)If the hypothesis of Lemma 4.1 are fulfilled,thanks to this lemma,the solution obtained in Proposition 4.1 is global in ΩT,

(ii)By Proposition 4.1 and Lemma 4.1 it follows that(u,v)is also a solution of(Ph)with h(u)=u(M−u).

The uniqueness of the solution is obtained in the following.

Theorem 4.1Letf ∈ L∞(ΩT)∩L2(0,T;H1(Ω))be a nonnegative function.Leth(u)=u(M − u)and suppose thatM ≥Then there exists a unique global solutionfor(Ph)which is nonnegative and such that

and

ProofWe skip here the proof of uniqueness since there is rigorously the same as in Theorem 5.1.

5 A Three Dimensional Domain

In order to prove the global existence of a solution of system(Ph),we supposed in the previous sections that Ω was a two dimensional domain and the initial conditions(u0,v0)∈L∞(Ω)×H1+ε0(Ω)were nonnegative and verifying some regularity conditions.These conditions are quite restrictive for a model of soil organic carbon and three dimensional domains are obviously more relevant in applications than bidimensional ones.

In this section we prove that if Ω is of dimension less than or equal to 3,if h=(4.1)and if both initial conditions and forcing term are nonnegative and less regular that in the previous section:(u0,v0) ∈ (L2(Ω))2and f ∈ L2(0,T;L2(Ω)),then the system(Ph)has a global nonnegative solution.Furthermore,if(u0,v0)∈ (L∞(Ω))2and f ∈ L∞(ΩT),then the solution is unique.

Here we use the following setting:

where(u0,v0)∈ (L2(Ω))2,f ∈ L2(0,T;L2(Ω))and u is a function in X=L2(ΩT).

For the sake of simplicity we take dim(Ω)=3,since all results remain the same if dim(Ω)＜ 3.

We will apply the Schauderfixed point theorem but let usfirst gather some more information.

First StepInvariant ball.

For any function u∈X the existence of a unique local solution of(P-S)(uu,vu)follows by direct application of Theorem 2.1.Additionally we have the following result.

Proposition 5.1Let(u0,v0)∈ (L2(Ω))2andf ∈ L2(0,T;L2(Ω)).

ProofTo prove that(uu,vu)is global in time,we multiply thefirst equation by uuand the second by vuand use Young inequality to get

Since u0,v0∈ L2(Ω)we deduce that uu,vuare bounded in L∞(0,T;L2(Ω))∩L2(0,T,H1(Ω)),and this bound does not depend on u.Using interpolation technique we obtain that uuis bounded in L4(0,T;(L3(Ω))and consequently|uu|uuis bounded in L2(0,T;(L32(Ω)),independent of u.

Combining Hö lder inequality,the boundedness of uu,vuin L2(0,T,H1(Ω))and L4(0,T;(L3(Ω))and the continous injection of L2(0,T;H1(Ω))into L2(0,T;(L3(Ω)),we obtain that∂tuu, ∂tvuare bounded in L2(0,T;(H1(Ω))′),independent of u.So we finish the proof.

We can then define the mapping Π :X → X such that uu= Π(u)is the unique solution of(P-S).From(5.1)the ball BR⊂ X is invariant by Π.

Second StepCompactness of Π(BR).The second statement of the previous proposition implies that Π(BR)⊂ {u ∈ W, ‖u‖W≤ C}.But the embedding of W into L2(0,T,L2(Ω))is compact thanks to the Aubin-Lions lemma.

Third StepΠ is a continuous mapping.Let zn∈ BRsuch that zn→ z in L2(ΩT)strong and let un= Π(zn).Then Un=(un,vn)satisfies the system(P-S)n:

Since the sequence(un,vn)n≥1is bounded in W2and(L∞([0,T];L2(Ω)))2,there exists by the Aubin-Lions lemma a subsequence(not relabeled)such that

To prove that∇u= ξ,we take a test function ϕ ∈ (D(ΩT))3,so that

where the last assertion is a straightforward consequence of the upper bound of sequence|un|unin L2(0,T,L32(Ω))and the a.e.convergence of the sequence(un)n≥1in ΩT.We obtain also a similar convergence for vntowards v as in(5.3).

Finally,thanks to suitable choices of test functions it follows that the limit function v is a solution of the following problem:

and we get u=Π(z).

By the uniqueness of the solution(u,v)of(5.5),we deduce that all the sequence converges.We conclude that Π is a continuous mapping.

We can now apply the Schauderfixed point theorem to prove the existence statement of the following result.

Proposition 5.2Letfbe a nonnegative function inL2(0,T;L2(Ω)).For each couple of nonnegative functions(u0,v0)∈ (L2(Ω))2,there exists a nonnegative solution for the problem(Ph)withh=

To prove the positivity of the solution,we proceed as in Section 4:We multiply thefirst equation by−u−and the second by−v−,we integrate in space and we add the two equations.Thanks to the identityh(u)∇v∇u−=0,a straightforward calculation gives:

with C ＞ 0.Wefinish the proof by applying the Gronwall lemma.

For the uniqueness of solution of problem(Ph)we have the following result.

Theorem 5.1Letf ∈ L∞(ΩT)be a nonnegative function.Consider(u0,v0)∈ (L∞(Ω))2such that0≤ u0≤ Mand0≤ v0(x)≤ vMalmost everywhere inΩ,wherevMis a positive constant.Then there exists a constantα≥0such that

and the solution of problem(Ph)is unique,whenh=~h.

Multiplying thefirst equation by~u+and the second by~v+and then adding the two equations gives

We sum up(5.8)—(5.9)multiplied by a constant σ ＞ 0 small enough,and we use(5.10)and(5.11)with a wise choise of ε,ε′and σ such that ε+ε′≤a and σ≤d.Thereby we prove the existence of a constant C＞0 such that

The Gronwall lemma entail that=0 for every t ∈ [0,T],which completes the proof.

6 Appendices

6.1 Appendix 1

Non-emergence of spatial patterns in(Ph)model without chemotaxis term(β=0)

Firstly we considere the PDEs system(Ph)without chemotaxis term(β=0).As in Lotka-Volterra systems(see[20]),also known as the predator-prey equations,diffusion alone cannot disturb a constant equilibrium,and so spatial heterogeneity cannot emerge.Using the following notation:

we obtain the following non-dimensional equations(we revoke the notation):

with the same initial conditions and boundary conditions as(Ph)system.Without diffusion,the system(6.1)has a unique positive steady state:

To assess the steady state stability,the system is linearised around(u∗,v∗).Setting

where 0＜ε≪1,gives the following linear system:

with no-flux boundary conditions.

As in[20—21],we look for a solution of the form:

Let k=|k|be the Euclidean norm of the wave vector.We obtain the following eigenvalue problem:

where A is the two by two matrix

The eigenvalue ρ depends on k.

Turing instability occurs(which means that spatial patterns appear)when ρ(k2) ＞ 0,for a given value of k.But the matrix A has a strictly negative trace and a positive determinant,and so ρ(k2)＜ 0 for all values of k.Hence no patterns will emerge in this case.

Emergence of spatial patterns in(Ph)model with β＞0

Finally,for the model with both diffusion and chemotaxis,it can be proven that the equilibrium solutions of the equation system(Ph)can be rendered non-stable under certain conditions,and thus produce patterns and spatial heterogeneity.As in the previous section,the system(Ph)was linearised around the steady state(u∗,v∗).We obtain the following system:

where

Looking for solutions like in(6.4),the following eigenvalue problem must be solved:

where B is the two by two matrix

In this case,the trace of matrix B is strictly negative while its determinant can be strictly negative for some values of k.Thus,taking chemotaxis into account in the system may lead to the emergence of spatial patterns.

6.2 Appendix 2

Numerical simulations

A set of validated parameters derived from studies published[7]was used to run numerical simulations.The data used came from an Andean Pramo site near Gavidia,Venezuela.As pattern geometries depend on the shape of the spatial domain(see[20]),two different forms of spatial domain were tested.Figures below show the numerical simulations of the soil microbial biomass compartment for the nearly rectangular and circular domains,using either h(u)=h1(u)=u which does not prevent explicitly any overcrowding(Figures 1—2),or h(u)=h2(u)=u(M − u)which explicitly does prevent overcrowding(Figures 3—4).These figures show the spatial variability and patterns obtained for soil microbial biomass after 60 days and for the two spatial domain shapes.The soil microbial biomass pattern agrees with the distribution within the soil matrix of the microbial hot spots at micron scale.Numerical simulations were performed using COMSOL Multiphysics 5.0.

Figure 1 Spatial microbial biomass distribution when h=h1after 60 days.

Figure 3 Spatial microbial biomass distribution when h=h2after 60 days.

Figure 4 Spatial microbial biomass distribution when h=h2after 60 days.

AcknowledgementWe thank DR Martial Bernoux-UMR Eco&Sols(INRA,SupAgro,CIRAD,IRD)for his expertise and advice on soil carbon dynamics.

[1]Vogel,C.,Mueller,C.W.,Hschen,C.,et al.,Submicron structures provide preferential spots for carbon and nitrogen sequestration in soils,Nature Communications,5(2947),2014.

[2]Manzoni,S.and Porporato,A.,Soil carbon and nitrogen mineralization:Theory and models across scales,Soil Biology and Biochemistry,41,2009,1355–1379.

[3]Elzein,A.and Balesdent,J.,Mechanistic simulation of vertical distribution of carbon concentrations and residence times in soils,Soil Science Society of America Journal,59,1995,1328–1335.

[4]Goudjo,C.,Leye,B.and Sy,M.,Weak solution to a parabolic nonlinear system arising in biological dynamic in the soil,International Journal of Differential Equations,2011,Article ID:831436,24 pages.

[5]Deckmyn,G.,Campioli,M.,Muys,B.,et al.,Simulating C cycles in forest soils:Including the active role of micro-organisms in the ANAFORE forest model,Ecological Modelling,222,2011,1972–1985.

[6]Pansu,M.,Bottner,P.,Sarmiento,L.,et al.,Comparison offive soil organic matter decomposition models using data from a 14C and 15N labeling field experiment,Global Biogeochemical Cycles,18,2004,1–11.

[7]Pansu,M.,Sarmiento,L.,Rujano,M.A.,et al.,Modeling organic transformations by microorganisms of soils in six contrasting ecosystems:Validation of the MOMOS model,Global Biogeochemical Cycles,24,2010.

[8]Hammoudi,A.,Iosifescu,O.and Bernoux,M.,Mathematical analysis of a nolinear model of soil carbon dynamics,Differ.Equ.Dyn.Syst.,23(4),2015,453–466.

[9]Hammoudi,A.,Iosifescu,O.and Bernoux,M.,Mathematical analysis of a spatially distributed soil carbon dynamics model,Analysis and Applications,2016,DOI:10.1142/SO219530516500081.

[10]Keller,E.F.and Segel,L.A.,Initiation of slime mold aggregation viewed as an instability,J.Theor.Biol.,26,1970,399–415.

[11]Bendahmane,M.,Mathematical analysis of reaction-diffusion system modeling predator-prey with preytaxis,Networks and Heterogeneous Media,3(4),2008,863–879.

[12]Ryu,S.-U.and Yagi,A.,Optimal Control of KellerSegel Equations,Journal of Mathematical Analysis and Applications,256(1),2001,45–66.

[13]Osaki,K.,Tsujikawa,T.,Yagi,A.,et al.,Exponential Attractor for a Chemotaxis-growth System of Equations,Nonlinear Anal.,51(1),2002,119–144.

[14]Hillen,T.and Painter,K.,Global existence for a parabolic chemotaxis model with prevention of overcrowding,Advances in Applied Mathematics,26(4),2001,280–301.

[15]Wrzosek,D.,Volume filling effect in modelling chemotaxis,Mathematical Modelling of Natural Phenomena,2010,5(1),123–147.

[16]Lions J.-L.,Quelques Methodes de Resolution des Problemes aux Limites Non Lineaires,Dunod,Gauthier-Villars,Paris,1969.

[17]Yagi,A.,AbstractParabolic Evolution Equations and Their Applications,Springer-Verlag,Berlin,Heidelberg,2010.

[18]Biler,P.,Hebisch,W.and Nadzieja T.,The Debye system:Existence and large time behavior of solutions,Nonlinear Anal.,1994,238,1189–1209.

[19]Temam,R.,Navier-Stokes Equations.Theory and Numerical Analysis,North-Holland,Amsterdam,New York,Oxford,1977.

[20]Murray,J.D.,Mathematical Biology I:An Introduction,Springer-Verlag,New York,2002.

[21]Turing,A.M.,The Chemical basis of Morphogenesis,Philosophical Transactions of the Royal Society of London,1952,237(641),37–72.