# Analysis of a Singular Convection Diffusion System Arising in Turbulence Modelling

Authors: P. Dreyfuss
Journal: International Journal of Partial Differential Equations
 (2013)

---

## Abstract

We shall study some singular stationary convection diffusion system governing the steady state of a turbulence model closely related to thek-ε one. We shall establish existence, positivity, and regularity results in a very general framework.


---

## Body

## 1. Introduction

We shall first recall some basic ideas concerning the statistical turbulence modelling for fluids. The reader can consult [1, 2] for a more detailed introduction.Letu, p, ρ, and T be the velocity, pressure, density, and temperature of a Newtonien compressible fluid. Let also Ω~⊂ℝ3 be a domain which is assumed to be bounded. Then the motion of the flow in Ω~ at a time t∈ℝ+ can be described by the compressible Navier Stokes equations (see system (C) page 8 in [3]). It is well known that direct simulation based on such a model is harder or even impossible at high Reynolds numbers. The reason is that too many points of discretization are necessary, and so only very simple configurations can be handled.Thus, engineers and physicists have proposed new sets of equations to describe the average of a turbulent flow. The most famous one is thek-ε model, introduced by Kolmogorov [4]. We shall briefly present its basic principles in the following. Let v denote a generic physical quantity subject to turbulent (i.e., unpredictable at the macroscopic scale); we introduce its mean part (or its esperance) 〈v〉 by setting:
(1)〈v(x,t)〉=∫𝒫‍v(ω,x,t)d𝒫(ω),
where the integral is taken in a probablistic context which we shall not detail any more here. Note, however, that the operation 〈·〉 is more generally called a filter. The probablistic meaning is one but not the only possible filter (see, for instance, [1] chapter 3). We shall then consider the decomposition: v=〈v〉+v′, where v′ is referred to the noncomputable or the nonrelevant part and 〈v〉 is called the mean part (i.e., the macroscopic part).The principle of thek-ε model is to describe the mean flow in terms of the mean quantities 〈u〉, 〈p〉, 〈ρ〉, and 〈T〉 together with two scalar functions k and ε, which contains relevant information about the small scales processes (or the turbulent processes). The variable k (SI: [m2/s2]) is called the turbulent kinetic energy, and ε [m2/s2] is the rate of dissipation of the kinetic energy. They are defined by
(2)k=12|u′|2,ε=ν2〈|∇u′+(∇u′)T|2〉,
where ν is the molecular viscosity of the fluid. The model is then constructed by averaging (i.e., by appling the operator 〈·〉 on) the Navier-Stokes equations. Under appropriate assumptions (i.e., the Reynolds hypothesis in the incompressible case and the Favre average in the compressible case), we obtain a closed system of equations for the variables 〈u〉, 〈p〉, 〈ρ〉, 〈T〉, k, and ε (see [1] pages 61-62 for the incompressible case and pages 116-117 in the compressible situation).Here, we shall focus on the equations fork and ε, and we consider that the others quantities 〈u〉, 〈p〉, 〈ρ〉, and 〈T〉 are known. Moreover, in order to simplify the readability, we do not use the notation 〈·〉; that is, in the sequel we will write u instead of 〈u〉 and ρ instead 〈ρ〉 to represent the mean velocity and density of the fluid. The equations for k and ε are of convection-diffusion-reaction type:
(3)∂tk+u·∇k-cνρdiv((ν+ρk2ε)∇k)=cνk2εF-23kD-ε,∂tε+u·∇ε-cερdiv((ν+ρk2ε)∇ε)=c1kF-2c13cνεD-c2ε2k,
where D(x,t):=divu(x,t), F(x,t):=(1/2)|∇u+(∇u)T|2-(2/3)D(x,t)2≥0 (see Appendix B in the appendices), and cν, cε, c1, and c2 are generally taken as positive constants (see (A.3) in the appendices).Note that (3) is only valid sufficiently far from the walls. In fact, in the vicinity of the walls of the domain Ω~, there is a thin domain Σ called logarithmic layer in which the modulus of the velocity goes from 0 to 𝒪(1). In this layer we can use some wall law or a one equation model (see [1] chap. 1 and [5]) instead of (3). However, (3) can be considered even in the logarithmic layer if we allow the coefficients cν,cε,c1, and c2 to depend appropriately on some local Reynolds numbers (see [1, pages 59-60 and page 115]). In this last situation the system is called Low-Reynolds number k-ε model.In the following, we focus on the study in the domainΩ:=Ω~∖Σ, and we assume that its boundary ∂Ω is Lipschitz (see [6, page 127]). We denote by n(x) the outward normal defined for almost all points x∈∂Ω. The boundary conditions for k and ε on ∂Ω are then on the form
(4)k=k0,ε=ε0on∂Ω,
where k0 and ε0 are strictly positive functions which can be calculated by using a wall law (see [1, page 59]) or a one equation model (see [5]). In the following, we assume that k0 and ε0 are given. Moreover, we can assume (see again [1] page 59) that
(5)u·n=0on∂Ω.We shall concentrate in this paper on a modified system obtained after introducing the new variablesθ [s] and φ [m−2] given by
(6)θ=kε,φ=ε2k3.
These variables have a physical meaning (see [2]): θ represents a characteristic time of turbulence and L=φ-1/2 is a characteristic turbulent, length scale. By using this change of variable in (3) and after considering some modelisation arguments for the diffusion processes (see the appendices), we obtain
(P)∂tθ+u·∇θ-1ρdiv((ν+cθρθφ)∇θ)=-c3θ2F+c4θD+c5,∂tφ+u·∇φ-1ρdiv((ν+cφρθφ)∇φ)=-φ(c6θF-c7D+c8θ-1),
where the coefficients cθ, cφ, and ci are all positive.Problem(P) differs from the k-ε one only by the diffusive parts, and it is attractive by some stronger mathematical properties. Another model closely related to these systems, and having some popularity, is the k-ω one (with ω=θ-1; see for instance [2]).In the papers [7, 8], the authors have established the existence of a weak solution for (P) and a property of positiveness. This last feature takes the model useful in practice: it can be used directly or also as an intermediary stabilization procedure to the k-ε one (see [9]). Another important property attempted for a turbulence model is its capability to predict the possible steady states. In the previous works (except in [10]) only the evolutive version of (P) was studied under very restrictive assumptions. In [10], however, the stationary problem is studied, but it is simplified by considering a perturbation of the viscosities that artificially cancel the singularity of the system.Hence, in this paper we shall study the stationary version of(P) on a bounded domain Ω⊂ℝN, N=2 or 3, on which we impose the boundary conditions θ=a, φ=b on ∂Ω. Remark that, by using (4) together with (6), we obtain
(7)a=k0ε0,b=ε02k03.
Hence, we can assume that a and b are strictly positive given functions.We shall establish existence, positivity, and regularity results in a very general framework.


## 2. Main Results

### 2.1. Assumptions and Notations

Let(Q) denote the stationary system associated with (P). For simplicity we introduce the new parameters Cind:=ρcind where the subscript “ind” takes the integer values 3, 4, 5, 6, 7, and 8 or the letters θ and φ. Then, our main model (Q) has the following form:
(Q)ρu·∇θ-div((ν+Cθθφ)∇θ)=-C3Fθ2+C4θD+C5inΩ,ρu·∇φ-div((ν+Cφθφ)∇φ)=-φ(C6θF-C7D+C8θ-1)inΩ,θ=a,φ=bon∂Ω.
For physical reasons we are only interested in positive solutions (θ,φ) for (Q). Note, however, that even with this restriction, the problem (Q) may be singular (i.e., the viscosities ν+(Cθ/θφ) and ν+(Cφ/θφ) may be unbounded). Moreover, because we allow ν≡0, the equations may degenerate (i.e., the viscosities may vanish). Hence, without additional restriction there may be various nonequivalent notions of weak solution (see for instance [11]).In fact a good compromise between respect of the physics, simplification of the mathematical study, and obtention of significative results is to restrictθ and φ to be within the class 𝒮 defined by
(8)𝒮={f∈H1(Ω)∩L∞(Ω),f-1∈L∞(Ω)f:Ω⟶ℝ+suchthatf∈H1(Ω)∩L∞(Ω),f-1∈L∞(Ω)}.
In particular, if the parameters appearing in (Q) are sufficiently regular and if we restrict θ and φ to be within the class 𝒮, then the notion of a weak solution for (Q) is univocally defined: it is a distributional solution (θ,φ) that satisfies the boundary conditions in the sense of the trace.In this last situation we will tell that (θ,φ) is a weak solution of (Q) in the class 𝒮.In order to be able to consider such a weak solution for(Q), we shall precise in the following some sufficient conditions of regularity for the data.LetN=2 or 3 denote the dimension of the domain Ω, and let r be a fixed number such that
(9)r>N2.
We then have the following continuous injection (see Lemma 5):
(10)Lr(Ω)⊂W-1,β(Ω),whereβ=r*>N.Recall thatD=div(u). We will consider the following assumptions.(i)
Assumptions onΩ:
(11)Ω⊂ℝNisopen,bounded,andithasaLipschitzboundary∂Ω.(ii)
Assumptions on the flow data (whenN=2, one assumption in (14) can be relaxed: u∈(Lp(Ω))2 with p>2 (instead of p=3) is sufficient, but this would not improve any result significantly.) u, F, D, ρ, and ν are as follows:
(12)ν≥0,(13)F,ρ:Ω⟶ℝ+,ρ,ρ-1∈L∞(Ω),F∈Lr(Ω),(14)u∈(L3(Ω))N,D∈Lr(Ω),div(ρu)=0,(15)u·n=0on∂Ω.(iii)
Assumptions on the turbulent quantities on the boundary are as follows:(16)a,b∈H1/2(∂Ω)∩L∞(∂Ω),a(x),b(x)≥δ>0a.e.x∈∂Ω,
where δ>0 is a fixed number.(iv)
Assumptions on the model coefficients are as follows:(17)Ci:Ω⟶ℝ+,C5,C8∈Lr(Ω),fori≠5,8:Ci∈L∞(Ω),(18)Cind:Ω×(ℝ+)2⟶ℝ+areCaratheodory,(19)Cind(x,v,w)≥αind>0∀v,w∈ℝ+,fora.a.x∈Ω,(20)∀v,w∈ℝ+,x⟶Cind(x,v,w)∈L∞(Ω),
where in (18)–(20) Cind means Cθ or Cφ. The assumption (18) signifies that for all v,w∈ℝ+, x→Cind(x,u,v) is measurable, and for a.a. x∈Ω:(v,w)→Cind(x,v,w) is continuous. This ensures that Cind(x,θ,φ) is measurable when θ and φ are measurable. The condition (19) means that Cind is uniformly positive, whereas (20) tells that Cind(x,θ,φ) remains bounded if θ and φ are bounded.We will study problem(Q) under themain assumption:
(21)(A0):(9)–(20)aresatisfied.Note that in the main situation(A0) the assumption (12) made for ν allows the possibility that ν≡0. In other words the molecular viscosity ν can be neglected in the model. This is often chosen in practice because the eddy viscosities Cθ/θφ and Cφ/θφ are dominant in the physical situations (see [1, 5]). Remark also that the coefficients Ci are allowed to depend on x, and the viscosity parameters Cθ, Cφ may depend on x, θ, φ.For a given functionf:Ω→ℝ, we shall use the notations f+ and f- to represent the positive and negative parts of f, that is
(22)f=f++f-,f+(x)=max(f(x),0)≥0,f-(x)=min(f(x),0)≤0.
We will also consider some assumptions of low compressibility of the form
(23)∥D+∥Lr(Ω)≤τ,
for some τ>0 that will be precised.This last kind of condition seems to be necessary (see the Appendices) in order to obtain a weak solution for(Q) in the three-dimensional case, whereas when N=2, we shall use some particularities of the situation to obtain a weak solution without any assumption of low compressibility. Nevertheless in this case we will assume that in addition to (A0) the following condition is satisfied:
(24)(A1):C42C3∈Lr,c6=0.014,c7=0.104,c8=0.84.
In this last condition the values for c6, c7, and c8 are in fact their classical constant values (see (A.8) in the Appendices)In the sequel we denote by DATA some quantity depending only on the data under the assumption(A0), that is,
(25)DATA=Const(Ω,a,b,∥u∥(L3)N,αθ,αφ,(∥Cj∥Lr)j=5,8,∥F∥Lr,(∥Ci∥L∞)i=3,4,6,7(∥Cj∥Lr)j=5,8).
Note that DATA does not depend on D and ν.The exact form of the dependency (i.e., the function Const) is allowed to change from one part of the text to another.


### 2.2. Main Results

We shall establish two theorems of existence. The first one applies ifN=2 or 3, and the second one is limited to the case N=2.Theorem 1.
Assume that (A0) holds. Then, there exists τ>0 such that, if  ∥D+∥Lr(Ω)≤τ, then problem (Q) admits at least one weak solution (θ,φ) in the class 𝒮.Theorem 2.
Assume thatN=2 and (A0), (A1) hold. Then, problem (Q) admits at least one weak solution (θ,φ) in the class 𝒮.In all the situations, we have the following regularity result.Theorem 3.
Let(θ,φ) be a weak solution of (Q) in the class 𝒮, and assume that (A0) is satisfied. We have the following. (i)
Ifu∈(L∞(Ω))N, and a,b are Hölder continuous, then (θ,φ)∈(𝒞0,α(Ω¯))2, for some α>0(ii)
Assume that in addition, the following conditions are satisfied:(26)∂Ω,aandbareofclass𝒞2,α,ρu∈(𝒞1,α(Ω¯))N,F∈𝒞0,α(Ω¯),∀i=3,…,8:Ci∈𝒞0,α(Ω¯),Cind(x,v,w)∈𝒞1,α(Ω¯×(ℝ+)2).
Then (θ,φ)∈(𝒞2,α(Ω¯))2, and it is a classical solution of (Q).


### 2.3. Discussion on the Results

Compared to previous works (see [1, 7, 12, 13]) our basic assumption (A0) made in Theorem 1 is very general. In particular, we do not artificially cancel the singularity in the model, and we only assume weak regularity on the data. For instance, the basic assumption we made on the mean flow is u∈(L3(Ω))N and divu∈Lr(Ω) with some r>3/2, whereas in the previous works it was assumed that u∈(L∞(Ω))N and divu∈L∞(Ω). Our condition is more interesting from a practical point of view because it is satisfied when u is a weak solution of the Navier Stokes equations. Hence, our work could be used for a future analysis of the full coupled system Navier-Stokes plus (Q).From a mathematical point of view the problem we study is a nonlinear, degenerate, and singular elliptic system. Several complications arise for its analysis. In particular, the balance between the increase/decrease of the source terms (i.e., the functionsCi, D, and F appearing in the second members of (Q), and the possible explosion/vanishing of the viscosities is difficult to control. The strategy followed here is to first carefully study some elliptic scalar problem (possibly degenerate and singular) of the form
(S)ρu·∇ζ-div((ν+κζ)∇ζ)=g(x,ζ)inΩ,ζ=ζ0on∂Ω.
By using and developing some techniques due to Stampacchia, we are able to establish existence, positivity, and regularity results for problem (S). These results (see Proposition 9), which also have an independent mathematical interest, are the key ingredients for proving Theorem 1. Under the additional assumption u∈(L∞(Ω))N, we give a Hölder continuity result for k and θ. Moreover, we establish an existence result for a classical solution under some smoothness assumptions on the data.


### 2.4. Organization of the Paper

(i)
In Section3, we shall recall some results concerning the truncature at a fixed level and the Stampacchia’s estimates. This last technique takes an important role in our analysis; moreover, we shall need a precise control of the estimates. Hence, we shall present it with some details and developments.
(ii)
In Section4, we introduce a sequence (Qn) of problems which approximate (Q). For n fixed (Qn) is a PDE system of two scalar equations of the form (S): one equation for the unknown θn+1 and one for φn+1. The point is that the unknowns θn+1 and φn+1 are only weakly coupled. The coupling of the two equations is essentially realized through the quantities θn and φn calculated at the previous step. Hence, we shall firstly study carefully the problem (S). The major tool used here is the Stampacchia’s estimates. We next prove that (Qn) is well posed. Hence, we obtain an approximate sequence of solutions (θn,φn) for problem (Q). Moreover, we prove that θn and φn are uniformly bounded from above and below, which are the key estimates.
(iii)
In Section5, we use the uniform bounds established in Section 4, in order to extract a converging subsequence from (θn,φn). We then prove that the limit (θ,φ) is a weak solution of (Q) in the class 𝒮.Under the additional assumption u∈(L∞(Ω))3, we are able to use the De Giorgi-Nash theorem, and we obtain an Hölder continuity result for θ,φ. By assuming in addition some smoothness properties for the data, we can iterate the Schauder estimates and prove Theorem 3.
(iv)
In the appendices, we present the derivation of theθ-φ model from the k-ε one. Moreover, we justify that the choice φ=ε2/k3 is valid even in the compressible situation. The justification uses in particular a property of positivity of the function F. We also discuss briefly the necessity of the low compressibility assumption when N=3. Finally, we recall a generalized version of the chain rule for G(u) where G is a Lipchitz function and u a Sobolev one.


## 2.1. Assumptions and Notations

Let(Q) denote the stationary system associated with (P). For simplicity we introduce the new parameters Cind:=ρcind where the subscript “ind” takes the integer values 3, 4, 5, 6, 7, and 8 or the letters θ and φ. Then, our main model (Q) has the following form:
(Q)ρu·∇θ-div((ν+Cθθφ)∇θ)=-C3Fθ2+C4θD+C5inΩ,ρu·∇φ-div((ν+Cφθφ)∇φ)=-φ(C6θF-C7D+C8θ-1)inΩ,θ=a,φ=bon∂Ω.
For physical reasons we are only interested in positive solutions (θ,φ) for (Q). Note, however, that even with this restriction, the problem (Q) may be singular (i.e., the viscosities ν+(Cθ/θφ) and ν+(Cφ/θφ) may be unbounded). Moreover, because we allow ν≡0, the equations may degenerate (i.e., the viscosities may vanish). Hence, without additional restriction there may be various nonequivalent notions of weak solution (see for instance [11]).In fact a good compromise between respect of the physics, simplification of the mathematical study, and obtention of significative results is to restrictθ and φ to be within the class 𝒮 defined by
(8)𝒮={f∈H1(Ω)∩L∞(Ω),f-1∈L∞(Ω)f:Ω⟶ℝ+suchthatf∈H1(Ω)∩L∞(Ω),f-1∈L∞(Ω)}.
In particular, if the parameters appearing in (Q) are sufficiently regular and if we restrict θ and φ to be within the class 𝒮, then the notion of a weak solution for (Q) is univocally defined: it is a distributional solution (θ,φ) that satisfies the boundary conditions in the sense of the trace.In this last situation we will tell that (θ,φ) is a weak solution of (Q) in the class 𝒮.In order to be able to consider such a weak solution for(Q), we shall precise in the following some sufficient conditions of regularity for the data.LetN=2 or 3 denote the dimension of the domain Ω, and let r be a fixed number such that
(9)r>N2.
We then have the following continuous injection (see Lemma 5):
(10)Lr(Ω)⊂W-1,β(Ω),whereβ=r*>N.Recall thatD=div(u). We will consider the following assumptions.(i)
Assumptions onΩ:
(11)Ω⊂ℝNisopen,bounded,andithasaLipschitzboundary∂Ω.(ii)
Assumptions on the flow data (whenN=2, one assumption in (14) can be relaxed: u∈(Lp(Ω))2 with p>2 (instead of p=3) is sufficient, but this would not improve any result significantly.) u, F, D, ρ, and ν are as follows:
(12)ν≥0,(13)F,ρ:Ω⟶ℝ+,ρ,ρ-1∈L∞(Ω),F∈Lr(Ω),(14)u∈(L3(Ω))N,D∈Lr(Ω),div(ρu)=0,(15)u·n=0on∂Ω.(iii)
Assumptions on the turbulent quantities on the boundary are as follows:(16)a,b∈H1/2(∂Ω)∩L∞(∂Ω),a(x),b(x)≥δ>0a.e.x∈∂Ω,
where δ>0 is a fixed number.(iv)
Assumptions on the model coefficients are as follows:(17)Ci:Ω⟶ℝ+,C5,C8∈Lr(Ω),fori≠5,8:Ci∈L∞(Ω),(18)Cind:Ω×(ℝ+)2⟶ℝ+areCaratheodory,(19)Cind(x,v,w)≥αind>0∀v,w∈ℝ+,fora.a.x∈Ω,(20)∀v,w∈ℝ+,x⟶Cind(x,v,w)∈L∞(Ω),
where in (18)–(20) Cind means Cθ or Cφ. The assumption (18) signifies that for all v,w∈ℝ+, x→Cind(x,u,v) is measurable, and for a.a. x∈Ω:(v,w)→Cind(x,v,w) is continuous. This ensures that Cind(x,θ,φ) is measurable when θ and φ are measurable. The condition (19) means that Cind is uniformly positive, whereas (20) tells that Cind(x,θ,φ) remains bounded if θ and φ are bounded.We will study problem(Q) under themain assumption:
(21)(A0):(9)–(20)aresatisfied.Note that in the main situation(A0) the assumption (12) made for ν allows the possibility that ν≡0. In other words the molecular viscosity ν can be neglected in the model. This is often chosen in practice because the eddy viscosities Cθ/θφ and Cφ/θφ are dominant in the physical situations (see [1, 5]). Remark also that the coefficients Ci are allowed to depend on x, and the viscosity parameters Cθ, Cφ may depend on x, θ, φ.For a given functionf:Ω→ℝ, we shall use the notations f+ and f- to represent the positive and negative parts of f, that is
(22)f=f++f-,f+(x)=max(f(x),0)≥0,f-(x)=min(f(x),0)≤0.
We will also consider some assumptions of low compressibility of the form
(23)∥D+∥Lr(Ω)≤τ,
for some τ>0 that will be precised.This last kind of condition seems to be necessary (see the Appendices) in order to obtain a weak solution for(Q) in the three-dimensional case, whereas when N=2, we shall use some particularities of the situation to obtain a weak solution without any assumption of low compressibility. Nevertheless in this case we will assume that in addition to (A0) the following condition is satisfied:
(24)(A1):C42C3∈Lr,c6=0.014,c7=0.104,c8=0.84.
In this last condition the values for c6, c7, and c8 are in fact their classical constant values (see (A.8) in the Appendices)In the sequel we denote by DATA some quantity depending only on the data under the assumption(A0), that is,
(25)DATA=Const(Ω,a,b,∥u∥(L3)N,αθ,αφ,(∥Cj∥Lr)j=5,8,∥F∥Lr,(∥Ci∥L∞)i=3,4,6,7(∥Cj∥Lr)j=5,8).
Note that DATA does not depend on D and ν.The exact form of the dependency (i.e., the function Const) is allowed to change from one part of the text to another.


## 2.2. Main Results

We shall establish two theorems of existence. The first one applies ifN=2 or 3, and the second one is limited to the case N=2.Theorem 1.
Assume that (A0) holds. Then, there exists τ>0 such that, if  ∥D+∥Lr(Ω)≤τ, then problem (Q) admits at least one weak solution (θ,φ) in the class 𝒮.Theorem 2.
Assume thatN=2 and (A0), (A1) hold. Then, problem (Q) admits at least one weak solution (θ,φ) in the class 𝒮.In all the situations, we have the following regularity result.Theorem 3.
Let(θ,φ) be a weak solution of (Q) in the class 𝒮, and assume that (A0) is satisfied. We have the following. (i)
Ifu∈(L∞(Ω))N, and a,b are Hölder continuous, then (θ,φ)∈(𝒞0,α(Ω¯))2, for some α>0(ii)
Assume that in addition, the following conditions are satisfied:(26)∂Ω,aandbareofclass𝒞2,α,ρu∈(𝒞1,α(Ω¯))N,F∈𝒞0,α(Ω¯),∀i=3,…,8:Ci∈𝒞0,α(Ω¯),Cind(x,v,w)∈𝒞1,α(Ω¯×(ℝ+)2).
Then (θ,φ)∈(𝒞2,α(Ω¯))2, and it is a classical solution of (Q).


## 2.3. Discussion on the Results

Compared to previous works (see [1, 7, 12, 13]) our basic assumption (A0) made in Theorem 1 is very general. In particular, we do not artificially cancel the singularity in the model, and we only assume weak regularity on the data. For instance, the basic assumption we made on the mean flow is u∈(L3(Ω))N and divu∈Lr(Ω) with some r>3/2, whereas in the previous works it was assumed that u∈(L∞(Ω))N and divu∈L∞(Ω). Our condition is more interesting from a practical point of view because it is satisfied when u is a weak solution of the Navier Stokes equations. Hence, our work could be used for a future analysis of the full coupled system Navier-Stokes plus (Q).From a mathematical point of view the problem we study is a nonlinear, degenerate, and singular elliptic system. Several complications arise for its analysis. In particular, the balance between the increase/decrease of the source terms (i.e., the functionsCi, D, and F appearing in the second members of (Q), and the possible explosion/vanishing of the viscosities is difficult to control. The strategy followed here is to first carefully study some elliptic scalar problem (possibly degenerate and singular) of the form
(S)ρu·∇ζ-div((ν+κζ)∇ζ)=g(x,ζ)inΩ,ζ=ζ0on∂Ω.
By using and developing some techniques due to Stampacchia, we are able to establish existence, positivity, and regularity results for problem (S). These results (see Proposition 9), which also have an independent mathematical interest, are the key ingredients for proving Theorem 1. Under the additional assumption u∈(L∞(Ω))N, we give a Hölder continuity result for k and θ. Moreover, we establish an existence result for a classical solution under some smoothness assumptions on the data.


## 2.4. Organization of the Paper

(i)
In Section3, we shall recall some results concerning the truncature at a fixed level and the Stampacchia’s estimates. This last technique takes an important role in our analysis; moreover, we shall need a precise control of the estimates. Hence, we shall present it with some details and developments.
(ii)
In Section4, we introduce a sequence (Qn) of problems which approximate (Q). For n fixed (Qn) is a PDE system of two scalar equations of the form (S): one equation for the unknown θn+1 and one for φn+1. The point is that the unknowns θn+1 and φn+1 are only weakly coupled. The coupling of the two equations is essentially realized through the quantities θn and φn calculated at the previous step. Hence, we shall firstly study carefully the problem (S). The major tool used here is the Stampacchia’s estimates. We next prove that (Qn) is well posed. Hence, we obtain an approximate sequence of solutions (θn,φn) for problem (Q). Moreover, we prove that θn and φn are uniformly bounded from above and below, which are the key estimates.
(iii)
In Section5, we use the uniform bounds established in Section 4, in order to extract a converging subsequence from (θn,φn). We then prove that the limit (θ,φ) is a weak solution of (Q) in the class 𝒮.Under the additional assumption u∈(L∞(Ω))3, we are able to use the De Giorgi-Nash theorem, and we obtain an Hölder continuity result for θ,φ. By assuming in addition some smoothness properties for the data, we can iterate the Schauder estimates and prove Theorem 3.
(iv)
In the appendices, we present the derivation of theθ-φ model from the k-ε one. Moreover, we justify that the choice φ=ε2/k3 is valid even in the compressible situation. The justification uses in particular a property of positivity of the function F. We also discuss briefly the necessity of the low compressibility assumption when N=3. Finally, we recall a generalized version of the chain rule for G(u) where G is a Lipchitz function and u a Sobolev one.


## 3. Mathematical Background

In this section, we shall recall some results concerning the truncature at a fixed level and the Stampacchia’s estimates. This last technique takes an important role in our analysis; moreover, we shall need a precise control of the estimates. Hence, we shall present here the technique with some details and developments. As in the rest of the paper we denote byΩ⊂ℝN a bounded open Lipschitz domain. These properties for Ω are always implicity assumed if they are not precised.


### 3.1. Truncatures and Related Properties

The technique of Stampacchia is based on the use of special test functions which are constructed by using some truncatures. We shall recall some basic properties of the truncatures used in the paper. An important tool is the generalized chain rule (see TheoremD.1 in the appendices).Letl>0; we denote by Tl the truncature function Tl:ℝ→ℝ defined by
(27)Tl(s)={lifs>lsif-l≤s≤l-lifs<-l.
Let v∈H1(Ω). By applying Theorem D.1, we see that Tl(v)∈H1(Ω)∩L∞(Ω). Moreover, if we denote by Ωv,l the set {x∈Ωs.t.|v(x)|≤l}, then we have
(28)∇Tl(v)={∇vinΩv,l0elsewhere.
Note that Tl(·) truncates both the positive and the negative large values. In some cases, we need only to truncate the positive or the neghative side. For this reason, we introduce the semitruncatures Tl,+ and Tl,- defined by
(29)Tl,+(t)={lift>ltelsewhere,Tl,-(t)={-lift<-ltelsewhere.
We then have the decomposition: Tl=Tl,+∘Tl,-=Tl,-∘Tl,+.For givenvl∈H1(Ω) and s>0, we shall also consider
(30)ψs(v)=vl-Ts(vl)={vl-sinAs,l+vl+sinAs,l-0elsewhere,
where we have used the notations:
(31)As,l+:={vl≥s},As,l-:={vl≤-s},As,l=As,l+∪As,l-.
Let also ψs,± be the functions defined above (30) while replacing Ts by Ts,+ or by Ts,-. It is easy to verify that ψs,+ (resp., ψs,+) is in fact the positive (resp., the negative) part of ψs. In other words, we have
(32)ψs,+(vl)=ψs+(vl)=(vl-s)1As,l+≥0,ψs,-(vl)=ψs-(vl)=(vl+s)1As,l-≤0.ψs=ψs++ψs-.
The function ψs has the following properties:Lemma 4.
(i)ψs∈H1(Ω) and ∇ψs=1As,l∇vl,
(ii) ifγvl∈L∞(∂Ω), then for all s>∥γvl∥L∞(∂Ω), we have ψs∈H01(Ω).
Whereγ:H1(Ω)→H1/2(∂Ω) denotes the trace function.Proof.
Point (i) is a direct consequence of TheoremD.1. Property (ii) is proved in [14, Lemma  3.3 p. 53].


### 3.2. The Stampacchia Estimates

The Stampacchia estimates is a general method which allows one to obtain anL∞-estimate for the weak solution of a large class of elliptic PDEs of the second order. The L∞-estimate presented in the original work [15] depends on various quantities related to the PDE problem studied, but the exact dependency is not established. In our analysis we need a precise control of the L∞-estimates with respect to some quantities (in particular with respect to the diffusion coefficient of the PDEs). Hence, in the following, we take over and detail the technique in order to obtain a more precise L∞-estimates.The Stampacchia estimates are established by using the test functionsψs (or ψs±) defined previously, where in this case v (resp., vl) is a weak solution of the problem (resp., the sequence of problems) considered.For technical reasons we need a classical result concerning some relationships betweenLr functions and linear form on Sobolev spaces.Lemma 5.
Let1≤r<∞ and E∈Lr(Ω). Then
(33)E∈W-1,β(Ω),withβ=r*=NrN-r,
and there exists E~∈(Lβ(Ω))3 such that
(34)∫Ω‍Eφ=∫Ω‍E~·∇φ∀φ∈𝒟(Ω),∥E~∥(Lβ(Ω))3≤∥E∥Lr(Ω).
Moreover, we have
(35)r>N2⟹β>N.Proof.
Property (33) is easy to prove by using the Sobolev injection Theorem; together with the Hölder inequality: φ→∫Ω‍Eφ is a linear form on W01,p if p*:=Np/(N-p)≥r/(r-1). This last condition holds for p=P:=Nr/(r(N+1)-N). Hence, P′=p/(p-1)=Nr/(N-r)=r*.
We next obtain (34) by using a classical result (see [16] Proposition IX.20).
Finally, if we assume thatr>N/2 then β=r*>(N2/2)/(N-N/2)=N.The Stampacchia technique works in two steps: the first one is dependent on the problem (or the sequence of problems) studied, and the second one is independent of it. Here, the purpose is to present the key ingredients of these two steps. Because the first one is dependent on the problem studied we cannot present it here in its entirety, but we will consider a simple problem which contains the main technical points (in fact this introductive presentation will be useful to treat a more complicated class of problems in Section4). Let (vl)⊂H01(Ω) be a sequence of functions satisfying
(36)∫Ω‍νl∇vl∇φ=∫Ω‍glφ,∀φ∈H01,
where (νl) is a given sequence of strictly positive bounded functions and (gl)⊂Lr(Ω), with r>N/2. Let also ml, Ml denote the bounds from above and below for νl, that is,
(37)0<ml≤νl≤Ml<∞,a.e.inΩ.Step 1.
By testing (36) with ψs=ψs(vl), we obtain
(38)ml∫Ω‍|∇ψs|2=∫Ω‍glψs=byLemma5∫Ω‍El∇ψs,
with El∈Lβ, ∥El∥Lβ≤∥gl∥Lr for some β=β(r)>N.
Recall also thatSuppψs(vl)⊂Al,s. Hence, by using the Hölder inequality we obtain
(39)∫Ω‍El∇ψs≤∥∇ψs∥L2(∫Al,s‍|El|2)1/2≤∥∇ψs∥L2∥El∥Lβ(∫Al,s‍1)1/2(β/2)′≤∥∇ψs∥L2∥El∥Lβ|Al,s|(β-2)/2β.
Consequently (38) leads to
(40)∥ψs∥H01(Ω)22≤C|Al,s|Φ,withC>0,Φ>22*=N-2N.
This is the key estimate needed to pass at the second step which is independent of the problem studied.
Note that with the particular sequence of problems (36) chosen here, the constants C and Φ are
(41)C=∥gl∥Lr2ml2,Φ=β-2β.
Hence, Φ does not depend on l, s. Moreover, if we assume that (gl) is uniformly bounded in the Lr-norm and that (ml) is uniformly bounded from above by a strictly positive constant, then C is also independent of l,s.
This is an important point because we will see it hereafter; an estimate (40) with C and Φ independent of l, s leads to a uniform L∞ bound for (vl).Step 2.
Assume that we have obtained (40). We can obtain an L∞-estimate for vl as follows.
Let2*=2N/(N-2) be the Sobolev exponent associated to 2 in dimension N. By using the Poincaré-Sobolev inequality we have
(42)(∫Al,s‍|ψs|2*)2/2*≤C1∥ψs∥H01(Ω)2,C1=C1(|Ω|).
Let now t>s. It is clear that Al,t⊂Al,s, and consequently
(43)(∫Al,s‍|ψs|2*)2/2*≥(∫Al,t|ψs|2*‍)2/2*≥(∫Al,t|t-s|2*‍)2/2*≥|t-s|2|Al,t|2/2*.
We set
(44)χ(s):=|Al,s|,∀s≥0.
For fixed l, χ is a decreasing function, and from the estimates (40)–(43), we obtain
(45)χ(t)≤C2|χ(s)|2*Φ/2(t-s)-2*∀t>s≥0,C2=(C1C)2*/2.
Recall that we have assumed in (40) that Φ>2/2*. Hence, 2*Φ/2>1, and by using Lemma 4.1 in [15] we obtain
(46)χ(d)=0,d=22*Φ/(2*Φ-2)C21/2*|Ω|(Φ/2)-(1/2*)<∞.
This property tells exactly that
(47)∥vl∥L∞(Ω)≤d,d=C3(Φ,N,|Ω|)C.
In particular, d does not depend on l if the constants C and Φ appearing in (40) are independent of l. For instance with the particular sequence of problems (36) the constants C and Φ are given by (41), and if we assume that ∥gl∥Lr≤K, ml≥m>0, we obtain
(48)d=C~Km,C~=C~(Ω,N,r).Remark 6.
(i) If you are only interested in obtaining a uniform majoration or minoration forvl then instead of (40) it suffices to have
(49)∥ψs±∥H01(Ω)22≤C±|Al,s±|Φ±,withC±>0,Φ±>N-2N.
In fact in this case we consider χ±(s):=|Al,s±| instead of (44). This function is decreasing, and we obtain again (46). But now this property tells that
(50)±vl(x)≤d±a.e.inΩ,d±=C~±(Φ±,N,|Ω|)C±.
(ii) Let again (vl) be a sequence of functions satisfying (36), and assume that ml≥m>0. Then, we have
(51)gl≤hl,H:=supl∥hl∥Lr<∞⟹vl≤d1a.e.inΩ,(52)gl≥kl,K:=supl∥kl∥Lr<∞⟹vl≥-d2a.e.inΩ,withd1,d2>0,d1=C~Hm,d2=C~Km,C~=C~(Φ,N,|Ω|).
The proof of (51) and (52) is obtained by taking over the first step of the technique of Stampacchia: we use the test function ψs± instead of ψs.
In fact, the functionψs+ is positive. Hence, instead of (38), we have the following. 									(53)m∫Ω‍|∇ψs+|2=∫Ω‍glψs+≤∫Ω‍hlψs+=∫Ω‍El,+∇ψs+,
with El,+∈Lβ, ∥El,+∥Lβ≤∥hl∥Lr≤H. Consequently, in this case, we obtain (49) for the function ψs+, with C+=H2/m2, and the uniform majoration (51) follows.
The relation (52) can be proven by using ψs- as test function in (36). In fact, we remark now that ψs- is negative. Hence, instead of (38) we obtain
(54)m∫Ω‍|∇ψs-|2=∫Ω‍glψs-≤∫Ω‍klψs-=∫Ω‍El,-∇ψs-,
with El,-∈Lβ, ∥El,-∥Lβ≤∥kl∥Lr≤K. Consequently, we now obtain (49) for the function ψs-, with C-=K2/m2. This implies -vl≤d2, and consequently vl≥-d2.


## 3.1. Truncatures and Related Properties

The technique of Stampacchia is based on the use of special test functions which are constructed by using some truncatures. We shall recall some basic properties of the truncatures used in the paper. An important tool is the generalized chain rule (see TheoremD.1 in the appendices).Letl>0; we denote by Tl the truncature function Tl:ℝ→ℝ defined by
(27)Tl(s)={lifs>lsif-l≤s≤l-lifs<-l.
Let v∈H1(Ω). By applying Theorem D.1, we see that Tl(v)∈H1(Ω)∩L∞(Ω). Moreover, if we denote by Ωv,l the set {x∈Ωs.t.|v(x)|≤l}, then we have
(28)∇Tl(v)={∇vinΩv,l0elsewhere.
Note that Tl(·) truncates both the positive and the negative large values. In some cases, we need only to truncate the positive or the neghative side. For this reason, we introduce the semitruncatures Tl,+ and Tl,- defined by
(29)Tl,+(t)={lift>ltelsewhere,Tl,-(t)={-lift<-ltelsewhere.
We then have the decomposition: Tl=Tl,+∘Tl,-=Tl,-∘Tl,+.For givenvl∈H1(Ω) and s>0, we shall also consider
(30)ψs(v)=vl-Ts(vl)={vl-sinAs,l+vl+sinAs,l-0elsewhere,
where we have used the notations:
(31)As,l+:={vl≥s},As,l-:={vl≤-s},As,l=As,l+∪As,l-.
Let also ψs,± be the functions defined above (30) while replacing Ts by Ts,+ or by Ts,-. It is easy to verify that ψs,+ (resp., ψs,+) is in fact the positive (resp., the negative) part of ψs. In other words, we have
(32)ψs,+(vl)=ψs+(vl)=(vl-s)1As,l+≥0,ψs,-(vl)=ψs-(vl)=(vl+s)1As,l-≤0.ψs=ψs++ψs-.
The function ψs has the following properties:Lemma 4.
(i)ψs∈H1(Ω) and ∇ψs=1As,l∇vl,
(ii) ifγvl∈L∞(∂Ω), then for all s>∥γvl∥L∞(∂Ω), we have ψs∈H01(Ω).
Whereγ:H1(Ω)→H1/2(∂Ω) denotes the trace function.Proof.
Point (i) is a direct consequence of TheoremD.1. Property (ii) is proved in [14, Lemma  3.3 p. 53].


## 3.2. The Stampacchia Estimates

The Stampacchia estimates is a general method which allows one to obtain anL∞-estimate for the weak solution of a large class of elliptic PDEs of the second order. The L∞-estimate presented in the original work [15] depends on various quantities related to the PDE problem studied, but the exact dependency is not established. In our analysis we need a precise control of the L∞-estimates with respect to some quantities (in particular with respect to the diffusion coefficient of the PDEs). Hence, in the following, we take over and detail the technique in order to obtain a more precise L∞-estimates.The Stampacchia estimates are established by using the test functionsψs (or ψs±) defined previously, where in this case v (resp., vl) is a weak solution of the problem (resp., the sequence of problems) considered.For technical reasons we need a classical result concerning some relationships betweenLr functions and linear form on Sobolev spaces.Lemma 5.
Let1≤r<∞ and E∈Lr(Ω). Then
(33)E∈W-1,β(Ω),withβ=r*=NrN-r,
and there exists E~∈(Lβ(Ω))3 such that
(34)∫Ω‍Eφ=∫Ω‍E~·∇φ∀φ∈𝒟(Ω),∥E~∥(Lβ(Ω))3≤∥E∥Lr(Ω).
Moreover, we have
(35)r>N2⟹β>N.Proof.
Property (33) is easy to prove by using the Sobolev injection Theorem; together with the Hölder inequality: φ→∫Ω‍Eφ is a linear form on W01,p if p*:=Np/(N-p)≥r/(r-1). This last condition holds for p=P:=Nr/(r(N+1)-N). Hence, P′=p/(p-1)=Nr/(N-r)=r*.
We next obtain (34) by using a classical result (see [16] Proposition IX.20).
Finally, if we assume thatr>N/2 then β=r*>(N2/2)/(N-N/2)=N.The Stampacchia technique works in two steps: the first one is dependent on the problem (or the sequence of problems) studied, and the second one is independent of it. Here, the purpose is to present the key ingredients of these two steps. Because the first one is dependent on the problem studied we cannot present it here in its entirety, but we will consider a simple problem which contains the main technical points (in fact this introductive presentation will be useful to treat a more complicated class of problems in Section4). Let (vl)⊂H01(Ω) be a sequence of functions satisfying
(36)∫Ω‍νl∇vl∇φ=∫Ω‍glφ,∀φ∈H01,
where (νl) is a given sequence of strictly positive bounded functions and (gl)⊂Lr(Ω), with r>N/2. Let also ml, Ml denote the bounds from above and below for νl, that is,
(37)0<ml≤νl≤Ml<∞,a.e.inΩ.Step 1.
By testing (36) with ψs=ψs(vl), we obtain
(38)ml∫Ω‍|∇ψs|2=∫Ω‍glψs=byLemma5∫Ω‍El∇ψs,
with El∈Lβ, ∥El∥Lβ≤∥gl∥Lr for some β=β(r)>N.
Recall also thatSuppψs(vl)⊂Al,s. Hence, by using the Hölder inequality we obtain
(39)∫Ω‍El∇ψs≤∥∇ψs∥L2(∫Al,s‍|El|2)1/2≤∥∇ψs∥L2∥El∥Lβ(∫Al,s‍1)1/2(β/2)′≤∥∇ψs∥L2∥El∥Lβ|Al,s|(β-2)/2β.
Consequently (38) leads to
(40)∥ψs∥H01(Ω)22≤C|Al,s|Φ,withC>0,Φ>22*=N-2N.
This is the key estimate needed to pass at the second step which is independent of the problem studied.
Note that with the particular sequence of problems (36) chosen here, the constants C and Φ are
(41)C=∥gl∥Lr2ml2,Φ=β-2β.
Hence, Φ does not depend on l, s. Moreover, if we assume that (gl) is uniformly bounded in the Lr-norm and that (ml) is uniformly bounded from above by a strictly positive constant, then C is also independent of l,s.
This is an important point because we will see it hereafter; an estimate (40) with C and Φ independent of l, s leads to a uniform L∞ bound for (vl).Step 2.
Assume that we have obtained (40). We can obtain an L∞-estimate for vl as follows.
Let2*=2N/(N-2) be the Sobolev exponent associated to 2 in dimension N. By using the Poincaré-Sobolev inequality we have
(42)(∫Al,s‍|ψs|2*)2/2*≤C1∥ψs∥H01(Ω)2,C1=C1(|Ω|).
Let now t>s. It is clear that Al,t⊂Al,s, and consequently
(43)(∫Al,s‍|ψs|2*)2/2*≥(∫Al,t|ψs|2*‍)2/2*≥(∫Al,t|t-s|2*‍)2/2*≥|t-s|2|Al,t|2/2*.
We set
(44)χ(s):=|Al,s|,∀s≥0.
For fixed l, χ is a decreasing function, and from the estimates (40)–(43), we obtain
(45)χ(t)≤C2|χ(s)|2*Φ/2(t-s)-2*∀t>s≥0,C2=(C1C)2*/2.
Recall that we have assumed in (40) that Φ>2/2*. Hence, 2*Φ/2>1, and by using Lemma 4.1 in [15] we obtain
(46)χ(d)=0,d=22*Φ/(2*Φ-2)C21/2*|Ω|(Φ/2)-(1/2*)<∞.
This property tells exactly that
(47)∥vl∥L∞(Ω)≤d,d=C3(Φ,N,|Ω|)C.
In particular, d does not depend on l if the constants C and Φ appearing in (40) are independent of l. For instance with the particular sequence of problems (36) the constants C and Φ are given by (41), and if we assume that ∥gl∥Lr≤K, ml≥m>0, we obtain
(48)d=C~Km,C~=C~(Ω,N,r).Remark 6.
(i) If you are only interested in obtaining a uniform majoration or minoration forvl then instead of (40) it suffices to have
(49)∥ψs±∥H01(Ω)22≤C±|Al,s±|Φ±,withC±>0,Φ±>N-2N.
In fact in this case we consider χ±(s):=|Al,s±| instead of (44). This function is decreasing, and we obtain again (46). But now this property tells that
(50)±vl(x)≤d±a.e.inΩ,d±=C~±(Φ±,N,|Ω|)C±.
(ii) Let again (vl) be a sequence of functions satisfying (36), and assume that ml≥m>0. Then, we have
(51)gl≤hl,H:=supl∥hl∥Lr<∞⟹vl≤d1a.e.inΩ,(52)gl≥kl,K:=supl∥kl∥Lr<∞⟹vl≥-d2a.e.inΩ,withd1,d2>0,d1=C~Hm,d2=C~Km,C~=C~(Φ,N,|Ω|).
The proof of (51) and (52) is obtained by taking over the first step of the technique of Stampacchia: we use the test function ψs± instead of ψs.
In fact, the functionψs+ is positive. Hence, instead of (38), we have the following. 									(53)m∫Ω‍|∇ψs+|2=∫Ω‍glψs+≤∫Ω‍hlψs+=∫Ω‍El,+∇ψs+,
with El,+∈Lβ, ∥El,+∥Lβ≤∥hl∥Lr≤H. Consequently, in this case, we obtain (49) for the function ψs+, with C+=H2/m2, and the uniform majoration (51) follows.
The relation (52) can be proven by using ψs- as test function in (36). In fact, we remark now that ψs- is negative. Hence, instead of (38) we obtain
(54)m∫Ω‍|∇ψs-|2=∫Ω‍glψs-≤∫Ω‍klψs-=∫Ω‍El,-∇ψs-,
with El,-∈Lβ, ∥El,-∥Lβ≤∥kl∥Lr≤K. Consequently, we now obtain (49) for the function ψs-, with C-=K2/m2. This implies -vl≤d2, and consequently vl≥-d2.


## 4. Approximate Sequence and Estimates

Letg∈H1/2(∂Ω); we denote by ℛg its harmonic lifting that is:
(55)ℛg∈H1(Ω),ℛg=gon∂Ω,Δℛg=0onΩ.
We define the functions θ0 and φ0 by the formula
(56)θ0:=ℛaφ0:=ℛb.
Hence, by the using the maximum principle (see [16, p. 189] and [17]) together with the condition (16), we obtain
(57)0<δ≤θ0≤∥a∥L∞(∂Ω),δ≤φ0≤∥b∥L∞(∂Ω).Let nown≥0, (θn,φn) be given, and
(58)Cθ(n)(·):=Cθ(·,θn(·),φn(·)),Cφ(n)(·):=Cφ(·,θn(·),φn(·)).
In order to construct an approximate solution (θn+1,φn+1) for problem (Q), we introduce the following system:
(Qn)ρu·∇θn+1-div((ν+Cθ(n)θn+1φn)∇θn+1)=gθ(n)inΩ,ρu·∇φn+1-div((ν+Cφ(n)θn+1φn+1)∇φn+1)=gφ(n)inΩ,θn+1=a,φn+1=bon∂Ω,
where we used the notations
(59)gθ(n):=C5-C3Fθn+12+C4θn+1D,gφ(n):=-φn(C8θn+1-1+C6θn+1F-C7D).
For n∈ℕ we denote by (Hn) the following condition:
(Hn)θn,φn,θn-1,φn-1∈L∞(Ω),θn,φn≥0.
Let φmax be a fixed real number such that
(60)φmax>∥b∥L∞(∂Ω).
We shall also consider the condition: (Kn):=(Hn)+(φn≤φmax).Note that (57) shows that the condition (Kn) is satisfied for n=0. We will prove in the sequel that, under condition (Kn), we can obtain a weak solution (θn+1,φn+1) for problem (Qn), with moreover (θn+1,φn+1) satisfying the condition (Kn+1). This last property ensures the right definition of an approximate sequence. More precisely, we haveProposition 7.
Letn∈ℕ be given, and assume that (A0) is satisfied. Let also (θn,φn) be given and satisfy condition (Kn). There exists τ>0 depending only on DATA such that if ∥D+∥Lr(Ω)≤τ, then problem (Qn) admits at least one weak solution (θn+1,φn+1)⊂(H1(Ω)×L∞(Ω))2.
Moreover,(θn+1,φn+1) satisfies condition (Kn+1) and the estimates
(61)0<θmin≤θn+1≤θmax,0<φmin≤φn+1≤φmax,
where φmax was fixed in (60) and φmin, θmin, and θmax are positive numbers depending on DATA, but not on n.Remark 8.
Proposition7 is the key result that will be used later on to prove Theorem 1, whereas for Theorem 2 we shall establish and use a more simple version of this proposition (see Section 5.2).In order to prove the proposition we establish intermediate results.


### 4.1. Auxiliary Results

Letn∈ℕ, and let (θn,φn) be given and satisfy (Kn). We want to obtain (θn+1,φn+1) by solving (Qn), and in order to iterate the algorithm we also want that (θn+1,φn+1) satisfies (Kn+1).Remark that the system(Qn) is composed of two coupled scalar elliptic equations in divergence form, with a possible singular and degenerate structure. Hence, the goal of this subsection is to study this last kind of scalar problem.In order to do this, we first introduce a weightκ:ℝ+→ℝ+ which is assumed to be measurable and satisfying
(62)0<κ0≤κ≤κ1a.e.inΩ,
where κ0 and κ1 are two given reals.Let alsog:Ω×ℝ+→ℝ be a Caratheodory function (i.e., for all u∈ℝ+:x→g(x,u) is measurable, and for a.a. x∈Ω:u→g(x,u) is continuous).Let us consider the scalar problem(S) defined in Section 2 after Theorem 1, where ζ0∈H1/2(∂Ω)∩L∞(∂Ω), ζ0≥δ>0 a.e. in δΩ, is given. We always assume that ρ, u, ν, Ω which appear in (S) satisfy their corresponding conditions in (A0).Recall that we allowν≡0 in (A0). Hence, problem (S) may degenerate (i.e., the viscosity ν+(κ/ζ) may vanish) when ζ→∞. Moreover, (S) is singular (i.e., the viscosity tends to infinity) when ζ→0.We want now to find sufficient additional conditions forg that guarantee the existence of a bounded positive weak solution for problem (S). Hence, we shall consider
(63)γ(x):=supu∈[0,1]|g-(x,u)|∈Lr(Ω),(64)g+(x,u)≤γ1(x)+γ2(x)h(u),
where γ1,γ2∈Lr(Ω) and h:ℝ+→ℝ+ is continuous. In fact, more than establishing only the existence of a bounded positive solution for (S), we are interested in obtaining some uniform (with respect to κ1) bounds from above and below and some regularity results. We have the following.Proposition 9.
(i) Letκ satisfy (62), and let g:Ω×ℝ+→ℝ be a Caratheodory function satisfying (63), (64). There exists a real τ>0 depending on κ0, h such that if ∥γ2∥Lr(Ω)≤τ, then there exists a weak solution ζ∈H1(Ω)∩L∞(Ω) for problem (S). Moreover, we have
(65)0<ζmin≤ζ≤ζmax,ζmin=e-C/κ0,ζmax=eC/κ0,
where C depends only on γ, γ1, Ω, r, N, ζ0. In particular ζmin and ζmax are independent of κ1.
In addition, the following extended (wheng+≡0, it is a maximum principle) maximum principle holds:
(66)∥ζ∥L∞(Ω)≤∥ζ0∥L∞(∂Ω)+C~∥g+∥Lrκ0∥ζ∥L∞(Ω),C~=C~(Ω,N,r,γ1).
(ii) Assume that in additionu∈(L∞(Ω))N and ζ0 is Hölder continuous. Then, ζ∈𝒞0,α(Ω¯) for some 0<α<1. Moreover, if ∂Ω is of class 𝒞2,α, g∈𝒞0,α(Ω¯×ℝ+), ρu∈(𝒞1,α(Ω¯))N, κ∈𝒞1,α(Ω¯), and ζ0∈𝒞2,α(∂Ω), then ζ∈𝒞2,α(Ω¯), and it is a classical solution of (S).Before proving Proposition9, we establish an intermediate result. In a first step we consider the change of variable v=lnζ in (S), and for l∈ℕ we introduce a truncated version (Sl) of the system obtained:
(Sl)ρu∇eTl(v)-div((νeTl(v)+κ)∇v)=g(x,eTl(v))inΩ,v=lnζ0on∂Ω.
We then establish the following.Lemma 10.
(i) Letκ satisfy (62), and let g:Ω×ℝ+→ℝ be a Caratheodory function satisfying (63), (64). Then, for any l∈ℕ, there exists a weak solution v=vl∈H1(Ω)∩L∞ for the problem (Sl).
(ii) Let(vl) be the sequence given in (i), and let l≥1 be a fixed integer. Then, there exists τl>0 such that if the function γ2 in (63) satisfies ∥γ2∥Lr≤τl then we have
(67)∥vl∥L∞(Ω)≤C∥γ+γ1∥Lrκ0,C=C(DATA,ζ0).
In particular,C is independent of κ, ν, and l.Proof .
(i) By using the divergence formula, we obtain, for allw∈H1(Ω),
(68)∫Ω‍ρu∇eTl(v)·w=-∫Ω‍eTl(v)div(ρuw)+∫∂Ω‍eTl(v)ρwu·ndσ︸=0=-∫Ω‍eTl(v)ρu·∇w-∫Ω‍eTl(v)div(ρu)︸=0w.
Letv0:=ln(ℛζ0), and consider the change of variable v~:=v-v0. Then, problem (Sl) is equivalent to find v~∈H01(Ω) such that
(69)-divσ(x,v~,∇v~)=f(x,v~)in𝒟′(Ω),
where σ:Ω×ℝ×ℝN→ℝN and f:Ω×ℝ→ℝ are defined by
(70)σ(x,w,G)=(νeTl(w+v0)+κ(x))G-ρ(x)u(x)eTl(w+v0(x)),f(x,w)=g(x,eTl(w+v0(x))).
We now remark that (69) is a quasilinear equation in divergence form. Moreover, it is easy to see that f and σ satisfy the classical growth assumptions and σ satisfies also the classical coercivity condition. Note that
(71)〈σ(x,w,G)-σ(x,w,G′),G-G′〉=(νeTl(w+v0)+κ)|G-G′|2≥κ0|G-G′|2.
Hence, σ is strictly monotonous in the third variable. We then conclude (see for instance [18, Theorem 1.5] or [14, Theorem  8.8 page 311]) that there exists a weak solution v~∈H01(Ω) for (69).
Consequently,vl:=v~+v0 is a weak solution for (Sl), that is, for all w∈H01(Ω),
(72)∫Ω‍(νeTl(vl)+κ)∇vl∇w+∫Ω‍ρu∇eTl(vl)w=∫Ω‍g(x,eTl(vl))w.
By applying Theorem  4.2 page 108 in [15], we obtain vl∈L∞(Ω).
(ii)(73)Letτl:=1h(el)>0,andassumenowthat∥γ2∥Lr≤τl.
With this additional assumption, we are able to obtain a useful estimation for ∥vl∥L∞. Technically, we will detail a method due to Stampacchia. (See Section 3.2 for the notations and for an introduction of the method. Here only the first step of the technique will be developed further). Let s>|ln∥ζ0∥L∞(∂Ω)|; we consider the function ψs=vl-Ts(vl). We have (see Lemma 4) ψs∈H01(Ω)∩L∞(Ω), and by testing (72) with ψs, we obtain
(74)∫Ω‍(νeTl(vl)+κ)|∇ψs|2+∫Ω‍ρu·∇eTl(vl)ψs︸=:I=∫Ω‍g(x,eTl(vl))ψs︸=:II.
(a) We will now evaluate the terms I and II.
The term I is simplified by writing one of its integrand factors, namely,∇eTl(vl)ψs, as a gradient. More precisely we have ∇eTl(vl)ψs=∇ζl, with ζl∈H1(Ω)∩L∞ (see Lemma D.2 in the Appendices). Hence, by applying the divergence formula, we see that I vanishes:
(75)I=∫Ω‍ρu·∇ζl=div.formula-∫Ω‍ζldivρu︸=0+∫∂Ω‍ρζlu·ndσ︸=0by(15)=0.
We next estimate the term II:
(76)II=∫Ω‍g(x,eTl(vl))ψs≤∫As,l+‍g+(x,eTl(vl))ψs︸:=II1+∫As,l-‍g-(x,eTl(vl))ψs︸:=II2.
Remark now that on As,l- we have vl≤-s≤0, which implies that eTl(vl)≤1. Consequently by using the assumption (63), we obtain
(77)II2≤∫Ω‍γψs1As,l-.
The term II1 is majorated as follows:
(78)II1≤by(64)∫As,l+‍(γ1+γ2h(el))ψs≤by(73)∫Ω‍(γ1+1)ψs1As,l+.
(b) At this point by using the estimates (75), (77), and (78) together with (74) and (62), we obtain
(79)κ0∫Ω‍|∇ψs|2≤∫Ω‍(γ1As,l-+(γ1+1)1As,l+)︸=:Eψs.
Note that E∈Lr(Ω) and ∥E∥Lr≤C0, C0=∥γ∥Lr+∥γ1+1∥Lr. On the other hand (see Lemma 5) there exists E~∈(Lβ(Ω))3 satisfying ∥E~∥(Lβ)N≤∥E∥Lr, and ∫Ω‍Eφ=∫Ω‍E~∇φ for all φ∈H01. Recall also that we have assumed in (9) that r>N/2 which implies β>N. By again using the Hölder inequality, we obtain for φ∈H01(80)∫Ω‍E~∇φ≤∥E~∥(Lβ)N∥φ∥H01|Suppφ|(β-2)/2β.
Consequently, (79) leads to:
(81)κ0∥ψs∥H012≤∫Ω‍E~∇ψs≤C0|As,l|(β-2)/2β∥ψs∥H01≤Youngineq.κ02∥ψs∥H012+C022κ0|As,l|(β-2)/β.
Let Φ:=(β-2)/β>(N-2)/N. We have obtained the estimate
(82)∥ψs∥H012≤C~0|As,l|Φ,C~0=C024κ02.
By now using the Stampacchia estimates (see Section 3.2), we obtain the existence of a real Λ independent of l such that |AΛ,l|=0. Hence, ∥v~l∥L∞(Ω)≤Λ, Λ=C(∥ζ0∥L∞(∂Ω),|Ω|,N,r)(C0/κ0).Proof of Proposition9.
We have the following.
(i) Existence and Estimates. Let (vl)l≥1 be the sequence given in Lemma 10. Let also l≥1 be given, and τl:=1/h(el). We assume that ∥γ2∥Lr≤τl. It follows from Lemma 10(ii) that ∥vl∥L∞≤K, where K=K(κ0,DATA) (K independent of l) is the integer defined by K=[Λ]+1.
Let nowτ:=τK, and assume that ∥γ2∥Lr≤τ. Then, we have ∥vK∥L∞≤K. Hence, TK(vK)=vK. On the other hand, vK satisfies (SK), that is
(83)ρu∇evK-div((νevK+κ)∇vK)=g(x,evK)inΩ,vK=lnζ0on∂Ω.
Let ζ:=evK. We have (see Theorem D.1 in the Appendices) ζ∈H1(Ω)∩L∞(Ω) and ∇ζ=ζ∇vK. Consequently, ζ is a solution of problem (S).
Moreover,∥ζ∥L∞(Ω)≤eK=C(κ0,DATA).
On the other hand,vK≥-K almost everywhere implies that ζ≥e-K a.e. in Ω. Hence, we obtain (65) by setting ζmin=e-K and ζmax=eK.
(a) The estimation (66) is obtained as follows.
By using the test functionψs+ instead of ψs, we obtain
(84)ζmax≤Λ1:=eC(∥γ1∥Lr/κ0).
This last estimation is only a first step in order to obtain the majoration for ζmax announced in (66).
In fact, letκ~:=ν+(κ/ζ). We have 0<ν+(κ0/Λ1)=:κ~0≤κ~≤ν+(κ0/ζmin)<∞, and
(85)ρu·∇ζ-div(κ~∇ζ)=g-+g+inΩ,ζ=ζ0on∂Ω
we can then consider the decomposition ζ=ζ1+ζ2, where ζ1 (resp., ζ2) ∈H1(Ω)∩L∞(Ω) satisfies the following problem (S1) (resp., (S2):
(S1)ρu·∇ζ1-div(κ~∇ζ1)=g-inΩζ1=ζ0on∂Ω,(S2)ρu·∇ζ2-div(κ~∇ζ2)=g+inΩζ2=0on∂Ω.
Note that the second member in the PDE in (S1) is negative. Hence, by the maximum principle (see [15, p. 80] or [16, p. 191] for a simplified situation), we obtain ζ1≤∥ζ0∥L∞(∂Ω) a.e. in Ω.
By using next the Stampacchia technique (see again Section3.2, Remark 6), we major the function ζ2 as follows:
(86)ζ2≤C′∥g+∥Lrκ~0≤by(84)eC(∥γ1∥Lr/κ0)∥g+∥Lrk0,C=C(DATA).
This leads to the majoration (66).
(ii) Regularity Results. (a) If we assume that u∈(L∞(Ω))N, then ρu∈(L∞(Ω))N. Moreover, by using the estimates (65) the diffusion coefficient ν+(κ/ζ) is bounded from above and below, and g(x,ζ)∈Lr(Ω) with r>3/2 fixed. Hence (see Lemma 5) g(x,ζ)∈W-1,β(Ω) with β>N. Consequently, by using the De Giorgi-Nash Theorem (see [19, Th. 8.22]) we obtain ζ∈𝒞0,α(Ω¯), for some α>0.
(b) Assume that in addition we have the following:(87)∂Ω,ζ0areofclass𝒞2,α,g∈𝒞0,α(Ω¯×ℝ+),ρu∈(𝒞1,α(Ω¯))N,κ∈𝒞1,α(Ω¯).
We have proved in the previous point that ζ∈𝒞0,α(Ω¯). We now iterate the Schauder estimates as follows. In a first step we see that g(x,ζ) and ν+(κ/ζ) are in 𝒞0,α(Ω¯), and by applying Theorem 2.7 in [20] we obtain ζ∈𝒞1,α(Ω¯). Consequently, we now obtain (see Appendix B in [11]) g(x,ζ),ν+(κ/ζ)∈𝒞1,α(Ω¯), and by using Theorem 2.8 in [20] we finally obtain ζ∈𝒞2,α(Ω¯). Hence, ζ is a classical solution of (S).


### 4.2. Proof of Proposition7

Letn∈ℕ, and let θn, φn be given. We assume that condition (Kn) is satisfied. Recall that this implies in particular φn≤φmax, where φmax>∥b∥L∞(∂Ω) was fixed in (60). Hence, let ε:=φmax-∥b∥L∞(∂Ω)>0.Step 1.
We introduce(88)κ(n):=Cθ(n)φn,gθ(x,u):=C5-C3Fu2+C4Du.
Hence, the first subproblem in (Qn) reads as
(Qn·1)ρu·∇θn+1-div((ν+κ(n)θn+1)∇θn+1)=gθ(x,θn+1)inΩθn+1=aon∂Ω.
Moreover, it is easy to verify that
(89)0<αθφmax︸=:κ0≤κ(n)≤∥Cθ(n)φn∥L∞<∞gθ+(x,u)=C5+C4D+u,gθ-(x,u)=-C3Fu2+C4D-u.
Note also that κ0 is independent of n. Hence, we can apply Proposition 9(i) (take ζ=θn+1,κ=κ(n), g=gθ, γ1=C5, γ2=C4D+, h(t)=t, γ=C3F+C4|D-|, ζ0=a). We obtain the existence of τ0>0 independent of n such that if ∥D+∥Lr≤τ0, then problem (Qn·1) admits at least one weak solution θn+1∈H1(Ω)∩L∞(Ω). Moreover, we have
(90)0<e-Cφmax≤θn+1≤eCφmax<∞,C=C(DATA).Step 2.
Let now(91)κ~(n):=Cφ(n)θn+1,gφ(n)(x,u)=gφ(n)(x)=-φn(C8θn+1-1+C6θn+1F-C7D).
With these notations, the second subproblem in (Qn) reads as
(Qn·2)ρu·∇φn+1-div((ν+κ~(n)φn+1)∇φn+1)=gφ(n)inΩ,φn+1=bon∂Ω.
We verify that
(92)0<αφθmax︸=:κ~0≤κ~(n)≤∥Cφ(n)∥L∞θmin<∞,gφ+=φnC7D+≤φmaxC7D+,gφ-=-φn(C8θn+1-1+C6θn+1F-C7D-),|gφ-|≤by(93)C~(DATA,φmax).
Hence, we can apply the Proposition 9(i) (take now ζ=φn+1, κ=κ~(n), g=gφ(n), γ1=φmaxC7D+, h=0, γ=C~, ζ0=b). Then, we obtain the existence of a weak solution (at this stage there is nomore additional condition needed for ∥D+∥Lr because h=0) φn+1∈H1(Ω)∩L∞(Ω) for problem (Qn·2). Moreover, we have
(93)0<e-C′eCφmax≤φn+1≤eC′eCφmax<∞,
where C, C′ depend on DATA but not on n. Moreover, by using (66), we have
(94)∥φn+1∥L∞(Ω)≤∥b∥L∞(∂Ω)+C′′∥gφ+∥LreCφmax∥φn+1∥L∞(Ω),C′′=C′′(DATA),≤∥b∥L∞(∂Ω)+C′′∥C7∥L∞eCφmax︸≔K(φmax,DATA)∥D+∥Lr∥φn+1∥L∞(Ω).
Assume now that
(95)∥D+∥Lr≤τ:=min(τ0,εφmaxK).
Then, (94) leads to
(96)∥φn+1∥L∞(Ω)≤∥b∥L∞(∂Ω)+εφmax∥φn+1∥L∞(Ω),
and it follows
(97)∥φn+1∥L∞(Ω)≤∥b∥L∞(∂Ω)1-(ε/φmax)≤φmax.Step 3.
If we assume that∥D+∥Lr≤τ, then by using the results established in the previous two steps, we conclude that there exists a solution (θn+1,φn+1)∈(H1∩L∞)2 for problem (Qn). Moreover, this solution satisfies (Kn+1), and (61) holds.


## 4.1. Auxiliary Results

Letn∈ℕ, and let (θn,φn) be given and satisfy (Kn). We want to obtain (θn+1,φn+1) by solving (Qn), and in order to iterate the algorithm we also want that (θn+1,φn+1) satisfies (Kn+1).Remark that the system(Qn) is composed of two coupled scalar elliptic equations in divergence form, with a possible singular and degenerate structure. Hence, the goal of this subsection is to study this last kind of scalar problem.In order to do this, we first introduce a weightκ:ℝ+→ℝ+ which is assumed to be measurable and satisfying
(62)0<κ0≤κ≤κ1a.e.inΩ,
where κ0 and κ1 are two given reals.Let alsog:Ω×ℝ+→ℝ be a Caratheodory function (i.e., for all u∈ℝ+:x→g(x,u) is measurable, and for a.a. x∈Ω:u→g(x,u) is continuous).Let us consider the scalar problem(S) defined in Section 2 after Theorem 1, where ζ0∈H1/2(∂Ω)∩L∞(∂Ω), ζ0≥δ>0 a.e. in δΩ, is given. We always assume that ρ, u, ν, Ω which appear in (S) satisfy their corresponding conditions in (A0).Recall that we allowν≡0 in (A0). Hence, problem (S) may degenerate (i.e., the viscosity ν+(κ/ζ) may vanish) when ζ→∞. Moreover, (S) is singular (i.e., the viscosity tends to infinity) when ζ→0.We want now to find sufficient additional conditions forg that guarantee the existence of a bounded positive weak solution for problem (S). Hence, we shall consider
(63)γ(x):=supu∈[0,1]|g-(x,u)|∈Lr(Ω),(64)g+(x,u)≤γ1(x)+γ2(x)h(u),
where γ1,γ2∈Lr(Ω) and h:ℝ+→ℝ+ is continuous. In fact, more than establishing only the existence of a bounded positive solution for (S), we are interested in obtaining some uniform (with respect to κ1) bounds from above and below and some regularity results. We have the following.Proposition 9.
(i) Letκ satisfy (62), and let g:Ω×ℝ+→ℝ be a Caratheodory function satisfying (63), (64). There exists a real τ>0 depending on κ0, h such that if ∥γ2∥Lr(Ω)≤τ, then there exists a weak solution ζ∈H1(Ω)∩L∞(Ω) for problem (S). Moreover, we have
(65)0<ζmin≤ζ≤ζmax,ζmin=e-C/κ0,ζmax=eC/κ0,
where C depends only on γ, γ1, Ω, r, N, ζ0. In particular ζmin and ζmax are independent of κ1.
In addition, the following extended (wheng+≡0, it is a maximum principle) maximum principle holds:
(66)∥ζ∥L∞(Ω)≤∥ζ0∥L∞(∂Ω)+C~∥g+∥Lrκ0∥ζ∥L∞(Ω),C~=C~(Ω,N,r,γ1).
(ii) Assume that in additionu∈(L∞(Ω))N and ζ0 is Hölder continuous. Then, ζ∈𝒞0,α(Ω¯) for some 0<α<1. Moreover, if ∂Ω is of class 𝒞2,α, g∈𝒞0,α(Ω¯×ℝ+), ρu∈(𝒞1,α(Ω¯))N, κ∈𝒞1,α(Ω¯), and ζ0∈𝒞2,α(∂Ω), then ζ∈𝒞2,α(Ω¯), and it is a classical solution of (S).Before proving Proposition9, we establish an intermediate result. In a first step we consider the change of variable v=lnζ in (S), and for l∈ℕ we introduce a truncated version (Sl) of the system obtained:
(Sl)ρu∇eTl(v)-div((νeTl(v)+κ)∇v)=g(x,eTl(v))inΩ,v=lnζ0on∂Ω.
We then establish the following.Lemma 10.
(i) Letκ satisfy (62), and let g:Ω×ℝ+→ℝ be a Caratheodory function satisfying (63), (64). Then, for any l∈ℕ, there exists a weak solution v=vl∈H1(Ω)∩L∞ for the problem (Sl).
(ii) Let(vl) be the sequence given in (i), and let l≥1 be a fixed integer. Then, there exists τl>0 such that if the function γ2 in (63) satisfies ∥γ2∥Lr≤τl then we have
(67)∥vl∥L∞(Ω)≤C∥γ+γ1∥Lrκ0,C=C(DATA,ζ0).
In particular,C is independent of κ, ν, and l.Proof .
(i) By using the divergence formula, we obtain, for allw∈H1(Ω),
(68)∫Ω‍ρu∇eTl(v)·w=-∫Ω‍eTl(v)div(ρuw)+∫∂Ω‍eTl(v)ρwu·ndσ︸=0=-∫Ω‍eTl(v)ρu·∇w-∫Ω‍eTl(v)div(ρu)︸=0w.
Letv0:=ln(ℛζ0), and consider the change of variable v~:=v-v0. Then, problem (Sl) is equivalent to find v~∈H01(Ω) such that
(69)-divσ(x,v~,∇v~)=f(x,v~)in𝒟′(Ω),
where σ:Ω×ℝ×ℝN→ℝN and f:Ω×ℝ→ℝ are defined by
(70)σ(x,w,G)=(νeTl(w+v0)+κ(x))G-ρ(x)u(x)eTl(w+v0(x)),f(x,w)=g(x,eTl(w+v0(x))).
We now remark that (69) is a quasilinear equation in divergence form. Moreover, it is easy to see that f and σ satisfy the classical growth assumptions and σ satisfies also the classical coercivity condition. Note that
(71)〈σ(x,w,G)-σ(x,w,G′),G-G′〉=(νeTl(w+v0)+κ)|G-G′|2≥κ0|G-G′|2.
Hence, σ is strictly monotonous in the third variable. We then conclude (see for instance [18, Theorem 1.5] or [14, Theorem  8.8 page 311]) that there exists a weak solution v~∈H01(Ω) for (69).
Consequently,vl:=v~+v0 is a weak solution for (Sl), that is, for all w∈H01(Ω),
(72)∫Ω‍(νeTl(vl)+κ)∇vl∇w+∫Ω‍ρu∇eTl(vl)w=∫Ω‍g(x,eTl(vl))w.
By applying Theorem  4.2 page 108 in [15], we obtain vl∈L∞(Ω).
(ii)(73)Letτl:=1h(el)>0,andassumenowthat∥γ2∥Lr≤τl.
With this additional assumption, we are able to obtain a useful estimation for ∥vl∥L∞. Technically, we will detail a method due to Stampacchia. (See Section 3.2 for the notations and for an introduction of the method. Here only the first step of the technique will be developed further). Let s>|ln∥ζ0∥L∞(∂Ω)|; we consider the function ψs=vl-Ts(vl). We have (see Lemma 4) ψs∈H01(Ω)∩L∞(Ω), and by testing (72) with ψs, we obtain
(74)∫Ω‍(νeTl(vl)+κ)|∇ψs|2+∫Ω‍ρu·∇eTl(vl)ψs︸=:I=∫Ω‍g(x,eTl(vl))ψs︸=:II.
(a) We will now evaluate the terms I and II.
The term I is simplified by writing one of its integrand factors, namely,∇eTl(vl)ψs, as a gradient. More precisely we have ∇eTl(vl)ψs=∇ζl, with ζl∈H1(Ω)∩L∞ (see Lemma D.2 in the Appendices). Hence, by applying the divergence formula, we see that I vanishes:
(75)I=∫Ω‍ρu·∇ζl=div.formula-∫Ω‍ζldivρu︸=0+∫∂Ω‍ρζlu·ndσ︸=0by(15)=0.
We next estimate the term II:
(76)II=∫Ω‍g(x,eTl(vl))ψs≤∫As,l+‍g+(x,eTl(vl))ψs︸:=II1+∫As,l-‍g-(x,eTl(vl))ψs︸:=II2.
Remark now that on As,l- we have vl≤-s≤0, which implies that eTl(vl)≤1. Consequently by using the assumption (63), we obtain
(77)II2≤∫Ω‍γψs1As,l-.
The term II1 is majorated as follows:
(78)II1≤by(64)∫As,l+‍(γ1+γ2h(el))ψs≤by(73)∫Ω‍(γ1+1)ψs1As,l+.
(b) At this point by using the estimates (75), (77), and (78) together with (74) and (62), we obtain
(79)κ0∫Ω‍|∇ψs|2≤∫Ω‍(γ1As,l-+(γ1+1)1As,l+)︸=:Eψs.
Note that E∈Lr(Ω) and ∥E∥Lr≤C0, C0=∥γ∥Lr+∥γ1+1∥Lr. On the other hand (see Lemma 5) there exists E~∈(Lβ(Ω))3 satisfying ∥E~∥(Lβ)N≤∥E∥Lr, and ∫Ω‍Eφ=∫Ω‍E~∇φ for all φ∈H01. Recall also that we have assumed in (9) that r>N/2 which implies β>N. By again using the Hölder inequality, we obtain for φ∈H01(80)∫Ω‍E~∇φ≤∥E~∥(Lβ)N∥φ∥H01|Suppφ|(β-2)/2β.
Consequently, (79) leads to:
(81)κ0∥ψs∥H012≤∫Ω‍E~∇ψs≤C0|As,l|(β-2)/2β∥ψs∥H01≤Youngineq.κ02∥ψs∥H012+C022κ0|As,l|(β-2)/β.
Let Φ:=(β-2)/β>(N-2)/N. We have obtained the estimate
(82)∥ψs∥H012≤C~0|As,l|Φ,C~0=C024κ02.
By now using the Stampacchia estimates (see Section 3.2), we obtain the existence of a real Λ independent of l such that |AΛ,l|=0. Hence, ∥v~l∥L∞(Ω)≤Λ, Λ=C(∥ζ0∥L∞(∂Ω),|Ω|,N,r)(C0/κ0).Proof of Proposition9.
We have the following.
(i) Existence and Estimates. Let (vl)l≥1 be the sequence given in Lemma 10. Let also l≥1 be given, and τl:=1/h(el). We assume that ∥γ2∥Lr≤τl. It follows from Lemma 10(ii) that ∥vl∥L∞≤K, where K=K(κ0,DATA) (K independent of l) is the integer defined by K=[Λ]+1.
Let nowτ:=τK, and assume that ∥γ2∥Lr≤τ. Then, we have ∥vK∥L∞≤K. Hence, TK(vK)=vK. On the other hand, vK satisfies (SK), that is
(83)ρu∇evK-div((νevK+κ)∇vK)=g(x,evK)inΩ,vK=lnζ0on∂Ω.
Let ζ:=evK. We have (see Theorem D.1 in the Appendices) ζ∈H1(Ω)∩L∞(Ω) and ∇ζ=ζ∇vK. Consequently, ζ is a solution of problem (S).
Moreover,∥ζ∥L∞(Ω)≤eK=C(κ0,DATA).
On the other hand,vK≥-K almost everywhere implies that ζ≥e-K a.e. in Ω. Hence, we obtain (65) by setting ζmin=e-K and ζmax=eK.
(a) The estimation (66) is obtained as follows.
By using the test functionψs+ instead of ψs, we obtain
(84)ζmax≤Λ1:=eC(∥γ1∥Lr/κ0).
This last estimation is only a first step in order to obtain the majoration for ζmax announced in (66).
In fact, letκ~:=ν+(κ/ζ). We have 0<ν+(κ0/Λ1)=:κ~0≤κ~≤ν+(κ0/ζmin)<∞, and
(85)ρu·∇ζ-div(κ~∇ζ)=g-+g+inΩ,ζ=ζ0on∂Ω
we can then consider the decomposition ζ=ζ1+ζ2, where ζ1 (resp., ζ2) ∈H1(Ω)∩L∞(Ω) satisfies the following problem (S1) (resp., (S2):
(S1)ρu·∇ζ1-div(κ~∇ζ1)=g-inΩζ1=ζ0on∂Ω,(S2)ρu·∇ζ2-div(κ~∇ζ2)=g+inΩζ2=0on∂Ω.
Note that the second member in the PDE in (S1) is negative. Hence, by the maximum principle (see [15, p. 80] or [16, p. 191] for a simplified situation), we obtain ζ1≤∥ζ0∥L∞(∂Ω) a.e. in Ω.
By using next the Stampacchia technique (see again Section3.2, Remark 6), we major the function ζ2 as follows:
(86)ζ2≤C′∥g+∥Lrκ~0≤by(84)eC(∥γ1∥Lr/κ0)∥g+∥Lrk0,C=C(DATA).
This leads to the majoration (66).
(ii) Regularity Results. (a) If we assume that u∈(L∞(Ω))N, then ρu∈(L∞(Ω))N. Moreover, by using the estimates (65) the diffusion coefficient ν+(κ/ζ) is bounded from above and below, and g(x,ζ)∈Lr(Ω) with r>3/2 fixed. Hence (see Lemma 5) g(x,ζ)∈W-1,β(Ω) with β>N. Consequently, by using the De Giorgi-Nash Theorem (see [19, Th. 8.22]) we obtain ζ∈𝒞0,α(Ω¯), for some α>0.
(b) Assume that in addition we have the following:(87)∂Ω,ζ0areofclass𝒞2,α,g∈𝒞0,α(Ω¯×ℝ+),ρu∈(𝒞1,α(Ω¯))N,κ∈𝒞1,α(Ω¯).
We have proved in the previous point that ζ∈𝒞0,α(Ω¯). We now iterate the Schauder estimates as follows. In a first step we see that g(x,ζ) and ν+(κ/ζ) are in 𝒞0,α(Ω¯), and by applying Theorem 2.7 in [20] we obtain ζ∈𝒞1,α(Ω¯). Consequently, we now obtain (see Appendix B in [11]) g(x,ζ),ν+(κ/ζ)∈𝒞1,α(Ω¯), and by using Theorem 2.8 in [20] we finally obtain ζ∈𝒞2,α(Ω¯). Hence, ζ is a classical solution of (S).


## 4.2. Proof of Proposition7

Letn∈ℕ, and let θn, φn be given. We assume that condition (Kn) is satisfied. Recall that this implies in particular φn≤φmax, where φmax>∥b∥L∞(∂Ω) was fixed in (60). Hence, let ε:=φmax-∥b∥L∞(∂Ω)>0.Step 1.
We introduce(88)κ(n):=Cθ(n)φn,gθ(x,u):=C5-C3Fu2+C4Du.
Hence, the first subproblem in (Qn) reads as
(Qn·1)ρu·∇θn+1-div((ν+κ(n)θn+1)∇θn+1)=gθ(x,θn+1)inΩθn+1=aon∂Ω.
Moreover, it is easy to verify that
(89)0<αθφmax︸=:κ0≤κ(n)≤∥Cθ(n)φn∥L∞<∞gθ+(x,u)=C5+C4D+u,gθ-(x,u)=-C3Fu2+C4D-u.
Note also that κ0 is independent of n. Hence, we can apply Proposition 9(i) (take ζ=θn+1,κ=κ(n), g=gθ, γ1=C5, γ2=C4D+, h(t)=t, γ=C3F+C4|D-|, ζ0=a). We obtain the existence of τ0>0 independent of n such that if ∥D+∥Lr≤τ0, then problem (Qn·1) admits at least one weak solution θn+1∈H1(Ω)∩L∞(Ω). Moreover, we have
(90)0<e-Cφmax≤θn+1≤eCφmax<∞,C=C(DATA).Step 2.
Let now(91)κ~(n):=Cφ(n)θn+1,gφ(n)(x,u)=gφ(n)(x)=-φn(C8θn+1-1+C6θn+1F-C7D).
With these notations, the second subproblem in (Qn) reads as
(Qn·2)ρu·∇φn+1-div((ν+κ~(n)φn+1)∇φn+1)=gφ(n)inΩ,φn+1=bon∂Ω.
We verify that
(92)0<αφθmax︸=:κ~0≤κ~(n)≤∥Cφ(n)∥L∞θmin<∞,gφ+=φnC7D+≤φmaxC7D+,gφ-=-φn(C8θn+1-1+C6θn+1F-C7D-),|gφ-|≤by(93)C~(DATA,φmax).
Hence, we can apply the Proposition 9(i) (take now ζ=φn+1, κ=κ~(n), g=gφ(n), γ1=φmaxC7D+, h=0, γ=C~, ζ0=b). Then, we obtain the existence of a weak solution (at this stage there is nomore additional condition needed for ∥D+∥Lr because h=0) φn+1∈H1(Ω)∩L∞(Ω) for problem (Qn·2). Moreover, we have
(93)0<e-C′eCφmax≤φn+1≤eC′eCφmax<∞,
where C, C′ depend on DATA but not on n. Moreover, by using (66), we have
(94)∥φn+1∥L∞(Ω)≤∥b∥L∞(∂Ω)+C′′∥gφ+∥LreCφmax∥φn+1∥L∞(Ω),C′′=C′′(DATA),≤∥b∥L∞(∂Ω)+C′′∥C7∥L∞eCφmax︸≔K(φmax,DATA)∥D+∥Lr∥φn+1∥L∞(Ω).
Assume now that
(95)∥D+∥Lr≤τ:=min(τ0,εφmaxK).
Then, (94) leads to
(96)∥φn+1∥L∞(Ω)≤∥b∥L∞(∂Ω)+εφmax∥φn+1∥L∞(Ω),
and it follows
(97)∥φn+1∥L∞(Ω)≤∥b∥L∞(∂Ω)1-(ε/φmax)≤φmax.Step 3.
If we assume that∥D+∥Lr≤τ, then by using the results established in the previous two steps, we conclude that there exists a solution (θn+1,φn+1)∈(H1∩L∞)2 for problem (Qn). Moreover, this solution satisfies (Kn+1), and (61) holds.


## 5. Proofs of the Theorems

We begin by a lemma.Lemma 11.
Under the assumptions of Proposition7, we can extract a subsequence (still denoted by (θn,φn)) such that
(98)θn⇀*θ,φn⇀*φinL∞(Ω),θn⇀θ,φn⇀φinH1(Ω),(99)1θnφn⟶1θφ,Cθ(n)⟶Cθ(x,θ,φ),Cφ(n)⟶Cφ(x,θ,φ)inLp(Ω),p<∞.Proof.
The first properties in (98) follow directly from Proposition 7. By next using θn-θ0 as test function in (Qn·1) and φn-φ0 as test function in (Qn·2) we obtain a uniform bound for (θn) and (φn) in the H1-norm. Hence, the second properties in (98) follow. Finally, property (99) is obtained by using the dominated convergence theorem. In fact, we have
(100)1θnφn⟶1θφ,Cθ(n)⟶Cθ(x,θ,φ),Cφ(n)⟶Cφ(x,θ,φ)a.einΩ,|1θnφn|≤1θminφmin,|Cθ(n)|≤sup(v,w)∈𝒦Cθ(x,v,w),|Cφ(n)|≤sup(v,w)∈𝒦Cφ(x,v,w),
where 𝒦=[0,θmax]×[0,φmax].


### 5.1. Proof of Theorem1

By using (98) together with (99) we obtain
(101)ρu·∇θn+1⇀ρu·∇θinL1(Ω),(ν+Cθ(n)θn+1φn)∇θn+1⇀(ν+Cθθφ)∇θinLq(Ω),q<2,ρu·∇φn+1⇀ρu·∇φinL1(Ω),(ν+Cφ(n)θn+1φn+1)∇φn+1⇀(ν+Cφθφ)∇φinLq(Ω),q<2.
Moreover by using (98) together with the property θn+1≥θmin>0, we obtain
(102)θn+12⟶θ2,θn+1⟶θ,θn+1-1⟶θ-1inLp(Ω),p<∞.
Hence we can pass to the limit in the approximate problems (Qn). We obtain a weak solution (θ,φ) for problem (Q). That is for all ψ∈𝒟(Ω):
(103)∫Ω‍ρu·∇θψ+∫Ω‍(ν+Cθθφ)∇θ·∇ψ=∫Ω‍(C5-C3Fθ2+C4θD)ψ,∫Ω‍ρu·∇φψ+∫Ω‍(ν+Cφθφ)∇φ·∇ψ=∫Ω‍-φ(C8θ-1+C6θF-C7D)ψ,θ=a,φ=bon∂Ω.
Moreover, this solution satisfies
(104)θ,φ∈H1(Ω)∩L∞(Ω),θ,φ≥min(θmin,φmin)>0a.e.inΩ.


### 5.2. Proof of Theorem2

WhenN=2, the function F has a stronger property of positivity (see Lemma B.1 in the Appendices):
(105)F≥D23.
We will see that this last property allows one to obtain a weak solution for problem (Q) under the assumptions (A0) and (A1) but without assuming a low compressibility condition of the form (23).In order to prove this result, we take over the proof of Proposition7 with slight modifications: if (θn,φn) is given and satisfies (Hn) (it is not useful to consider (Kn) here), then problem (Qn) has at least one solution (θn+1,φn+1) satisfying in addition (Hn+1) and the estimates (61).Step 1.
By using property (105) we major gθ+ as follows:
(106)gθ+(x,u)=C5+C4D+u-C3Fu2≤C5+C4D+u-C33(D+u)2≤C5+D+u︸≥0(C4-C33(D+u))︸≤0ifD+u≥(3C4/C3)≤C5+3C42C3.
Hence, we have here estimated gθ+ independently of the second variable. We then apply Proposition 9(i), but now we take γ1=C5+(3C42/C3) instead of C5 and h(t)=0 instead of h(t)=t. It follows that there exists (without any condition on ∥D+∥Lr because h=0) a weak solution θn+1∈H1(Ω)∩L∞(Ω) for problem (Qn·1), with the estimate
(107)0<e-C∥φn∥L∞≤θn+1≤eC∥φn∥L∞<∞,C=C(DATA).Step 2.
By taking over the arguments presented in the proof of Proposition7 we see that problem (Qn·2) has at least one positive solution φn+1∈H1(Ω)∩L∞(Ω).
Hence, at this point we have obtained a weak solution(θn+1,φn+1) for (Qn) satisfying in addition (Hn+1). It remains to prove that the estimates (61) hold. We have made a first step in this direction by proving (107). We will now prove
(108)0<e-C∥θn+1∥L∞≤φn+1≤∥b∥L∞(∂Ω),C=C(DATA).
In fact, by using the additional assumption (A1) we major the function gφ(n) as follows:
(109)gφ(n)(x)=-φn(C8θn+1-1+C6θn+1F-C7D)≤by(106)-φn(C8θn+1-1+C63θn+1D2-C7D)≤-ρφnθn+1(c8+c63(θn+1D)2-c7(θn+1D))=-ρφnθn+1𝒫(θn+1D),
with 𝒫(X):=(c6/3)X2-c7X+c8. We remark that the discriminant Δ of 𝒫 is negative: Δ=c72-(4/3)c6c8=-4.864*10-3<0. It follows that 𝒫 is positive and consequently gφ(n) is negative. Hence, by applying Proposition 9(i) with now g+≡0 we obtain (108).Step 3.
By using (107) together with (108) we obtain the estimates (61). Hence, we have recovered the conclusions of Proposition 7. The remainder of the proof for Theorem 3 is exactly the same as for Theorem 1: we can extract a subsequence with the properties (98)-(99). These properties are sufficient to pass to the limit n→∞ in (Qn), and Theorem 2 follows.


### 5.3. Proof of Theorem3

Let (θ,φ) be a weak solution for (Q) in the class 𝒮, and consider the notations
(110)gθ(x,u):=C5-C3Fu2+C4Du,gφ(x,u):=-u(C8θ-1+C6θF-C7D).(i)
It suffices to remark that the coefficientsCθ/φ and Cφ/θ are bounded from above and below, and gθ, gφ are Caratheodory functions satisfying (63) and (64). Hence, we can apply the first point in Proposition 9(ii) in each equation of (Q). We obtain θ, φ∈𝒞0,α(Ω¯), for some α>0.(ii)
Assume that in addition we have the following:(111)∂Ω,a,bareofclass𝒞2,α,F∈𝒞0,α(Ω¯),ρu∈(𝒞1,α(Ω¯))N,Cθ,Cφ∈𝒞1,α(Ω¯×(ℝ+)2).
We remark now that the conditions in the second part of Proposition 9(ii) are satisfied for each equation of (Q). Hence, θ,φ∈𝒞2,α(Ω¯), and it is a classical solution of (Q).


## 5.1. Proof of Theorem1

By using (98) together with (99) we obtain
(101)ρu·∇θn+1⇀ρu·∇θinL1(Ω),(ν+Cθ(n)θn+1φn)∇θn+1⇀(ν+Cθθφ)∇θinLq(Ω),q<2,ρu·∇φn+1⇀ρu·∇φinL1(Ω),(ν+Cφ(n)θn+1φn+1)∇φn+1⇀(ν+Cφθφ)∇φinLq(Ω),q<2.
Moreover by using (98) together with the property θn+1≥θmin>0, we obtain
(102)θn+12⟶θ2,θn+1⟶θ,θn+1-1⟶θ-1inLp(Ω),p<∞.
Hence we can pass to the limit in the approximate problems (Qn). We obtain a weak solution (θ,φ) for problem (Q). That is for all ψ∈𝒟(Ω):
(103)∫Ω‍ρu·∇θψ+∫Ω‍(ν+Cθθφ)∇θ·∇ψ=∫Ω‍(C5-C3Fθ2+C4θD)ψ,∫Ω‍ρu·∇φψ+∫Ω‍(ν+Cφθφ)∇φ·∇ψ=∫Ω‍-φ(C8θ-1+C6θF-C7D)ψ,θ=a,φ=bon∂Ω.
Moreover, this solution satisfies
(104)θ,φ∈H1(Ω)∩L∞(Ω),θ,φ≥min(θmin,φmin)>0a.e.inΩ.


## 5.2. Proof of Theorem2

WhenN=2, the function F has a stronger property of positivity (see Lemma B.1 in the Appendices):
(105)F≥D23.
We will see that this last property allows one to obtain a weak solution for problem (Q) under the assumptions (A0) and (A1) but without assuming a low compressibility condition of the form (23).In order to prove this result, we take over the proof of Proposition7 with slight modifications: if (θn,φn) is given and satisfies (Hn) (it is not useful to consider (Kn) here), then problem (Qn) has at least one solution (θn+1,φn+1) satisfying in addition (Hn+1) and the estimates (61).Step 1.
By using property (105) we major gθ+ as follows:
(106)gθ+(x,u)=C5+C4D+u-C3Fu2≤C5+C4D+u-C33(D+u)2≤C5+D+u︸≥0(C4-C33(D+u))︸≤0ifD+u≥(3C4/C3)≤C5+3C42C3.
Hence, we have here estimated gθ+ independently of the second variable. We then apply Proposition 9(i), but now we take γ1=C5+(3C42/C3) instead of C5 and h(t)=0 instead of h(t)=t. It follows that there exists (without any condition on ∥D+∥Lr because h=0) a weak solution θn+1∈H1(Ω)∩L∞(Ω) for problem (Qn·1), with the estimate
(107)0<e-C∥φn∥L∞≤θn+1≤eC∥φn∥L∞<∞,C=C(DATA).Step 2.
By taking over the arguments presented in the proof of Proposition7 we see that problem (Qn·2) has at least one positive solution φn+1∈H1(Ω)∩L∞(Ω).
Hence, at this point we have obtained a weak solution(θn+1,φn+1) for (Qn) satisfying in addition (Hn+1). It remains to prove that the estimates (61) hold. We have made a first step in this direction by proving (107). We will now prove
(108)0<e-C∥θn+1∥L∞≤φn+1≤∥b∥L∞(∂Ω),C=C(DATA).
In fact, by using the additional assumption (A1) we major the function gφ(n) as follows:
(109)gφ(n)(x)=-φn(C8θn+1-1+C6θn+1F-C7D)≤by(106)-φn(C8θn+1-1+C63θn+1D2-C7D)≤-ρφnθn+1(c8+c63(θn+1D)2-c7(θn+1D))=-ρφnθn+1𝒫(θn+1D),
with 𝒫(X):=(c6/3)X2-c7X+c8. We remark that the discriminant Δ of 𝒫 is negative: Δ=c72-(4/3)c6c8=-4.864*10-3<0. It follows that 𝒫 is positive and consequently gφ(n) is negative. Hence, by applying Proposition 9(i) with now g+≡0 we obtain (108).Step 3.
By using (107) together with (108) we obtain the estimates (61). Hence, we have recovered the conclusions of Proposition 7. The remainder of the proof for Theorem 3 is exactly the same as for Theorem 1: we can extract a subsequence with the properties (98)-(99). These properties are sufficient to pass to the limit n→∞ in (Qn), and Theorem 2 follows.


## 5.3. Proof of Theorem3

Let (θ,φ) be a weak solution for (Q) in the class 𝒮, and consider the notations
(110)gθ(x,u):=C5-C3Fu2+C4Du,gφ(x,u):=-u(C8θ-1+C6θF-C7D).(i)
It suffices to remark that the coefficientsCθ/φ and Cφ/θ are bounded from above and below, and gθ, gφ are Caratheodory functions satisfying (63) and (64). Hence, we can apply the first point in Proposition 9(ii) in each equation of (Q). We obtain θ, φ∈𝒞0,α(Ω¯), for some α>0.(ii)
Assume that in addition we have the following:(111)∂Ω,a,bareofclass𝒞2,α,F∈𝒞0,α(Ω¯),ρu∈(𝒞1,α(Ω¯))N,Cθ,Cφ∈𝒞1,α(Ω¯×(ℝ+)2).
We remark now that the conditions in the second part of Proposition 9(ii) are satisfied for each equation of (Q). Hence, θ,φ∈𝒞2,α(Ω¯), and it is a classical solution of (Q).



---

*Source: 940924-2013-08-26.xml*