Singular limit of a two-phase flow problem in porous medium as the air viscosity tends to zero

In this paper we consider a two-phase flow problem in porous media and study its singular limit as the viscosity of the air tends to zero; more precisely, we prove the convergence of subsequences to solutions of a generalized Richards model.


Introduction
Hydrologists have studied air-water flow in soils, mainly using the so-called Richards approximation. At least two hypotheses are physically required for this model to be applicable: the water pressure in the saturated region must be larger than the atmospheric pressure and all the unsaturated regions must have a boundary connected to the surface. However, in many situations, these hypotheses are not satisfied and a more general twophase flow model must be considered. This work explores the limit of this general model as the viscosity of the air tends to zero, which is one of the hypotheses required in the Richards model. To that purpose we prove the existence of a weak solution of the twophase flow problem and prove estimates which are uniform in the air viscosity. In this paper, we assume that the air and water phases are incompressible and immiscible. The geometric domain is supposed to be horizontal, homogeneous and isotropic. Our starting point is the following two-phase flow model, which one can deduce from Darcy's law where u and p are respectively the saturation and the pressure of the water phase, k w and k a are respectively the relative permeabilities of the water and the air phase, µ is the ratio between the viscosity of the air phase and that of the water phase, p c is the capillary pressure, s w is an internal source term for the water phase and s a is an internal source term for the air phase; these source terms are used to represent exchanges with the outside. We suppose in particular that the physical functions k w , k a and p c only depend on the saturation u of the water phase, and that k w (1) = k a (0) = 1. The aim of this paper is the study of the limit of the two-phase flow problem as µ ↓ 0.
The classical Richards model as formulated by the engineers is given by where the properties of capillary pressure p c = p c (u) are describes in hypothesis (H 8 ) below. For the existence and uniqueness of the solution of Richards model together with suitable initial and boundary conditions as well as qualitative properties of the solution and methods for numerical approximations we refer to [1], [6], [10], [11]. In this article, we will show that the singular limit as µ ↓ 0 of the two phase flow problem (T P) has the form (F BP) u t − div(k w (u)∇p) = s w u = 1 or ∇(p + p c (u)) = 0 a.e. in Ω × (0, T ).
We remark that a solution of (R) with u > 0 satisfies (F BP).
This paper is organized as follows. In Section 2 we present a complete mathematical formulation of the problem, and state the main mathematical results, which include a precise formulation of the singular limit problem. We give a sequence of regularized problems in Section 3, and prove the existence of a classical solution. In Section 4 we present a priori estimates, which are uniform in an extra regularization parameter δ and in the air viscosity µ. In Section 5, we let δ ↓ 0 and prove that the solution converges to a solution of the two phase flow problem. We study its limiting behavior as the air viscosity µ tends to zero in Section 6. Finally in Section 7 we propose a finite volume algorithm in a one dimensional context and present a variety of numerical solutions.

Mathematical formulation and main results
We consider the two-phase flow problem where T is a positive constant, Q T := Ω × (0, T ) and where we suppose that In this model, u and p are respectively the saturation and the pressure of the water phase, k w and k a are respectively the mobilities of the water phase and the mobility of the non-water phase and p c is the capillary pressure. We assume in particular that the permeability functions k w , k a and the capillary pressure p c only depend on the saturation u of the water phase. Here, we suppose that the flow of the water phase in the reservoir is driven by an injection term f µ (c)s and an extraction term f µ (u)s where s and s are given space dependent functions, c is the saturation of the injected fluid; if c = 1, only water will be injected, if c = 0, only air will be injected, whereas a mixture of water and air will be injected if 0 < c < 1. The function f µ is the fractional flow of the water phase, namely In particular, we remark that f µ (s) is non decreasing. (2.8) Next we introduce a set of notations, which will be useful in the sequel. and for all s ∈ [0, 1]. This implies in particular that 14) for all ϕ in C := {w ∈ W 2,1 2 (Q T ), w(., T ) = 0 in Ω}.
Our first result, which we prove in Section 3, is the following Theorem 2.2 Suppose that the hypotheses (H 1 ) − (H 9 ) are satisfied, then there exists a weak solution (u µ , p µ ) of Problem (S µ ).
Next we define the discontinuous function χ by as well as the graph The main goal of this paper is to prove the following convergence result, Theorem 2.3 Suppose that the hypotheses (H 1 ) − (H 9 ) are satisfied, then there exists a subsequence ((u µn , p µn )) n∈N of weak solutions of Problem (S µn ) and functions u, p,f such that and (u µn ) n∈N tends to u strongly in L 2 (Q T ), (p µn ) n∈N tends to p weakly in L 2 (0, T ; H 1 (Ω)), as µ n tends to zero and T 0 Ω and Ω p(x, t)dx = 0, for almost every t ∈ (0, T ). (2.18) Formally, u satisfies the following limit problem More precisely the following corollary holds where τ > 0 and Ω t , for t ∈ [τ, T ], are smooth subdomains of Ω and O is a smooth domain of Ω × [τ, T ] and that u = 1 in Finally we remark that another form of the limit problem involves a parabolic equation, which is close to the standard Richards equation. Indeed if we set φ(s) := p c (0) − p c (s) and denote by β the inverse function of φ, the function v := φ(u) is a weak solution of the problem withf ∈ H(β(v)).

Existence of a solution of an approximate problem
Let δ be an arbitrary positive constant. In order to prove the existence of a solution of Problem (S µ ) we introduce a sequence of regularized problems (S µ δ ), namely where u δ 0 , c δ , s δ and s δ are smooth functions such that u δ 0 tends to u 0 in L 2 (Ω) and c δ , s δ and s δ tend respectively to c, s and s in L 2 (Q T ), as δ ↓ 0. In particular we suppose that there exists a positive constant C such that Moreover we suppose that u δ 0 , c δ satisfy Adding up (3.1) and (3.2) we deduce the equation We formulate below an equivalent form of Problem (S µ δ ). To that purpose we define the global pressure, P, by We rewrite the equation (3.1) of Problem (S µ δ ) as is continuous on [0, 1] and differentiable on [0, 1). Multiplying (3.11) by f µ (u µ δ ) and adding the result to (3.12) we deduce that This yields a problem equivalent to (S µ δ ), namely (3.14) In order to prove the existence of a smooth solution of (S µ δ ), we introduce the set where α ∈ (0, 1), and we prove the following result Proof: Let T 1 be the map defined for all V ∈ K by T 1 (V ) = W , where W is the unique solution of the elliptic problem, By standard theory of elliptic system (see [7] Theorem 3.2 p 137), we have that For W solution of (Q 1 V ) we consider T 2 defined by T 2 (W ) =V , whereV is the solution of the parabolic problem, ¿From the standard theory of parabolic equations, we have that Moreover defining by L the parabolic operator arising in (Q 2 W ), namely we remark that (2.8), the property (3.9) of c δ and the fact that, by (3.7), s δ is positive imply that , as m → ∞, which ensures the continuity of the map T . It follows from the Schauder fixed point theorem that there exists a solution (u µ δ , P µ δ ) of (S µ δ ) such that This concludes the proof of Lemma 3.1. Moreover we deduce from Lemma 3.1 the existence of a solution of (S µ δ ), namely

A priori Estimates
In view of (2.8) and (3.23) we deduce the following bounds for all (x, t) ∈ Ω × (0, T ). Next we state some essential a priori estimates.
Lemma 4.1 Let (u µ δ , p µ δ ) be a solution of Problem (S µ δ ). There exists a positive constant C, which only depends on Ω, k w , k a and T such that and Proof: We first prove (4.9). Multiplying (3.11) by P = p µ δ + R µ (u µ δ ) and integrating the result on Q T = Ω × (0, T ) we obtain for all h > 0. Moreover we have by Poincaré-Wirtinger inequality that

Using (3.3) and (4.7), it follows that
which we substitute into (4.14) with h = k w (u m ) 2C 1 to deduce, also in view of (3.7) and (4.5), that Furthermore multiplying (3.1) by p µ δ and (3.2) by p µ δ + p c (u µ δ ), adding up both results and integrating on Q T we obtain where We check below that first term on the left-hand-side of (4.16) and I are bounded. Denoting by P c a primitive of p c we have that Since P c is continuous and u µ δ is bounded this gives Moreover we have using (3.3) and (2.13) that In view of (H 5 ), (3.7), (4.1), (4.2), (4.7) and (4.8) we obtain This together with (4.15) yields I ≤ C 5 C 3 + C 6 . Substituting this into (4.16) and also using (4.17) we obtain that which implies (4.9). In view of (4.3), we also deduce from (4.19) the estimate (4.10). Next we prove (4.11). By the definition (2.9) of g, we obtain from (3.10) that Multiplying (4.20) by f µ (u µ δ ) and subtracting the result from (3.1) we deduce that which we substitute into (4.21) to obtain for all a ∈ [0, 1], so that by the definition (2.11) of Q µ we have ∇D µ (u µ δ ) = p c (u µ δ )∇(f µ (u µ δ )). Substituting this into (4.22), which we have multiplied by p c (u µ δ ), we deduce that ]. (4.24) Multiplying (4.20) by D µ (u µ δ ), adding the result to (4.24) and also using the fact that we deduce that Integrating (4.25) on Q T and using the fact that the definition (4.23) of D µ implies where It follows from (4.1), (4.2), (4.6), (4.8) and (3.7) that |J| ≤ C 8 . Substituting this into (4.26) and also using (4.17) we obtain that Furthermore we remark that which together with (4.27) and the fact that By the definition (2.10) of ζ, we have −∇p c (u µ δ )∇g(u µ δ ) = |∇ζ(u µ δ )| 2 . This together with (4.28) implies (4.11) and (4.12), which in view of (4.4) gives (4.13). This completes the proof of Lemma 4.1.
In what follows we give estimates of differences of space translates of p µ δ and g(u µ δ ). We set for r ∈ I R + sufficiently small: and T 0 Ωr where ξ ∈ I R N and |ξ| ≤ 2r.
Next we estimate differences of time translates of g(u µ δ ). Lemma 4.3 Let (u µ δ , p µ δ ) be a solution of Problem (S µ δ ) then there exists a positive constant C such that for all τ ∈ (0, T ).

Proof: We set
Since g is a non decreasing Lipschitz continuous function with the Lipschitz constant C g we have that where we have used (3.1). Integrating by parts this gives (4.33) Next we estimate the right hand side of (4.32). Using (4.3) we have that Similarly we have that Moreover using (4.1) and (4.2) we obtain from the definition (4.33) of K that This together with (3.7) and the fact that the function g(u µ δ ) is bounded uniformly on µ and δ yields Substituting (4.34), (4.35) and (4.36) into (4.32) we deduce that In view of (4.10) and (4.13) we deduce (4.31), which completes the proof of Lemma 4.3.
Letting δ tend to 0, we deduce from the estimates given in Lemmas 4.1 and 4.2 the existence of a weak solution of Problem (S µ ). More precisely, we have the following result,

4)
where ξ ∈ I R N and |ξ| ≤ 2r. Moreover the following estimate of differences of time translates holds for all τ ∈ (0, T ).
The goal of this section is to prove Theorem 2.3. We first deduce from the estimates (5.2), (5.4) and (5.5) that there exists a couple of functions (u, p) and a subsequence ((u µn , p µn )) n∈N such that (u µn ) n∈N tends to u strongly in L 2 (Q T ) and almost everywhere in Q T , (p µn ) n∈N tends to p weakly in L 2 (0, T ; H 1 (Ω)), as µ n tends to zero. Moreover since there exists a functionf ∈ L 2 (Q T ) with 0 ≤f ≤ 1 and a subsequence (f µn m (u µn m )) nm∈N of (f µn (u µn )) n∈N such that (f µn m (u µn m )) nm∈N tends tof weakly in L 2 (Q T ) as µ nm tends to zero. Moreover we deduce respectively from (5.10) and (5.11) that 0 ≤ u ≤ 1 and that Ω p(x, t)dx = 0, for almost every t ∈ (0, T ), which gives (2.18). As it is done in Section 6 in the proof of (5.1), one can first check that k a (u µn m )(∇p µn m + ∇p c (u µn m )) tends to k a (u)(∇p + ∇p c (u)) weakly in L 2 (Q T ), as µ nm ↓ 0 and then deduce from (5.1) the estimate (2.17). Furthermore letting µ nm tends to zero into (5.12) we obtain, since lim µn m ↓0 f µn m (s) = χ(s) for all s ∈ [0, 1], that which coincides with (2.16) and concludes the proof of Theorem 2.3.

Numerical simulations 7.1 The saturation equation and the numerical algorithm
In this section we present numerical simulations in one space dimension. To that purpose we apply the finite volume method, which we present below. To begin with, we rewrite the equations (2.1) and (2.2) in the case that Ω = (0, 1); this gives for (x, t) ∈ (0, 1) × (0, T ) Adding up both equations and using the boundary conditions (2.4) and (2.5) we obtain Substituting (7.3) into (7.1) yields Moreover we deduce from (7.3) and the definition (2.12) of R µ that ∂ x p µ = −∂ x (R µ (u µ )), so that in view of (2.3) we have In the sequel, we compare numerically the solution u µ of (7.4) with the solution u of the limit equation in the case that u < 1, namely .
Furthermore u µ is given by the line with crosses, p µ g is given by the lines with diam and the limit function u corresponds to the continuous line. We note that, for µ small, the functions u and u µ are very close. Here we only start with water and inject a mixture of water and air. The air immediately invades the whole domain. Figure 1 illustrates the result which we proved in this paper, namely that u µ tends to the solution u of the limit equation (7.6) as µ tends to 0 and moreover that the pressure p µ a is constant. This is indeed the case since u < 1.
Second test case: The case that c = 0, 7 and u 0 (x) = 0, 1 on [0, 1/3] 0, 7 on (1/3, 1] . We obtain the following pictures for t = 0, 01 and for t = 0, 1 respectively The injection of a mixture of water and air (c = 0.7) takes place in a region of low water saturation. We first remark that both functions u µ and u evolve very slowly. Here again we have that u(x, t) < 1 for all (x, t) ∈ (0, 1) × (0, T ) and we remark that the graphs of the two functions u µ and u nearly coincide.