Finite volume scheme for two-phase flows in heterogeneous porous media involving capillary pressure discontinuities
aa r X i v : . [ m a t h . A P ] M a r Finite volume scheme for two-phase flows in heterogeneous porousmedia involving capillary pressure discontinuities
Cl´ement Canc`es ∗† Abstract
We study a one-dimensional model for two-phase flows in heterogeneous media, in which the capillary pressurefunctions can be discontinuous with respect to space. We first give a model, leading to a system of degeneratednonlinear parabolic equations spatially coupled by nonlinear transmission conditions. We approximate thesolution of our problem thanks to a monotonous finite volume scheme. The convergence of the underlyingdiscrete solution to a weak solution when the discretization step tends to 0 is then proven. We also show, underassumptions on the initial data, a uniform estimate on the flux, which is then used during the uniqueness proof.A density argument allows us to relax the assumptions on the initial data and to extend the existence-uniquenessframe to a family of solution obtained as limit of approximations. A numerical example is then given to illustratethe behavior of the model.
MSC subject classification. keywords. capillarity discontinuities, degenerate parabolic equation, finite volume scheme
Introduction
The models of immiscible two-phase flows in porous media are widely used in petroleum engineering in order topredict the positions where oil could be collected. The discontinuities of the physical characteristics due to brutalchange of lithology play a crucial role in the phenomenon of oil trapping, preventing the light hydrocarbons fromreaching the surface. It seems that the discontinuities with respect to the space variable of a particular function,called the capillary pressure, are responsible of the phenomenon of oil-trapping [38, 10].In this paper, we consider one-dimensional two-phase flows in heterogeneous porous media, which are madeof several homogeneous submedia. A simplified model of two-phase flow within this rock is described in the firstsection, leading to the definition of the weak solution. The transmission conditions at the interface between thedifferent submedia are written using the graph formalism introduced in [19] for the connection of the capillarypressures, which is simple to manipulate and allows to deal with any type of discontinuity of the domain, withoutany compatibility constraint, contrary to what occurs in [14] and to a lesser extent in [10, 23].The graph way to connect the capillary pressures at the interfaces is well suited to be discretized by a monotonousFinite Volume scheme. A discretization is proposed in the second section of the paper. Adapting the material fromthe book of Eymard, Gallou¨et & Herbin [25] to our case, it is shown that the discrete solution provided by the schemeconverges, up to a subsequence, to a weak solution as the step of the discretization tends to 0. The monotonicityof the transmission conditions is fundamental for proving the convergence of the scheme.Unfortunately, we are not able to show the uniqueness of the weak solution to the problem, because of the lackof regularity. As it will be shown in the fourth section, supposing that the fluxes are uniformly bounded with regardto space and time is sufficient to claim the uniqueness of the solution. The uniqueness proof is an adaptation of theone given in [19] to the case where the convection is not neglected. Here again, the monotonicity of the transmissionconditions at the interfaces is strongly used.The existence of a bounded flux solution is the topic of Section 3. It is shown that if the initial data is regularenough to ensure that the initial flux is bounded with respect to space, then the flux will remain bounded withrespect to space and to time. Such a result has already been obtained in [19], where a parabolic regularizationof the problem had been introduced. A maximum principle on the flux follows. We also quote [10], in which a BV -estimate is shown on the flux. Since the monotonous schemes introduce some numerical diffusion, a stronganalogy can be done between a uniformly parabolic regularization of the problem and the numerical approximation ∗ ENS Cachan Bretagne, IRMAR, UEB, av. Robert Schuman, 35170 Bruz – France, [email protected] † The author is partially supported by the Groupement MoMaS whole sequence of discrete solutions built usingthe finite volume scheme introduced in Section 2 converges towards the unique SOLA to the problem.Finally, a numerical example is given in Section 6. This example gives an evidence of the entrapment of a certainquantity of oil under the interface.
We consider a one-dimensional heterogeneous porous medium, which is an apposition of homogeneous porous media,representing the different geological layers. The physical properties of the medium only depend on the rock typeand are piece-wise constant.For the sake of simplicity, we only deal with two geological layers of same size. A generalization to an arbitraryfinite number of geological layers would only lead to notation difficulties. In the sequel, we denote by Ω = ( − , = ( − , = (0 ,
1) the two homogeneous layers. The interfacebetween the layers is thus { x = 0 } . T is a positive real value.We consider an incompressible and immiscible oil-water flow through Ω. Writing the conservation of each phase,and using Darcy’s law leads to: for all ( x, t ) ∈ Ω i × (0 , T ), φ i ∂ t u − ∂ x [ µ o,i ( u ) ( ∂ x P o,i − ρ o g )] = 0 , (1) − φ i ∂ t u − ∂ x [ µ w,i ( u ) ( ∂ x P w,i − ρ w g )] = 0 , (2)where φ i ∈ (0 ,
1) is the porosity of the porous media Ω i , u is the oil-saturation (then (1 − u ) is the water-saturation), µ β,i is the mobility of the phase β = w, o , where w stands for water, and o for oil. We denote by P β,i the pressureof the phase β , by ρ β its density, and by g the gravity.Adding (1) and (2) shows that : ∂ x q = 0 , where q = − µ w,i ( u ) ( ∂ x P w,i − ρ w g ) − µ o,i ( u ) ( ∂ x P o,i − ρ o g ) (3)is the total flow-rate. For the sake of simplicity, we suppose that q does not depend on time, even if all the resultspresented below still hold for q ∈ BV (0 , T ), as it is shown in [13, chapter 4].Using (3) in (1) and (2) yields: φ i ∂ t u + ∂ x (cid:18) µ o,i ( u ) µ o,i ( u ) + µ w,i ( u ) q + λ i ( u )( ρ o − ρ w ) g − λ i ( u ) ∂ x ( P o,i − P w,i ) (cid:19) = 0 , (4)where λ i ( u ) = µ o,i ( u ) µ w,i ( u ) µ o,i ( u ) + µ w,i ( u ) . One assumes that the capillary pressure ( P o,i − P w,i ) depends only on the saturation and of the rock type. Moreprecisely, ( P o,i − P w,i ) = π i ( u ), where π i ( u ) is supposed to be an increasing Lipschitz continuous function. Theequation (4) becomes φ i ∂ t u + ∂ x ( f i ( u ) − λ i ( u ) ∂ x π i ( u )) = 0 , (5)where f i ( u ) = µ o,i ( u ) µ o,i ( u ) + µ w,i ( u ) q ( t ) + λ i ( u )( ρ o − ρ w ) g. We do the following assumptions on the functions appearing in the equation.
Assumptions 1
For i = 1 , , one has: . π i is an increasing Lipschitz continuous function;2. µ o,i is an increasing Lipschitz continuous function on [0 , , with µ o,i (0) = 0 ;3. µ w,i is a decreasing Lipschitz continuous function on [0 , , with µ w,i (1) = 0 . Remark 1.1
It is often supposed for such problems that the functions µ β,i are monotonous in a large sense, andthat there exist irreducible saturations s i , S i ∈ (0 , , with s i < S i , such that µ o,i ( u ) = 0 if u ∈ [0 , s i ] , µ w,i ( u ) = 0 if u ∈ [ S i , . If we assume that the functions µ β,i are strictly monotonous on their support, a convenient scaling would allow usto suppose that assumptions 1 are fulfilled. We denote by ϕ i ( s ) = R s λ i ( a ) π ′ i ( a ) da, then (5) can be rewritten φ i ∂ t u + ∂ x ( f i ( u ) − ∂ x ϕ i ( u )) = 0 . (6) Properties 1.1
It follows directly from assumptions 1 that for i = 1 , :1. f i is Lipschitz continuous and f i (0) = 0 , f i (1) = q ;2. λ i is Lipschitz continuous, and λ i (0) = λ i (1) = 0 , λ i ( u ) > if u > ;3. ϕ i is an increasing Lipschitz continuous fulfilling ϕ i (0) = 0 , ϕ ′ i (0) = ϕ ′ i (1) = 0 . We deduce from the properties 1.1 that (6) is a degenerated nonlinear parabolic equation.Let us now focus on the transmission conditions through the interface { x = 0 } . We denote by α i = lim s → π i ( s )and β i = lim s → π i ( s ). We define the monotonous graphs ˜ π i by:˜ π i ( s ) = π i ( s ) if s ∈ (0 , , ( −∞ , α i ] if s = 0 , [ β i , + ∞ ) if s = 1 . (7)Let u i denote the trace of u | Ω i on { x = 0 } (which is supposed to exist for the moment). The trace on { x = 0 } from Ω i of the pressure P β,i of the phase β is still denoted by P β,i . As it is exposed in [23] (see also [19]), thepressure of the phase β can be discontinuous through the interface { x = 0 } in the case where it is missing in theupstream side. This can be written µ β, ( u )( P β, − P β, ) + − µ β, ( u )( P β, − P β, ) + = 0 , β ∈ { o, w } . (8)The conditions (8) have direct consequences on the connection of the capillary pressures through { x = 0 } .Indeed, if 0 < u , u <
1, then the partial pressures P o and P w have both to be continuous, thus the connection ofthe capillary pressures π ( u ) = π ( u ) is satisfied. If u = 0 and 0 < u ≤
1, then P o, ≥ P o, and P w, ≤ P w, ,thus π ( u ) ≤ π (0). The same way, u = 1 and 0 ≤ u < π ( u ) ≥ π (1). Checking that the definitionof the graphs ˜ π and ˜ π implies ˜ π (0) ∩ ˜ π (0) = ∅ , ˜ π (1) ∩ ˜ π (1) = ∅ , we can claim that (8) implies:˜ π ( u ) ∩ ˜ π ( u ) = ∅ . (9)The conservation of each phase leads to the connection of the fluxes on { x = 0 } . Denoting by F i the flux in Ω i , i.e.for all x ∈ Ω i , F i ( x, t ) = f i ( u )( x, t ) − ∂ x ϕ i ( u )( x, t ) , the connection of the fluxes through the interface can be written F (0 , · ) = F (0 , · ) , (10)where (10) has to be understood in a weak sense.We now turn to the problem of the boundary conditions. Because of technical difficulties occurring duringsection 4, we want that the solution to the flow admits bounded fluxes, at least for regular initial data. This willforce us to consider specific boundary conditions, which will involve bounded fluxes.Let G i : ( a, b ) G i ( a, b ) ( i = 1 ,
2) be a function fulfilling the following properties :3 G i is Lipschitz continuous, non-decreasing w.r.t. its first argument, and non-increasing w.r.t. the second. • for all a ∈ [0 , G i ( a, a ) = f i ( a ).Let u, u ∈ L ∞ (0 , T ), 0 ≤ u, u ≤ F ( − , t ) = G ( u ( t ) , u ( − , t )) , F (1 , t ) = G ( u (1 , t ) , u ( t )) . (11)The way in which we approximate the boundary condition shall be judiciously compared with the discretizationof the boundary conditions for scalar hyperbolic conservation laws using monotonous Finite Volume schemes (see[39]).We consider an initial data u ∈ L ∞ (Ω), with 0 ≤ u ≤
1, then we can write the initial-boundary-value problem: φ i ∂ t u + ∂ x [ f i ( u ) − ∂ x ϕ i ( u )] = 0 in Ω i × (0 , T ) ,F (0 , · ) = F (0 , · ) on (0 , T ) , ˜ π ( u ) ∩ ˜ π ( u ) = ∅ on (0 , T ) ,u ( t = 0) = u in Ω ,F ( − , t ) = G ( u ( t ) , u ( − , t )) on (0 , T ) ,F (1 , t ) = G ( u (1 , t ) , u ( t )) on (0 , T ) . ( P )We now define the notion of weak-solution Definition 1.1
A function u is said to be a weak solution to the problem ( P ) if it fulfills:1. u ∈ L ∞ (Ω × (0 , T )) , with ≤ u ≤ ;2. for i = 1 , , ϕ i ( u ) ∈ L (0 , T ; H (Ω i )) ;3. for a.e. t ∈ (0 , T ) , ˜ π ( u ( t )) ∩ ˜ π ( u ( t )) = ∅ , where u i denotes the trace of u | Ω i on { x = 0 } ;4. for all ψ ∈ D (Ω × [0 , T [) , denoting by u (1 , · ) and u ( − , · ) the traces of u on the boundary, Z T X i =1 , Z Ω i φ i u ( x, t ) ∂ t ψ ( x, t ) dxdt + X i =1 , Z Ω i φ i u ( x ) ψ ( x, dx + Z T X i =1 , Z Ω i [ f i ( u )( x, t ) − ∂ x ϕ i ( u )( x, t )] ∂ x ψ ( x, t ) dxdt + Z T G ( u ( t ) , u ( − , t )) ψ ( − , t ) dt − Z T G ( u (1 , t ) , u ( t )) ψ (1 , t ) dt = 0 . (12) In this section, we build an implicit finite volume scheme in order to approximate a solution of ( P ). We will adaptthe convergence proofs stated in [21, 25, 23], which are based on monotonicity properties of the scheme. This willallow us to claim the convergence in L p (Ω × (0 , T )), up to a subsequence, of the discrete solutions built using thefinite volume scheme towards a weak solution to the problem as step of the the discretization tends to 0. We first need to discretize all the data, so that we can define an approximate problem through the finite volumescheme.
Discretization of Ω : for the sake of simplicity, we will only deal with uniform spatial discretizations. Let N ∈ N ⋆ , one defines: ( x j = j/N, ∀ j ∈ [[ − N , N]] ,x j +1 / = j + 1 / N , ∀ j ∈ [[ − N , N − . One denotes by δx = 1 /N . Discretization of (0 , T ) : once again, we will only deal with uniform discretizations. Let M ∈ N ⋆ , one defines:for all n ∈ [[0 , M]], t n = nT /M . One denotes by δt = T /M . We denote by D the discretization of Ω × (0 , T ) deducedof those of Ω and (0 , T ). 4 iscretization of u : ∀ j ∈ [[ − N , N − u , D ( x j +1 / ) = u j +1 / = 1 δx Z x j +1 x j u ( x ) dx. (13) Discretization of the boundary conditions: ∀ n ∈ [[0 , M]], u n +1 = 1 δt Z t n +1 t n u ( t ) dt, u n +1 = 1 δt Z t n +1 t n u ( t ) dt. The Finite Volume scheme: the first equation of ( P ) can be rewritten: φ i ∂ t u + ∂ x F i ( x, t ) = 0 , in Ω × (0 , T )with F i ( x, t ) = f i ( u ) − ∂ x ϕ i ( u ). We consider the following implicit scheme: ∀ j ∈ [[ − N , N − ∀ n ∈ [[0 , M − φ i u n +1 j +1 / − u nj +1 / δt δx + F n +1 j +1 − F n +1 j = 0 (14)where F n +1 j is an approximation of the mean flux through x j on ( t n , t n +1 ), and i is chosen such that ( x j , x j +1 ) ⊂ Ω i .This notation will hold all along the paper. We choose a monotonous discretization of the flux: ∀ j ∈ [[ − N + 1 , − ∪ [[1 , N − ∀ n ∈ [[0 , M − F n +1 j = G i ( u n +1 j − / , u n +1 j +1 / ) − ϕ i ( u n +1 j +1 / ) − ϕ i ( u n +1 j − / ) δx , (15)where G i is the same function as the one defined in (11). We also define F n +1 − N = G ( u n +1 , u − N +1 / ) , F n +1 N = G ( u N − / , u n +1 ) , (16) F n +10 = G ( u n +1 − / , u n +10 , ) − ϕ ( u n +10 , ) − ϕ ( u n +1 − / )) δx (17)= G ( u n +10 , , u n +11 / ) − ϕ ( u n +11 / ) − ϕ ( u n +10 , )) δx , (18)where u n +10 , , u n +10 , moreover satisfy ˜ π ( u n +10 , ) ∩ ˜ π ( u n +10 , ) = ∅ . (19) Remark 2.1
The choice of the boundary conditions F n +1 ± N has been done in order to ensure (cid:12)(cid:12) F n +1 − N (cid:12)(cid:12) ≤ k G k ∞ < ∞ , (cid:12)(cid:12) F n +1 N (cid:12)(cid:12) ≤ k G k ∞ < ∞ . Thanks to the following lemma, such a couple ( u n +10 , , u n +10 , ) is unique in [0 , , thus the discrete transmissionconditions system (17)-(18)-(19) is well posed. Lemma 2.1
For all ( a, b ) ∈ [0 , , there exists a unique couple ( c, d ) ∈ [0 , such that: ( G ( a, c ) − ϕ ( c ) − ϕ ( a )) δx = G ( d, b ) − ϕ ( b ) − ϕ ( d )) δx , ˜ π ( c ) ∩ ˜ π ( d ) = ∅ . (20) Furthermore, ( a, b ) c and ( a, b ) d are continuous and nondecreasing w.r.t. each one of their arguments. Proof
For i = 1 ,
2, ˜ π − i are continuous non-decreasing functions, increasing on [ π i (0) , π i (1)] and constant otherwise. Thenwe can build the continuous non-decreasing function Λ, defined byΛ : ( R → R p G (˜ π − ( p ) , b ) − G ( a, ˜ π − ( p )) + 2 δx (cid:0) ϕ ◦ ˜ π − ( p ) − ϕ ( a ) + ϕ ◦ ˜ π − ( p ) − ϕ ( b ) (cid:1) . p such that Λ( p ) = 0, the couple (˜ π − ( p ) , ˜ π − ( p )) is a solution to the discrete transmission conditions system(17)-(18)-(19). It is easy to check, using the monotonicity of the functions G i that for all p ≤ min π i (0), Λ( p ) ≤ p ≥ max π i (1), Λ( p ) ≥
0. Thus there exists p ⋆ such that Λ( p ⋆ ) = 0.Suppose that there exists i such that p ⋆ ∈ ( π i (0) , π i (1)), then since ϕ i is increasing, Λ is increasing on aneighborhood of p ⋆ , and then the solution to the system (20) is unique.Suppose now that p ⋆ / ∈ S i ( π i (0) , π i (1)). Either p ⋆ ≤ min i π i (0), then c = d = 0, or p ⋆ ≥ max i π i (1), then c = d = 1, or p ⋆ ∈ [ π k (1) , π l (0)] for k = l . We can suppose without any loss of generality that p ⋆ ∈ [ π (1) , π (0)],then the unique solution to the system (20) is given by c = 1, d = 0.To conclude the proof of the lemma, it only remains to check that ( a, b ) Λ is decreasing w.r.t. each one of itsarguments, then the monotonicity of Λ and ˜ π − i ensures that ( a, b ) c and ( a, b ) d are non-decreasing. (cid:3) We will now work on the implicit finite volume scheme given by (13)-(19) to show that this approximate problemis well-posed.
Definition 2.1
Let
N, M be two positive integers and D be the associated discretization of Ω × (0 , T ) . One defines: X D ,i = (cid:26) z ∈ L ∞ (Ω i × (0 , T )) / ∀ ( x j , x j +1 ) ⊂ Ω i , ∀ n ∈ [[0 , M − ,z | ( x j ,x j +1 ) × ( t n ,t n +1 ] is a constant (cid:27) , and X D = (cid:8) z ∈ L ∞ (Ω × (0 , T )) / ∀ i = 1 , , z | Ω i × (0 ,T ) ∈ X D ,i (cid:9) . One defines u D ( x, t ) ∈ X D , called discrete solution , given almost everywhere in ( − , × (0 , T ) by: for all j ∈ [[ − N , N − , for all n ∈ [[0 , M − , ( u D ( x,
0) = u , D ( x ) = u j +1 / if ( x, t ) ∈ ( x j , x j +1 ) ,u D ( x, t ) = u n +1 j +1 / if ( x, t ) ∈ ( x j , x j +1 ) × ( t n , t n +1 ] , where (cid:16) u n +1 j +1 / (cid:17) j,n are given by the scheme (14) . The monotonicity of the flux F n +1 j w.r.t. (cid:16) u n +1 k +1 / (cid:17) k allows us to rewrite the scheme (14) under the form H j +1 / (cid:18) u n +1 j +1 / , u nj +1 / (cid:16) u n +1 k +1 / (cid:17) k = j (cid:19) = 0 , (21)where H j +1 / is continuous, increasing w.r.t. its first argument, and non-increasing w.r.t. all the others. Definition 2.2
A function v D is said to be a discrete supersolution (resp. w D is a discrete subsolution) if itbelongs to X ( D ) , and if it satisfies: ∀ j ∈ [[ − N , N − ,H j +1 / (cid:18) v n +1 j +1 / , v nj +1 / (cid:16) v n +1 k +1 / (cid:17) k = j (cid:19) ≥ , ( resp. H j +1 / (cid:18) w n +1 j +1 / , w nj +1 / (cid:16) w n +1 k +1 / (cid:17) k = j (cid:19) ≤ ). Remark 2.2
A function u D is a discrete solution to the scheme if and only if it is both a supersolution and asubsolution. Remark 2.3
It follows from the definition of the scheme, particularly from the definitions of the discrete boundaryconditions (16) and of the discrete fluxes at the interface (17) - (18) , that the constant function equal to is a discretesupersolution, and that the constant function equal to is a discrete subsolution. We now focus on the existence and the uniqueness of the discrete solution to the scheme. In order to prove theexistence of a discrete solution, we first need an a priori estimate on it.6 emma 2.2
Let u D be a discrete solution to the scheme associated to the initial data u , D , let v D be a discretesupersolution associated to the initial data v , D , then for all t ∈ [0 , T ] , X i =1 , Z Ω i φ i ( u D ( x, t ) − v D ( x, t )) + dx ≤ X i =1 , Z Ω i φ i ( u , D ( x ) − v , D ( x )) + dx. Symmetrically, if w D is a subsolution associated to the initial w , D , X i =1 , Z Ω i φ i ( u D ( x, t ) − w D ( x, t )) − dx ≤ X i =1 , Z Ω i φ i ( u , D ( x ) − w , D ( x )) − dx. Proof
Denoting by a ⊤ b = max( a, b ), and a ⊥ b = min( a, b ), it follows from the monotonicity of the functions H j +1 / impliesthat H j +1 / (cid:18) u n +1 j +1 / , u nj +1 / ⊤ w nj +1 / (cid:16) u n +1 k +1 / ⊤ w n +1 k +1 / (cid:17) k = j (cid:19) ≤ ,H j +1 / (cid:18) w n +1 j +1 / , u nj +1 / ⊤ w nj +1 / (cid:16) u n +1 k +1 / ⊤ w n +1 k +1 / (cid:17) k = j (cid:19) ≤ , where w D is a subsolution. Since u n +1 j +1 / ⊤ w n +1 j +1 / is either equal to u n +1 j +1 / or to w n +1 j +1 / , H j +1 / (cid:18) u n +1 j +1 / ⊤ w n +1 j +1 / , u nj +1 / ⊤ w nj +1 / (cid:16) u n +1 k +1 / ⊤ w n +1 k +1 / (cid:17) k = j (cid:19) ≤ . (22)Thanks to the conservativity of the scheme, subtracting (21) to (22), and summing on j ∈ [[ − N , N − X i =1 , Z Ω i φ i (cid:0) u D ( x, t n +1 ) − w D ( x, t n +1 ) (cid:1) − dx ≤ X i =1 , Z Ω i φ i ( u D ( x, t n ) − w D ( x, t n )) − dx. Since this inequality holds for any n ∈ [[0 , M − ∀ t ∈ [0 , T ], X i =1 , Z Ω i φ i ( u D ( x, t ) − w D ( x, t )) − dx ≤ X i =1 , Z Ω i φ i ( u D ( x, − w D ( x, − dx. (23)The proof of the discrete comparison principle between a discrete solution and a discrete supersolution can beperformed similarly. (cid:3) Let us now state the existence and the uniqueness of the discrete solution.
Proposition 2.3
Let u ∈ L ∞ (Ω) , ≤ u ≤ a.e., then there exists a unique discrete solution u D to the scheme,which furthermore fulfills ≤ u D ≤ a.e.. Moreover, if v stands for another initial data, ≤ v ≤ , approximatedby v , D following (13) , and if we denote by v D the corresponding discrete solution, then the following L -contractionprinciple holds; X i =1 , Z Ω i φ i ( u D ( x, t ) − v D ( x, t )) ± dx ≤ X i =1 , Z Ω i φ i ( u , D ( x ) − v , D ( x )) ± dx, ∀ t ∈ [0 , T ] . Proof
It follows from Remark 2.3 and from Lemma 2.2 that the following L ∞ a priori estimate holds:0 ≤ u D ( x, t ) ≤ , for all t ∈ [0 , T ] , for almost all x ∈ Ω . Thanks to this estimate, mimicking the proof given in [24], we can claim the existence of a discrete solution u D .Suppose that u D and v D are two solutions associated to the initial data u , D and v , D . As it was stressed in theremark 2.2, both u D and v D are both discrete sub- and supersolutions. Then, Lemma 2.2 ensures that the following L -contraction principle holds: X i =1 , Z Ω i φ i ( u D ( x, t ) − v D ( x, t )) ± dx ≤ X i =1 , Z Ω i φ i ( u , D ( x ) − v , D ( x )) ± dx, ∀ t ∈ [0 , T ] . The uniqueness of the discrete solution u D corresponding to the initial data u follows. (cid:3) .3 The L ((0 , T ); H (Ω i )) estimates The current subsection is devoted to the proof of the discrete energy estimate stated in Proposition 2.4. Sincethe discrete solutions are only piecewise constant, we need to introduce discrete semi-norms, which are discreteanalogues to the L ((0 , T ); H (Ω i )) semi-norms. Definition 2.3
Let i = 1 , , one defines the discrete L (0 , T ; H (Ω i )) semi-norms | · | , D ,i on X D ,i by: ∀ z ∈ X D ,i , | z | , D ,i = M − X n =0 δt X j ∈ J int ,i δx (cid:18) z ( x j +1 / , t n +1 ) − z ( x j − / , t n +1 ) δx (cid:19) , where J int , = [[ − N + 1 , − and J int , = [[1 , N − . Proposition 2.4
For i = 1 , , one defines the Lipschitz continuous increasing functions ξ i : s Z s p λ i ( a ) π ′ i ( a ) da. There exists
C > only depending on π i , φ i , T, G i such that: X i =1 , | ξ i ( u D ) | , D ,i ≤ C. This estimate is the discrete analogue to: X i =1 , Z T Z Ω i | ∂ x ξ i ( u )( x, t ) | dxdt ≤ C. In order to prove Proposition 2.4, we will need the following technical lemma. To understand this lemma, firstsuppose that the total flow rate q is 0. Then, roughly speaking, it claims that, in the case where the capillarypressure is discontinuous at the interface, the discrete flux is oriented from the high capillary pressure to the lowcapillary pressure. Suppose now that q = 0. In order to respect the conservation of mass, some fluid will have togo through the interface, but we keep a control on the energy. Lemma 2.5
Let ( a, b ) ∈ [0 , , and let ( c, d ) ∈ [0 , be the unique solution to the system (20), as stated inLemma 2.1, then the following inequality holds: ( π ( c ) − π ( d )) (cid:18) G ( a, c ) + ϕ ( a ) − ϕ ( c ) δx/ (cid:19) = ( π ( c ) − π ( d )) (cid:18) G ( d, b ) + ϕ ( d ) − ϕ ( b ) δx/ (cid:19) ≥ −| q || π ( c ) − π ( d ) | . Proof
In this proof, we suppose that π (0) ≥ π (0) and π (1) ≥ π (1), the other cases do not bring any other difficulties.One has ˜ π ( c ) ∩ ˜ π ( d ) = ∅ , so there are three different cases: • π ( c ) = π ( d ): in this case, one has directly:( π ( c ) − π ( d )) (cid:18) G ( a, c ) + ϕ ( a ) − ϕ ( c ) δx/ (cid:19) = 0 . • π ( d ) < π (0): the relation ˜ π ( c ) ∩ ˜ π ( d ) = ∅ ensures that c = 0, thus it follows from the monotonicity of ϕ and G that ϕ ( a ) ≥ ϕ (0) = 0, and G ( a, ≥ G (0 ,
0) = f (0) = 0. This gives:( π ( c ) − π ( d )) (cid:18) G ( a, c ) + ϕ ( a ) − ϕ ( c ) δx/ (cid:19) ≥ . • π ( c ) > π (1): this implies d = 1. From the monotonicity of ϕ and G , we deduce that ϕ ( b ) ≤ ϕ (1) and G (1 , b ) ≥ G (1 ,
1) = q . This yields( π ( c ) − π ( d )) (cid:18) G ( d, b ) + ϕ ( d ) − ϕ ( b ) δx/ (cid:19) ≥ q | π ( c ) − π ( d ) | . Proof of Proposition 2.4.
First check that the scheme (14) can be rewritten u n +1 j +1 / − u nj +1 / δt δx + (cid:16) F n +1 j +1 − f i ( u n +1 j +1 / ) (cid:17) − (cid:16) F n +1 j − f i ( u n +1 j +1 / ) (cid:17) = 0 . We multiply the previous equation by δtπ i ( u n +1 j +1 / ) and sum on j = − N, N −
1. This leads to A n +1 + B n +1 + C n +1 + D n +1 + E n +1 = 0 , (24)where A n +1 = N − X j = − N φ i π i ( u n +1 j +1 / ) (cid:16) u n +1 j +1 / − u nj +1 / (cid:17) δx ; B n +1 = X j / ∈{− N, ,N } δt π i ( u n +1 j − / ) (cid:16) G i ( u n +1 j − / , u n +1 j +1 / ) − G i ( u n +1 j − / , u n +1 j − / ) (cid:17) − π i ( u n +1 j +1 / ) (cid:16) G i ( u n +1 j − / , u n +1 j +1 / ) − G i ( u n +1 j +1 / , u n +1 j +1 / ) (cid:17) ; C n +1 = δt X j / ∈{− N, ,N } (cid:16) π i ( u n +1 j +1 / ) − π i ( u n +1 j − / ) (cid:17) ϕ i ( u n +1 j +1 / ) − ϕ i ( u n +1 j − / ) δx ; D n +1 = δtF n +10 (cid:16) π ( u n +1 − / ) − π ( u n +11 / ) (cid:17) δx − δtπ ( u n +1 − / ) f ( u n +1 − / ) + δtπ ( u n +11 / ) f ( u n +11 / ) ; E n +1 = δtπ ( u n +1 − N +1 / ) (cid:16) G (cid:16) u n +1 , u n +1 − N +1 / (cid:17) − G (cid:16) u n +1 − N +1 / , u n +1 − N +1 / (cid:17)(cid:17) + δtπ ( u n +1 N − / ) (cid:16) G (cid:16) u n +1 N − / , u n +1 (cid:17) − G (cid:16) u n +1 N − / , u n +1 N − / (cid:17)(cid:17) . Denoting by L G a Lipschitz constant of both G i , E n +1 ≥ − δtL G ( k π k ∞ + k π k ∞ ) . (25)One has F n +10 (cid:16) π ( u n +1 − / ) − π ( u n +11 / ) (cid:17) = F n +10 (cid:16) π ( u n +1 − / ) − π ( u n +10 , ) (cid:17) + F n +10 (cid:0) π ( u n +10 , ) − π ( u n +10 , ) (cid:1) + F n +10 (cid:16) π ( u n +10 , ) − π ( u n +11 / ) (cid:17) . It has been proven in Lemma 2.5 that there exists C depending only on q and π i such that F n +10 (cid:0) π ( u n +10 , ) − π ( u n +10 , ) (cid:1) ≥ C . (26)Using the definition of F n +10 , it is then easy to check that there exists C only depending on G i , q , π i , D n +1 ≥ δtC + δt (cid:16) π ( u n +1 − / ) − π ( u n +10 , ) (cid:17) ϕ ( u n +1 − / ) − ϕ ( u n +10 , ) δx/ δt (cid:16) π ( u n +11 / ) − π ( u n +10 , ) (cid:17) ϕ ( u n +11 / ) − ϕ ( u n +10 , ) δx/ . Since π i is a non-decreasing function, G i : s Z s φ i π i ( a ) da is convex, then: ∀ n ∈ [[0 , M − A n +1 ≥ N − X j = − N (cid:16) G i ( u n +1 j +1 / ) − G i ( u nj +1 / ) (cid:17) δx. (27)We denote by Ψ i ( s ) = Z s π i ( τ ) f ′ i ( τ ) dτ , then an integration by parts leads toΨ i ( b ) − Ψ i ( a ) = π i ( a )( G i ( a, b ) − f i ( a )) − π i ( b ) ( G i ( a, b ) − f i ( b )) − Z ba π ′ i ( s )( f i ( s ) − G i ( a, b )) ds. f i ( s ) = G i ( s, s ), it follows from the monotonicity of G i and π i that Z ba π ′ i ( s )( f i ( s ) − G i ( a, b )) ds ≥ . Thus B n +1 ≥ δt X j / ∈{− N, ,N } (cid:16) Ψ i ( u n +1 j +1 / ) − Ψ i ( u n +1 j − / ) (cid:17) ≥ δt (cid:16) Ψ ( u n +1 − / ) − Ψ ( u n +1 − N +1 / ) + Ψ ( u n +1 N − / ) − Ψ ( u n +11 / ) (cid:17) . So there exists C , only depending on π i , f i such that B n +1 ≥ δtC . (28)Let ξ i : s Z s p λ i ( a ) π ′ i ( a ) da , Cauchy-Schwarz inequality yields: ∀ ( a, b ) ∈ [0 , ,( π i ( a ) − π i ( b ))( ϕ i ( a ) − ϕ i ( b )) ≥ ( ξ i ( a ) − ξ i ( b )) . This ensures that C n +1 ≥ δt X j / ∈{− N, ,N } (cid:16) ξ i ( u n +1 j +1 / ) − ξ i ( u n +1 j − / ) (cid:17) δx ; (29) D n +1 ≥ δtC + δt (cid:16) ξ ( u n +10 , ) − ξ ( u n +1 − / ) (cid:17) δx/ δt (cid:16) ξ ( u n +11 / ) − ξ ( u n +10 , ) (cid:17) δx/ . (30)Summing (24) on n ∈ [[0 , M − C , depending only on T , π i , G i , φ i such that X i =1 , | ξ i ( u D ) | , D ,i + M − X n =0 δt (cid:16) ξ ( u n +10 , ) − ξ ( u n +1 − / ) (cid:17) δx/ (cid:16) ξ ( u n +11 / ) − ξ ( u n +10 , ) (cid:17) δx/ ≤ C. (31) (cid:3) Remark 2.4
The estimate (31) is stronger than the one stated in Proposition 2.4, since it lets also appear somecontributions coming from the interface. They will be useful in the sequel. Indeed, if we denote by u D ,i the trace of ( u D ) | Ω i on the interface { x = 0 } , and if we denote by γ D ,i ( t ) = u n +10 ,i if t ∈ ( nδt, ( n + 1) δt ] , then it follows from (31) that lim δt,δx → k u D ,i − γ D ,i k L p (0 ,T ) = 0 , ∀ p ∈ [1 , ∞ ) . Suppose that u D i converges in L p (0 , T ) towards a function u i , as it will be proven later. Then, we directly obtainthat γ D ,i also converges towards u i . Moreover, for all t > , ˜ π ( γ Dd, ( t )) ∩ ˜ π ( γ Dd, ( t )) = ∅ . Since F = { ( a, b ) ∈ [0 , | ˜ π ( a ) ∩ ˜ π ( b ) = ∅} is a closed set of [0 , , we can claim that ˜ π ( u ) ∩ ˜ π ( u ) = ∅ a.e. in (0 , T ) . Let ( M p ) p ∈ N , ( N p ) p ∈ N be two sequences of positive integers tending to + ∞ . We denote D p the discretizationof Ω × (0 , T ) associated to M p , and N p . The L ∞ -estimate stated in Proposition 2.3 shows that there exists u ∈ L ∞ (Ω × (0 , T )), 0 ≤ u ≤
1, such that, up to a subsequence, u D p → u in the L ∞ (Ω × (0 , T )) weak- ⋆ sense as p → + ∞ .We just need to prove that u D p → u almost everywhere in Ω × (0 , T ) to get the convergence of ( u D p ) towards u in L r (Ω × (0 , T )) for any 1 ≤ r < + ∞ . To apply Kolmogorov criterion (see e.g. [12]) we need some estimates onthe space and time translates of ξ i ( u D ). 10 emma 2.6 (space and time translates estimates) For all η ∈ R , for i = 1 , , one denotes Ω i,η = { x ∈ Ω i / ( x + η ) ∈ Ω i } , then the following estimate holds: k ξ i ( u D )( · + η, · ) − ξ i ( u D )( · , · ) k L (Ω i,η × (0 ,T )) ≤ | ξ i ( u D ) | , D ,i | η | ( | η | + 2 δx ) . (32) One denotes w i, D the function defined almost everywhere by: w i, D ( x, t ) = (cid:26) ξ i ( u D )( x, t ) in Ω i × (0 , T ) , in R \ (Ω i × (0 , T )) . There exists C depending only on π i , φ i , T, G i and C only depending on π i , φ i , T, λ i , G i such that: ∀ η ∈ R , k w i, D ( · + η, · ) − w i, D ( · , · ) k L ( R ) ≤ C η, (33) ∀ τ ∈ (0 , T ) , k w i, D ( · , · + τ ) − w i, D ( · , · ) k L (Ω i × (0 ,T − τ )) ≤ C τ. (34)The previous lemma is in fact a compilation of Lemmata 4.2, 4.3 and 4.6 of [25] adapted to our framework. Theestimates (33) and (34) allows us to use the Kolmogorov compactness criterion on the sequence ( w i, D p ) p ∈ N , andthus, there exists w i ∈ L (Ω i × (0 , T )) such that for almost every ( x, t ) ∈ (Ω i × (0 , T )), ξ i ( u D p )( x, t ) → w i ( x, t ),and then thanks to the L ∞ -estimate 0 ≤ u D p ( x, t ) ≤
1, one can claim that ξ i ( u D p ) → w i in L r (Ω i × (0 , T )), for all r ∈ [1 , + ∞ [. Letting p tend to + ∞ in (32) insures that w i belongs to L (0 , T ; H (Ω i )). Since ξ − i is a continuousfunction, we can identify the limit: w i = ξ i ( u ) . Thus ξ i ( u ) ∈ L (0 , T ; H (Ω i )), and since ϕ i ◦ ξ − i is a Lipschitz function, there exists C depending only on T, π i , φ i , G i , λ i such that: k ϕ i ( u ) k L (0 ,T ; H (Ω i )) ≤ C, (35)and that ξ i ( u D p ) → ξ i ( u ), up to a subsequence, in L r (Ω i × (0 , T )) as p → + ∞ for any r ∈ [1 , + ∞ ). Since ξ i , i = 1 , u D p converges a.e. in Ω × (0 , T ) towards u , and then: u D p → u in the L ∞ (Ω × (0 , T ))-weak- ⋆ sense, (36) u D p → u in L r (Ω i × (0 , T )) , ∀ r ∈ [1 , + ∞ [ . (37)Roughly speaking, the approximation u D obtained via a monotonous finite volume scheme, which introducesnumerical diffusion, is “close” to the approximation u ǫ obtained by adding additional diffusion − ǫ ∆ u ǫ to theproblem. For such a continuous problem, we would have an estimate of type Z T Z Ω i ( ∂ x ξ i ( u ǫ )) dxdt ≤ C ′ , which would lead to the relative compactness of the family ( ξ i ( u ǫ )) ǫ> in L (Ω i × (0 , T )). Then the family ( ξ i ( u ǫ )) ǫ> is also relatively compact in L ((0 , T ); H s (Ω i )) for all s ∈ (1 / , ξ i ( u ǫ )) converge in L (0 , T ). The continuity of ξ − i , and the L ∞ -estimateensure that the traces of u ǫ on the boundary and on the interface converge in L r (0 , T ), for all r ∈ [1 , ∞ ).This sketch has to be modified in order to deal with discrete solutions, which do not belong to L ((0 , T ); H s (Ω i ))for s > /
2. Nevertheless, a convenient estimate on the translates at the boundary, based on the discrete L ((0 , T ); H (Ω i )) estimate stated in Proposition 2.4, leads to the following convergence result, which is provenin the multidimensional case in [17, Proposition 3.10]. Lemma 2.7
Let i = 1 , , and let α ∈ ∂ Ω i . We denote by u α, D p ,i the trace of ( u D p ) | Ω i on { x = α } . Then, one has:for all r ∈ [1 , ∞ ) , u α, D p ,i → u | Ω i ( α, · ) in L r (0 , T ) as p → + ∞ . If we denote by u i ( t ) = u | Ω i (0 , t ), it follows from the remark 2.4 that ˜ π ( u ) ∩ ˜ π ( u ) = ∅ a.e. in (0 , T ) . We cansummarize all the results of this subsection in the following proposition :
Proposition 2.8
Let ( M p ) p , ( N p ) p tend to ∞ as p → ∞ , and let ( D p ) p be the corresponding sequence of discretiza-tions of Ω × (0 , T ) . Let (cid:0) u D p (cid:1) p be the sequence of corresponding discrete solutions to the scheme, then, up to a sub-sequence (still denoted by (cid:0) u D p (cid:1) p ), there exists u ∈ L ∞ (Ω × (0 , T )) , ≤ u ≤ a.e., with ξ i ( u ) ∈ L ((0 , T ); H (Ω i )) ( i = 1 , ) such that: u D p → u a.e. in Ω × (0 , T ) as p → ∞ . oreover, keeping the notations of Lemma 2.7, u − , D p , ( t ) → u ( − , t ) for a.e. t ∈ (0 , T ) as p → ∞ ,u , D p , ( t ) → u (1 , t ) for a.e. t ∈ (0 , T ) as p → ∞ ,u , D p ,i ( t ) → u | Ω i (0 , t ) = u i ( t ) for a.e. t ∈ (0 , T ) as p → ∞ , i = 1 , , and ˜ π ( u ) ∩ ˜ π ( u ) = ∅ almost everywhere in (0 , T ) . We will now achieve the proof of the following result.
Theorem 2.9
Let ( M p ) p ∈ N , ( N p ) p ∈ N be two sequences of positive integers tending to + ∞ , and ( D p ) p ∈ N the asso-ciated sequence of discretizations of Ω × (0 , T ) . Then, up to a subsequence, the sequence (cid:0) u D p (cid:1) p of the discretesolutions converges in L r (Ω × (0 , T )) for all r ∈ [1 , ∞ ) to a weak solution to the problem ( P ) in the sense ofDefinition 1.1. Proof
As it has been seen in Proposition 2.8, the discrete solution u D p converges, up to a subsequence, towards a function u fulfilling all the regularity criteria to be a weak solution. In order to prove the convergence of the subsequence toa weak solution, it only remains to show that the weak formulation (12) is satisfied by the limit u .In order to simplify the proof of convergence of the scheme towards a weak solution, we will use a density result,which is a simple particular case of those stated in [22]. Lemma 2.10
Let a, b ∈ R , a < b , then: { ψ ∈ C ∞ c ([ a, b ]) /ψ ′ ∈ C ∞ c (( a, b )) } is dense in W ,q ( a, b ) , q ∈ [1 , + ∞ [ . This lemma particularly allows us, thanks to a straightforward generalization, to restrict the set of test functions ψ for the weak formulation (12) to T = { ψ ∈ D (Ω × [0 , T )) /∂ x ψ ∈ D (( ∪ i =1 , Ω i ) × [0 , T )) } . Let ψ ∈ T . For j ∈ [[ − N p , N p − n ∈ [[0 , M p − ψ nj +1 / = ψ ( x j +1 / , t n ). Assume that p is largeenough to ensure: ψ n − / = ψ n / , ∀ n ∈ [[0 , M p − , (38) ∀ n ∈ [[0 , M p − , ( ψ n − N +1 / = ψ ( − , t n ) ,ψ nN − / = ψ (1 , t n ) . (39)One has also ψ M p j +1 / = 0 , ∀ j ∈ [[ − N p , N p − . (40)For j ∈ [[ − N p , N p − n ∈ [[0 , M p − ψ nj +1 / δt , and sum on j ∈ [[ − N p , N p − n ∈ [[0 , M p − M p − X n =0 N p − X j = − N p φ i ( u n +1 j +1 / − u nj +1 / ) ψ nj +1 / δx + M p − X n =0 δt N p − X j = − N p ( F n +1 j +1 − F n +1 j ) ψ nj +1 / = 0 , which can be rewritten thanks to (38), (39), (40): M p − X n =0 N p − X j = − N p φ i u n +1 j +1 / ( ψ nj +1 / − ψ n +1 j +1 / ) δx − N p − X j = − N p φ i u j +1 / ψ j +1 / δx + M p − X n =0 δt X j / ∈{− N p , ,N p } F n +1 j ( ψ nj − / − ψ nj +1 / ) − M p − X n =0 δtF n +1 − N ψ ( − , t n ) + M p − X n =0 δtF n +1 N ψ (1 , t n ) = 0 . (41)12sing the definition of F n +1 j , we obtain A p + B p + C p + D p + E p = 0 , (42)with, using the notation ψ n − N − / = ψ n − N +1 / , and ψ nN +1 / = ψ nN − / , A p = M p − X n =0 N p − X j = − N p φ i u n +1 j +1 / ( ψ nj +1 / − ψ n +1 j +1 / ) δx ; B p = − N p − X j = − N p φ i u j +1 / ψ j +1 / δx ; C p = M p − X n =0 δt X j / ∈{− N p , ,N p } G i ( u n +1 j − / , u n +1 j +1 / )( ψ nj − / − ψ nj +1 / ) ; D p = − M p − X n =0 δt N p − X j = − N p ϕ i ( u n +1 j +1 / ) ψ nj +3 / − ψ nj +1 / + ψ nj − / δx ; E p = − M p − X n =0 δtG ( u n +1 , u n +1 − N +1 / ) ψ ( − , t n ) + M p − X n =0 δtG ( u n +1 N − / , u n +1 ) ψ (1 , t n ) . Since h A p : ( x, t ) ψ nj +1 / − ψ n +1 j +1 / δt if ( x, t ) ∈ ( x j , x j +1 ) × ( t n , t n +1 ), converges uniformly towards − ∂ t ψ as p → ∞ ,and since u D converges in L (Ω × (0 , T )) towards u ,lim p →∞ A p = − Z T X i =1 , Z Ω i φ i u ( x, t ) ∂ t ψ ( x, t ) dxdt. (43)Thanks to the convergence in L (Ω) of u D ( x,
0) towards u , we havelim p →∞ B p = − X i =1 , Z Ω i φ i u ( x ) ψ ( x, dx. (44)We have to rewrite C p = C p + C p , with C p = M p − X n =0 δt X j / ∈{− N p , ,N p } f i ( u n +1 j − / )( ψ nj − / − ψ nj +1 / ) ,C p = M p − X n =0 δt X j / ∈{− N p , ,N p } (cid:16) G i ( u n +1 j − / , u n +1 j +1 / ) − f i ( u n +1 j − / ) (cid:17) ( ψ nj − / − ψ nj +1 / ) . Thanks to Proposition 2.8, the quantity C p converges towards − Z T X i =1 , Z Ω i f i ( u )( x, t ) ∂ x ψ ( x, t ) dxdt as p → ∞ .Concerning C p , since (cid:12)(cid:12)(cid:12) G i ( u n +1 j − / , u n +1 j +1 / ) − f i ( u n +1 j − / (cid:12)(cid:12)(cid:12) = (cid:12)(cid:12)(cid:12) G i ( u n +1 j − / , u n +1 j +1 / ) − G i ( u n +1 j − / , u n +1 j − / ) (cid:12)(cid:12)(cid:12) ≤ L G | u n +1 j +1 / − u n +1 j − / | , and since ( x, t ) ψ nj − / − ψ nj +1 / δx on ( x j − / , x j +1 / ) × ( t n , t n +1 ) is uniformly bounded by k ∂ x ψ k ∞ , | C p | ≤ L G k ∂ x ψ k ∞ Z T Z Ω | δu D p ( x, t ) | dxdt, where δu D p ( x, t ) = u n +1 j +1 / − u n +1 j − / on ( x j − / , x j +1 / ) × ( t n , t n +1 ), ( j ∈ [[ − N p + 1 , N p − L ((0 , T ); H (Ω i )) estimates stated in Proposition 2.4, that δu D p tends to 0in L (Ω × (0 , T )) as p → ∞ . Thenlim p →∞ C p = − Z T X i =1 , Z Ω i f i ( u )( x, t ) ∂ x ψ ( x, t ) dxdt. (45)13ince, using Proposition 2.8, ϕ i ( u D p ) tends to ϕ i ( u ) ∈ L ((0 , T ); H (Ω i )) in the L (Ω i × (0 , T ))-topology, one has:lim p →∞ D p = − Z T X i =1 , Z Ω i ϕ i ( u )( x, t ) ∂ xx ψ ( x, t ) dxdt = Z T X i =1 , Z Ω i ∂ x ϕ i ( u )( x, t ) ∂ x ψ ( x, t ) dxdt. (46)The strong convergence of the traces, stated in Proposition 2.8 allows us to claim thatlim p →∞ E p = − Z T G ( u ( t ) , u ( − , t )) ψ ( − , t ) dt + Z T G ( u (1 , t ) , u ( t )) ψ (1 , t ) dt. (47)We can thus take the limit for p → ∞ in (42), and it follows from (43)-(44)-(45)-(46)-(47) that u fulfills the weakformulation (12). (cid:3) In this section, we show that, under some regularity assumptions on the initial data, there exists a solution withbounded fluxes. This existence result is the consequence of some additional estimates on the discrete solution, andwill be necessary to get the uniqueness result of Theorem 4.1.
Definition 3.1
A function u is said to be a bounded-flux solution to the problem ( P ) if:1. u is a weak solution to the problem ( P ) in the sense of Definition 1.1;2. ∂ x ϕ i ( u ) belongs to L ∞ (Ω i × (0 , T )) . In order to get an existence result, we need more regularity on the initial data, as stated below.
Assumptions 2
We assume that:1. ∂ x ϕ i ( u ) ∈ L ∞ (Ω i ) , ≤ u ≤ ;2. ˜ π ( u , ) ∩ ˜ π ( u , ) = ∅ , where u ,i is the trace of u | Ω i on { x = 0 } , Theorem 3.1
Suppose that assumptions 2 are fulfilled. Let ( M p ) p ∈ N , ( N p ) p ∈ N be two sequences of positive integerstending to + ∞ . Let ( u D p ) p ∈ N be the sequence of the associated discrete solutions obtained via the finite volumescheme (14), and let u be an adherence value of the sequence ( u D p ) p ∈ N . Then u is a bounded flux solution to theproblem ( P ) in the sense of Definition 3.1. This particularly ensures the existence of such a bounded-flux solution. All the section 3 will be devoted to the proof of Theorem 3.1. We only need to verify the second point in Defini-tion 3.1, because we have already proven in Theorem 2.9 that u is a weak solution. So the aim of this section isto get the uniform bound on the fluxes. Such an estimate can be found in [19] in the case where the convection isneglected. It is obtained using a thin regular transition layer between Ω and Ω , and a regularization of the initialdata u . This technique was also used in [10] to get a BV -estimate on the fluxes in the case of a non-boundeddomain Ω, and for particular values of the data (which are supposed to be more regular). In this paper, we onlydeal with the discrete solution, which can be seen as a regularization of the solution to the continuous problem ( P ).We extend the definitions of the discrete internal fluxes (15)-(18) to the case n = −
1, i.e. in the time t = 0. Forall j ∈ [[ − N + 1 , N − j = 0, F j = G i ( u j − / , u j +1 / ) − ϕ i ( u j +1 / ) − ϕ i ( u j − / ) δx . (48)Thanks to Lemma 2.1, there exists a unique couple ( u , , u , ) solution to the system F = G ( u − / , u , ) − ϕ ( u , ) − ϕ ( u − / ) δx/ G ( u , , u / ) − ϕ ( u / ) − ϕ ( u , ) δx/ , (49)˜ π ( u , ) ∩ ˜ π ( u , ) = ∅ . (50) Remark 3.1 u , and u , are given by Lemma 2.1, and so they are different of u , and u , . emma 3.2 There exists
C > depending only on u , ϕ i , q such that max j ∈ [[ − N+1 , N − | F j | ≤ C. Proof
Since ϕ i ( u ) is a Lipschitz continuous function, and ϕ − i is continuous, u | Ω i is a continuous function, and thereexists y j +1 / ∈ ( x j , x j +1 ) such that u j +1 / = u ( y j +1 / ). Then (48) can be rewritten F j = G i ( u ( y j − / ) , u ( y j +1 / )) − ϕ i ( u ( y j +1 / )) − ϕ i ( u ( y j − / )) δx . Using the fact that ∂ x ϕ i ( u ) ∈ L ∞ (Ω i ) gives directly: ∀ j ∈ [[ − N + 1 , N − \ { } , (cid:12)(cid:12) F j (cid:12)(cid:12) ≤ max i =1 , k G i k ∞ + 2 max i =1 , (cid:0) k ∂ x ϕ i ( u ) k L ∞ (Ω i ) (cid:1) . (51)The monotony of the transmission conditions ˜ π ( u , ) ∩ ˜ π ( u , ) = ∅ and ˜ π ( u , ) ∩ ˜ π ( u , ) = ∅ implies that either u , ≥ u , and u , ≥ u , , or u , ≤ u , and u , ≤ u , . Assume for example that u , ≥ u , and u , ≥ u , —the other case could be treated similarly— then one deduce from (49) that: G ( u , , u ( y / )) − ϕ ( u ( y / )) − ϕ ( u , ) δx/ ≤ F ≤ G ( u ( y − / ) , u , ) − ϕ ( u , ) − ϕ ( u ( y − / )) δx/ , and so since ϕ i ( u ) is a Lipschitz continuous function, | F | ≤ max i =1 , k G i k ∞ + 2 max i =1 , (cid:0) k ∂ x ϕ i ( u ) k L ∞ (Ω i ) (cid:1) . (cid:3) Proposition 3.3
There exists
C > depending only on u , ϕ i , G i , such that max j ∈ [[ − N+1 , N − (cid:18) max n ∈ [[0 , M]] | F nj | (cid:19) ≤ C. Proof
For all j ∈ [[ − N + 1 , N − \ { } , for all n ∈ [[0 , M − F n +1 j − F nj = (cid:16) G i ( u n +1 j − / , u n +1 j +1 / ) − G i ( u nj − / , u n +1 j +1 / ) (cid:17) + (cid:16) G i ( u nj − / , u n +1 j +1 / ) − G i ( u nj − / , u nj +1 / ) (cid:17) + ϕ i ( u n +1 j − / ) − ϕ i ( u nj − / ) δx ! − ϕ i ( u n +1 j +1 / ) − ϕ i ( u nj +1 / ) δx ! . Thus, using (14) yields F n +1 j − F nj = δtφ i δx G i ( u n +1 j − / , u n +1 j +1 / ) − G i ( u nj − / , u n +1 j +1 / ) u n +1 j − / − u nj − / (cid:0) F n +1 j − F n +1 j − (cid:1) + δtφ i δx G i ( u nj − / , u n +1 j +1 / ) − G i ( u nj − / , u nj +1 / ) u n +1 j +1 / − u nj +1 / (cid:0) F n +1 j +1 − F nj (cid:1) + δtφ i δx ϕ i ( u n +1 j − / ) − ϕ i ( u nj − / ) u n +1 j − / − u nj − / (cid:0) F n +1 j − F n +1 j − (cid:1) − δtφ i δx ϕ i ( u n +1 j +1 / ) − ϕ i ( u nj − / ) u n +1 j +1 / − u nj +1 / (cid:0) F n +1 j +1 − F n +1 j (cid:1) . The monotonicity of the scheme is once again crucial, since it implies that there exist two non-negative values a n +1 j,j +1 , a n +1 j,j − such that (cid:0) a n +1 j,j − + a n +1 j,j +1 (cid:1) F n +1 j − a n +1 j,j − F n +1 j − − a n +1 j,j +1 F n +1 j +1 = F nj . (52)15he monotonicity of the graph transmission condition (19) ensures that either u n +10 , ≥ u n , and u n +10 , ≥ u n , , or u n +10 , ≤ u n , and u n +10 , ≤ u n , . Suppose for example that u n +10 , ≥ u n , and u n +10 , ≥ u n , , the other case beingcompletely symmetrical. F n +10 − F n = (cid:16) G ( u n +1 − / , u n +10 , ) − G ( u n − / , u n , ) (cid:17) + ϕ ( u n +1 − / ) − ϕ ( u n − / ) δx/ ! − ϕ ( u n +10 , ) − ϕ ( u n , ) δx/ ! (53)= (cid:16) G ( u n +10 , , u n +11 / ) − G ( u n , , u n / ) (cid:17) + ϕ ( u n +10 , ) − ϕ ( u n / ) δx/ ! − ϕ ( u n +10 , ) − ϕ ( u n / ) δx/ ! . (54)It follows from (53) and from the monotony of G , ϕ that F n +10 − F n ≤ (cid:16) G ( u n +1 − / , u n , ) − G ( u n − / , u n , ) (cid:17) + ϕ ( u n +1 − / ) − ϕ ( u n − / ) δx/ ! . Similar computations as those done to obtain (52) provide the existence of a non-negative value a n +10 , − such that (cid:0) a n +10 , − (cid:1) F n +10 − a n +10 , − F n +1 − ≤ F n . (55)Considering (54) instead of (53) shows the existence of a non-negative value b n +10 , such that (cid:0) b n +10 , (cid:1) F n +10 − b n +10 , F n +11 ≥ F n . (56)We denote by j n +1max (resp. j n +1min ) the integer such that F n +1 j n +1max = max j ∈ [[ − N , N]] F n +1 j (resp. F n +1 j n +1min = min j ∈ [[ − N , N]] F n +1 j ) . Either j n +1max ∈ {− N, N } , then it follows from the remark 2.1 that max j ∈ [[ − N , N]] F n +1 j ≤ max i =1 , k G i k ∞ , or j n +1max ∈ [[ − N + 1 , N − j F n +1 j = F n +1 j n +1max ≤ F nj n +1max ≤ max j F nj . Similarly, (52) and (55) yield min j F n +1 j = F n +1 j n +1min ≥ F nj n +1min ≥ min j F nj . We obtain a kind of discrete maximum principle on the discrete fluxes, which corresponds to the uniform bound onthe continuous fluxes proven in [19]. It follows from Lemma 3.2 thatmax n ∈ [[0 , M]] (cid:18) max j ∈ [[ − N+1 , N − (cid:12)(cid:12) F n +1 j (cid:12)(cid:12)(cid:19) ≤ max i =1 , k G i k ∞ + 2 max i =1 , (cid:0) k ∂ x ϕ i ( u ) k L ∞ (Ω i ) (cid:1) . (cid:3) Conclusion of proof of Theorem 3.1
Let ( N p ) p ∈ N , ( M p ) p ∈ N be two sequences of positive integers tending to+ ∞ , and let ( u D p ) p ∈ N the sequences of associated discrete solutions. It has been seen in theorem 2.9 that ( u D p ) p tends to a weak solution u in L r (Ω × (0 , T )), for all r ∈ [1 , + ∞ ).Let i = 1 ,
2, let ( x, y ) ∈ Ω i , let t ∈ (0 , T ]. For p large enough,there exists j , j ∈ J int such that x j ≤ x ≤ x j +1 and x j ≤ y ≤ x j +1 , and there exists n such that t ∈ ( t n , t n +1 ]. (cid:12)(cid:12) ϕ i ( u D p )( x, t ) − ϕ i ( u D p )( y, t ) (cid:12)(cid:12) = (cid:12)(cid:12)(cid:12) ϕ i ( u n +1 j +1 / ) − ϕ i ( u n +1 j +1 / ) (cid:12)(cid:12)(cid:12) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) j X j = j +1 ϕ i ( u n +1 j +1 / ) − ϕ i ( u n +1 j − / ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ≤ j X j = j +1 (cid:12)(cid:12)(cid:12) ϕ i ( u n +1 j +1 / ) − ϕ i ( u n +1 j − / ) (cid:12)(cid:12)(cid:12) . (cid:12)(cid:12) ϕ i ( u D p )( x, t ) − ϕ i ( u D p )( y, t ) (cid:12)(cid:12) ≤ j X j = j +1 δx (cid:12)(cid:12)(cid:12) F n +1 j − G i ( u n +1 j − / , u n +1 j +1 / ) (cid:12)(cid:12)(cid:12) . We deduce from Proposition 3.3 that there exists
C >
0, depending only on u , ϕ i , G i such that: (cid:12)(cid:12) ϕ i ( u D p )( x, t ) − ϕ i ( u D p )( y, t ) (cid:12)(cid:12) ≤ j X j = j +1 δxC ≤ C ( | x − y | + 2 δx ) . Letting p tend towards + ∞ , i.e. δx and δt towards 0 gives | ϕ i ( u )( x, t ) − ϕ i ( u )( y, t ) | ≤ C | x − y | . (57)So we deduce from (57) that ∂ x ϕ i ( u ) ∈ L ∞ (Ω i × (0 , T )). (cid:3) This section is devoted to the proof of Theorem 4.1, which is an adaptation of [19, Theorem 5.1] to the case wherethe convection is taken into account.
Theorem 4.1 If u , v are bounded-flux solutions in the sense of Definition 3.1 associated to the initial data u , v , then for all p ∈ [1 , + ∞ [ , u and v belong to C ([0 , T ]; L p (Ω)) , and the following L -contraction principle holds: ∀ t ∈ [0 , T ] , Z Ω i φ i ( u ( x, t ) − v ( x, t )) ± dx ≤ Z Ω i φ i ( u ( x ) − v ( x )) ± dx. This particularly implies the uniqueness of the bounded flux solution to the problem ( P ) Obtaining a L -contraction principle for a nonlinear parabolic equation is classical. We refer for example to [5, 26,34, 20, 31, 33, 11] for the case of homogeneous domains, and for boundary conditions of Dirichlet or Neumann type.We have to adapt the proof of the L -contraction principle to our problem, and thus particularly to the boundaryconditions and to the transmission conditions at the interface.We need to introduce the cut-off functions ρ εα ∈ C , (Ω , R + ) defined by ρ εα ( x ) = (cid:18) ε − | x − α | ε (cid:19) + . Lemma 4.2
For all θ ∈ D + ([0 , T )) , lim inf ε → Z T θ ( t ) X i =1 , Z Ω i (cid:16) sign ± ( u − v )( f i ( u ) − f i ( v )) − ∂ x ( ϕ i ( u ) − ϕ i ( v )) ± (cid:17) ∂ x ρ ε ( x ) dxdt ≥ . Proof
We define the subsets of (0 , T ) E u>v = n t ∈ (0 , T ) (cid:12)(cid:12)(cid:12) u ( t ) > v ( t ) or u ( t ) > v ( t ) o ,E u ≤ v = ( E u>v ) c = n t ∈ (0 , T ) (cid:12)(cid:12)(cid:12) u ( t ) ≤ v ( t ) and u ( t ) ≤ v ( t ) o . Since the trace on { x = 0 } of the function sign + ( u − v )( f i ( u ) − f i ( v )) is equal to 0 for all t ∈ E u ≤ v , it is easy tocheck that lim ε → Z E u ≤ v θ X i =1 , Z Ω i sign + ( u − v )( f i ( u ) − f i ( v )) ∂ x ρ ε dxdt = 0 . Thanks to the fact that the trace of ( ϕ i ( u ) − ϕ i ( v )) + is equal to 0 on the interface, one has also,lim inf ε → Z E u ≤ v θ X i =1 , Z Ω i ∂ x ( ϕ i ( u ) − ϕ i ( v )) + ∂ x ρ ε dxdt ≥ . ε → Z E u ≤ v θ X i =1 , Z Ω i h sign + ( u − v )( f i ( u ) − f i ( v )) − ∂ x ( ϕ i ( u ) − ϕ i ( v )) + i ∂ x ρ ε dxdt ≥ . (58)Since u, v are two weak solutions, subtracting their corresponding weak formulation (12) for the test function ψ ( x, t ) = θ ( t ) ρ ε ( x ) leads to Z T X i =1 , Z Ω i φ i ( u − v ) ρ ε ∂ t θdxdt + X i =1 , Z Ω i φ i ( u − v ) ρ ε θ (0) dx + Z T θ X i =1 , Z Ω i (( f i ( u ) − f i ( v )) − ∂ x ( ϕ i ( u ) − ϕ i ( v ))) ∂ x ρ ε dxdt = 0 . Since ρ ε tends to 0 in L (Ω) as ε →
0, one has:lim ε → Z T X i =1 , Z Ω i φ i ( u − v ) ρ ε ∂ t θdxdt + X i =1 , Z Ω i φ i ( u − v ) ρ ε θ (0) dx = 0 , thus lim ε → Z T θ X i =1 , Z Ω i (( f i ( u ) − f i ( v )) − ∂ x ( ϕ i ( u ) − ϕ i ( v ))) ∂ x ρ ε dxdt = 0 . (59)Thanks to the L ∞ (Ω × (0 , T )) bound on the fluxes, one has (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)Z T θ X i =1 , Z Ω i (( f i ( u ) − f i ( v )) − ∂ x ( ϕ i ( u ) − ϕ i ( v ))) ∂ x ρ ε dxdt (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ≤ C k θ k L (0 ,T ) , then, using a density argument, (59) holds for all θ ∈ L (0 , T ).Replacing θ by θχ E u>v in (59), and splitting the positive and the negative parts a = a + − a − , we obtain Z E u>v θ X i =1 , Z Ω i (cid:0) sign + ( u − v )( f i ( u ) − f i ( v )) − ∂ x ( ϕ i ( u ) − ϕ i ( v )) + (cid:1) ∂ x ρ ε dxdt = Z E u>v θ X i =1 , Z Ω i (cid:0) sign − ( u − v )( f i ( u ) − f i ( v )) − ∂ x ( ϕ i ( u ) − ϕ i ( v )) − (cid:1) ∂ x ρ ε dxdt + r ( ε ) , (60)with lim ε → r ( ε ) = 0 . For almost every t ∈ E u>v , it follows from the monotonicity of the graph relations for the capillary pressure˜ π ( u ) ∩ ˜ π ( u ) = ∅ , ˜ π ( v ) ∩ ˜ π ( v ) = ∅ , that t belongs to E u ≥ v = n t ∈ (0 , T ) (cid:12)(cid:12)(cid:12) u ( t ) ≥ v ( t ) and u ( t ) ≥ v ( t ) o . So we obtain exactly in the same way thatfor (58), thatlim inf ε → Z E u>v θ X i =1 , Z Ω i h sign − ( u − v )( f i ( u ) − f i ( v )) − ∂ x ( ϕ i ( u ) − ϕ i ( v )) − i ∂ x ρ ε dxdt ≥ . It follows directly from (60) thatlim inf ε → Z E u>v θ X i =1 , Z Ω i h sign + ( u − v )( f i ( u ) − f i ( v )) − ∂ x ( ϕ i ( u ) − ϕ i ( v )) + i ∂ x ρ ε dxdt ≥ . (61)Adding (58) and (61) achieves the proof of Lemma 4.2. (cid:3) emma 4.3 For all θ ∈ D + ([0 , T )) , lim inf ε → Z T θ Z Ω (cid:0) sign ± ( u − v ) ( f ( u ) − f ( v )) − ∂ x ( ϕ ( u ) − ϕ ( v )) ± (cid:1) ∂ x ρ ε − dxdt ≥ , (62)lim inf ε → Z T θ Z Ω (cid:0) sign ± ( u − v ) ( f ( u ) − f ( v )) − ∂ x ( ϕ ( u ) − ϕ ( v )) ± (cid:1) ∂ x ρ ε dxdt ≥ . (63) Proof
For the sake of simplicity, we will only provelim inf ε → Z T θ Z Ω (cid:0) sign + ( u − v ) ( f ( u ) − f ( v )) − ∂ x ( ϕ ( u ) − ϕ ( v )) + (cid:1) ∂ x ρ ε dxdt ≥ , but all the steps of the proof can be extended to the other cases. We denote by F u>v and F u ≥ v the subsets of (0 , T )given by F u>v = { t ∈ (0 , T ) | u (1 , t ) > v (1 , t ) } , F u ≤ v = ( F u>v ) c = { t ∈ (0 , T ) | u (1 , t ) ≤ v (1 , t ) } . Let ε >
0. For almost every t ∈ F u ≤ v , one has Z Ω ∂ x ( ϕ ( u )( x, t ) − ϕ ( v )( x, t )) + ∂ x ρ ε ( x ) dx ≤ . Then, using the fact that for almost every t ∈ F u ≤ v , the trace of sign + ( u ( · , t ) − v ( · , t )) ( f ( u )( · , t ) − f ( v )( · , t )) on { x = 1 } is equal to 0,lim inf ε → Z F u ≤ v θ Z Ω (cid:0) sign + ( u − v ) ( f ( u ) − f ( v )) − ∂ x ( ϕ ( u ) − ϕ ( v )) + (cid:1) ∂ x ρ ε dx ≥ . (64)We deduce from the weak formulation that for all θ ∈ D ([0 , T )),lim ε → Z T θ (cid:18)Z Ω ( f ( u ) − ∂ x ϕ ( u )) ∂ x ρ ε dx − G ( u (1 , t ) , u ( t )) (cid:19) dt = 0 . (65)Since the fluxes f ( u ) − ∂ x ϕ ( u ) and f ( v ) − ∂ x ϕ ( v ) belong to L ∞ (Ω × (0 , T )), a density argument, which hasalready been used during the proof of Lemma 4.2, allows us to claim that (65) still holds for any θ ∈ L (0 , T ). So,it particularly holds if we replace θ by θχ F u>v This leads tolim ε → Z F u
Proof of the Theorem 4.1.
First, since u and v are weak solutions to a parabolic equation, they are also entropysolutions (see [26], [20]), and it has been proven in [18] that u and v belong to C ([0 , T ] , L p (Ω)), in the sense thatthere exists ˜ u, ˜ v ∈ C ([0 , T ] , L p (Ω)) such that u = ˜ u , v = ˜ v almost everywhere in Ω × (0 , T ).Let u and v be two weak solutions, then some classical computations, based on the doubling variable techniqueapplied on both the time and the space variable (see e.g. [26], [20]) yields yields that for any ψ ∈ D + (Ω i × [0 , T )) Z T Z Ω i φ i ( u ( x, t ) − v ( x, t )) ± ∂ t ψ ( x, t ) dxdt + Z Ω i φ i ( u ( x ) − v ( x )) ± ψ ( x, dx + Z T Z Ω i sign ± ( u ( x, t ) − v ( x, t )) ( f i ( u )( x, t ) − f i ( v )( x, t )) ∂ x ψ ( x, t ) dxdt − Z T Z Ω i ∂ x ( ϕ i ( u )( x, t ) − ϕ i ( v )( x, t )) ± ∂ x ψ ( x, t ) dxdt ≥ . (68)Let θ ∈ D + ([0 , T )), then summing (68) with respect to i = 1 ,
2, choosing ψ ( x, t ) = θ ( t ) (cid:0) − ρ ε − ( x ) − ρ ε ( x ) − ρ ε ( x ) (cid:1) as test function, and letting ε tend to 0 leads to, thanks to Lemmata 4.2 and 4.3 : Z T ∂ t θ ( t ) X i =1 , Z Ω i φ i ( u ( x, t ) − v ( x, t )) ± dxdt + X i =1 , Z Ω i φ i ( u ( x ) − v ( x )) ± θ (0) dx ≥ . (69)Since u, v belong to C ([0 , T ]; L (Ω)), the relation (69) still holds for any θ ∈ BV (0 , T ) with θ ( T + ) = 0. Let t ∈ [0 , T ],we choose θ = χ [0 ,t ) in (69), obtaining this way the L -contraction and comparison principle stated in the Theorem4.1. (cid:3) We aim in this section to extend the existence-uniqueness result obtained in Theorems 3.1 and 4.1 for any initialdata u ∈ L ∞ (Ω), 0 ≤ u ≤ P ) in such a general case, but we are able to prove the existence and the uniqueness of the solutionobtained as limit of approximation by bounded flux solution. Moreover, this limit is the weak solution obtained viathe convergence of the implicit scheme (14) studied previously. Definition 5.1
A function u is said to be a SOLA (solution obtained as limit of approximation) to the problem ( P )if it fulfils: • u is a weak solution to the problem ( P ), • there exists a sequence ( u ν ) ν ∈ N of bounded flux solutions such that u n → u in C ([0 , T ]; L (Ω)) , as n → + ∞ . heorem 5.1 Let u ∈ L ∞ (Ω) , ≤ u ≤ a.e., then there exists a unique SOLA u to the problem ( P ) in thesense of Definition 5.1.Furthermore, if ( M p ) p ∈ N , ( N p ) p ∈ N are to sequences of positive integers tending to + ∞ , and if ( u D p ) p ∈ N is thecorresponding sequence of discrete solutions, then u D p → u in L r (Ω × (0 , T )) , r ∈ [1 , + ∞ ) . Proof
The set E = n u ∈ L ∞ (Ω) (cid:12)(cid:12)(cid:12) ≤ u ≤ , ∂ x ϕ i ( u ) ∈ L ∞ (Ω i ) , ˜ π ( u , ) ∩ ˜ π ( u , ) = ∅ o is dense in { u ∈ L ∞ (Ω) | ≤ u ≤ } for the L (Ω)-topology. Then we can build a sequence ( u ,ν ) ν ∈ N such thatlim ν →∞ k u ,ν − u k L (Ω) = 0 . Let ( u ν ) ν be the corresponding sequence of bounded flux solutions, then we deduce from the Theorem 4.1 that forall ν, µ ∈ N , ∀ t ∈ [0 , T ] , X i =1 , Z Ω i φ i ( u ν ( x, t ) − u µ ( x, t )) ± dx ≤ X i =1 , Z Ω i φ i ( u ,ν ( x ) − u ,µ ( x )) ± dx. (70)Then ( u ν ) ν is a Cauchy sequence in C ([0 , T ]; L (Ω)), thus it converges towards u ∈ C ([0 , T ]; L (Ω)), and ∀ t ∈ [0 , T ] , X i =1 , Z Ω i φ i ( u ν ( x, t ) − u ( x, t )) ± dx ≤ X i =1 , Z Ω i φ i ( u ,ν ( x ) − u ( x )) ± dx. (71)Let us now check that u is a weak solution. Since ϕ i is continuous, and since 0 ≤ u ν ≤ ϕ i ( u ν )converges in L (Ω i × (0 , T )) towards ϕ i ( u ). The L ((0 , T ); H (Ω i )) estimate (35) does not depend on u , thus,up to a subsequence, ( ϕ i ( u ν )) ν converges weakly to ϕ i ( u ) in L ((0 , T ); H (Ω i )). It also converges strongly in L ((0 , T ); H s (Ω i )) for all s ∈ (0 , ϕ i ( u ν )) ν onthe interface. Since ϕ − i is continuous, we obtain the strong convergence of the traces of ( u ν ) ν . Checking that theset F = { ( a, b ) ∈ [0 , | ˜ π ( a ) ∩ ˜ π ( b ) = ∅} is closed in [0 , , the limits u i fulfill ˜ π ( u ) ∩ ˜ π ( u ) = ∅ , and so u is a weak solution, then it is a SOLA.If u and v are two SOLAs associated to the initial data u and v , we can easily prove, using the Theorem 4.1that ∀ t ∈ [0 , T ] , X i =1 , Z Ω i φ i ( u ( x, t ) − v ( x, t )) ± dx ≤ X i =1 , Z Ω i φ i ( u ( x ) − v ( x )) ± dx. (72)The uniqueness particularly follows.Let u ∈ L ∞ (Ω), 0 ≤ u ≤
1, and let ( u ,ν ) ν ⊂ E a sequence of approximate initial data tending to u in L (Ω).We denote by u the unique SOLA associated to u , and by ( u ν ) ν the bounded flux solutions associated to ( u ,ν ) ν .Let ( M p ) p ∈ N , ( N p ) p ∈ N be two sequences of positive integers tending to + ∞ . Let p ∈ N , ν ∈ N , let u D p the discretesolution corresponding to u , and let u ν, D p the discrete solution corresponding to u ,ν . k u D p − u k L (Ω × (0 ,T )) ≤ k u D − u ν, D p k L (Ω × (0 ,T )) + k u ν, D p − u ν k L (Ω × (0 ,T )) + k u ν − u k L (Ω × (0 ,T )) . From the discrete L -contraction principle (23), and from the continuous one (71), we have k u D p − u k L (Ω × (0 ,T )) ≤ T max φ i min φ i k u , D − u ,ν, D p k L (Ω) + k u ν, D p − u ν k L (Ω × (0 ,T )) + T max φ i min φ i k u ,ν − u k L (Ω) . Letting p tend to ∞ , it follows from the definition of ( u ,ν, D ) (adapted from (13)) thatlim p →∞ k u , D − u ,ν, D p k L (Ω) = k u ,ν − u k L (Ω) . We have proven in the Theorem 3.1 that the sequence of discrete solutions converges, under assumption on theinitial data to the unique bounded flux solution, thuslim p →∞ k u ν, D p − u ν k L (Ω × (0 ,T )) = 0 . Figure 1: Saturation profiles for t = 20 , t = 100 , t = 200This implies lim sup p →∞ k u D p − u k L (Ω × (0 ,T )) ≤ T max φ i min φ i k u ,ν − u k L (Ω) . Letting ν tend to ∞ provides lim p →∞ k u D p − u k L (Ω × (0 ,T )) = 0 . The convergence occurs in L (Ω × (0 , T )), but the uniform bound on the sequence (cid:0) u D p (cid:1) in L ∞ (Ω × (0 , T )) ensuresthat the convergence also take place in all the L p (Ω × (0 , T )), for p ∈ [1 , ∞ ). (cid:3) In order to illustrate this model, we use a test case developed by Anthony Michel [32]. The porous medium Ω = (0 , x ∈ (0 , . ∪ (0 . , x ∈ (0 . , . First case:
The total flow rate is equal to 0, since f sand (1) = f shale (1) = 0, and the convection is the exclusive of the volumemass difference between the oil, which is lighter, and the water. The convection functions are given by: f sand ( u ) = 100 ∗ f shale ( u ) = 50 ∗ u (1 − u )1 − u + 2 u . The capillary pressures are first given by π sand ( u ) = u , π shale ( u ) = 0 . u . The function ϕ sand and ϕ shale , given by ϕ sand ( u ) = 10 ∗ Z u s (1 − s )1 − s + 2 s π ′ sand ( s ) ds, ϕ shale ( u ) = 0 . ∗ Z u s (1 − s )1 − s + 2 s π ′ shale ( s ) ds, are computed using an approximate integration formula. The initial data u is equal to 0, and u = 0 . u = 0.The convection is approximated by a Godunov scheme, defined by G i ( a, b ) = min s ∈ [ a,b ] f i ( s ) if a ≤ b, max s ∈ [ b,a ] f i ( s ) otherwise.A little quantity of oil enters the domain from the left boundary condition, and it moves forward in the firstpart made of sand. The discontinuity of the capillary pressure (figure 2) stops the migration of oil, which beginsto collect at the left of the interface, as shown on the figure 1. One can check on the figure 3 that for t smallenough, the oil-flux through the interface { x = 0 . } is equal to 0. The accumulation of oil at the left of { x = 0 . } implies an increase of the capillary pressure. As soon as the capillary pressure connects at { x = 0 . } , the oil canflow through the shale. The next discontinuity at { x = 0 . } does not impede the progression of the oil, since thecapillary pressure force, oriented from the large pressure to the small pressure (here from the left to the right),works in the same direction that the buoyancy, which drives the migration of oil.22 Figure 2: Capillary pressure profiles for t = 20 , t = 100 , t = 200 −5 flux t=20 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1012345x 10 −5 flux t=100 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1012345x 10 −5 flux t=200 Figure 3: Oil-flux profiles for t = 20 , t = 100 , t = 200For t = 200, the presented solution is a steady solution, with constant flux (figure 3). Some oil remains blockedin the first subdomain (0 , . u ( t ) = 0 for t ≥ { x = 0 . } and leave the porous medium (0 , u s ( x ) = (cid:26) x / ∈ (0 . , . ,π − (5( x − . x ∈ (0 . , . u = 0. It is easy to check that u ( · , ≥ u s , thus the comparison principlestated in the Theorem 4.1 ensures that for all t ≥ u ( · , t ) ≥ u s . Thus for all t ≥ Z . u ( x, t ) dx ≥ Z . u s ( x ) dx > . This quantity is said to be trapped by the geology change. Further illustrations, and a scheme comparison will begiven in [32].
Second case:
We only change the values of the capillary pressure functions (and also the linked functions ϕ sand and ϕ shale ). Theamplitude of the variation of each function is reduced from 1 to 0 .
2, i.e. π sand ( u ) = 0 . ∗ u , π shale ( u ) = 0 . . ∗ u . The graph transmission condition for the capillary pressure turns to(1 − u sand ) u shale = 0 , where u sand (resp u shale ) denotes the trace of the oil saturation at the interfaces { x = 0 . } and { x = 0 . } . In thiscase, no oil can overpass the first interface, which is thus impermeable for oil. The only steady solution is u s ( x ) = (cid:26) x < . , x > . . An asymptotic study for capillary pressures tending to functions depending only of space, and not on the saturationhas been performed in [13, Chapter 5&6] (see also [15, 16]). It has been proven that either the limit solution for thesaturation is an entropy solution for th e hyperbolic scalar conservation law with discontinuous fluxes in the senseof [36, 37, 35, 1, 2, 3, 4, 6, 9, 8, 7, 27] (see also [28, 29, 30]), mainly when the capillary forces at the interface areoriented in the same direction that the gravity forces, or that non-classical shocks can occur at the interfaces whenthe capillary forces and the gravity are oriented in opposite directions.23
Figure 4: Saturation profiles for t = 100 , t = 500 , t = 900 Figure 5: Capillary pressure profiles for t = 100 , t = 500 , t = 900 Acknowledgements.
The author would like to acknowledge the Professor Thierry Gallou¨et for his numerousrecommendations and Anthony Michel from IFP for the fruitful discussions on the models. He also thanks AlicePivan for her help with the English language.
References [1] Adimurthi and G. D. Veerappa Gowda. Conservation law with discontinuous flux.
J. Math. Kyoto Univ. ,43(1):27–70, 2003.[2] Adimurthi, J. Jaffr´e, and G. D. Veerappa Gowda. Godunov-type methods for conservation laws with a fluxfunction discontinuous in space.
SIAM J. Numer. Anal. , 42(1):179–208 (electronic), 2004.[3] Adimurthi, Siddhartha Mishra, and G. D. Veerappa Gowda. Optimal entropy solutions for conservation lawswith discontinuous flux-functions.
J. Hyperbolic Differ. Equ. , 2(4):783–837, 2005.[4] Adimurthi, Siddhartha Mishra, and G. D. Veerappa Gowda. Existence and stability of entropy solutions for aconservation law with discontinuous non-convex fluxes.
Netw. Heterog. Media , 2(1):127–157 (electronic), 2007.[5] H. W. Alt and S. Luckhaus. Quasilinear elliptic-parabolic differential equations.
Math. Z. , 183(3):311–341,1983. −5 flux t=100 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 100.511.522.533.544.55x 10 −5 flux t=500 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 100.511.522.533.544.55x 10 −5 flux t=900 Figure 6: Oil-flux profiles for t = 100 , t = 500 , t = 900246] F. Bachmann. Analysis of a scalar conservation law with a flux function with discontinuous coefficients. Adv.Differential Equations , 9(11-12):1317–1338, 2004.[7] F. Bachmann.
Equations hyperboliques scalaires `a flux discontinu . PhD thesis, Universit´e Aix-Marseille I, 2005.[8] F. Bachmann. Finite volume schemes for a non linear hyperbolic conservation law with a flux function involvingdiscontinuous coefficients.
Int. J. Finite Volumes , 3, 2006.[9] F. Bachmann and J. Vovelle. Existence and uniqueness of entropy solution of scalar conservation laws witha flux function involving discontinuous coefficients.
Comm. Partial Differential Equations , 31(1-3):371–395,2006.[10] M. Bertsch, R. Dal Passo, and C. J. van Duijn. Analysis of oil trapping in porous media flow.
SIAM J. Math.Anal. , 35(1):245–267 (electronic), 2003.[11] D. Blanchard and A. Porretta. Stefan problems with nonlinear diffusion and convection.
J. DifferentialEquations , 210(2):383–428, 2005.[12] H. Br´ezis.
Analyse Fonctionnelle: Th´eorie et applications . Masson, 1983.[13] C. Canc`es. ´Ecoulements diphasiques en milieux poreux h´et´erog`enes : mod´elisation et analyse des effets li´es auxdiscontinuit´es de la pression capillaire . PhD thesis, Universit´e de Provence, 2008.[14] C. Canc`es. Nonlinear parabolic equations with spatial discontinuities.
NoDEA Nonlinear Differential EquationsAppl. , 15(4-5):427–456, 2008.[15] C. Canc`es. Asymptotic behavior of two-phase flows in heterogeneous porous media for capillarity dependingonly of the space. I. Convergence to an entropy solution. submitted, arXiv:0902.1877, 2009.[16] C. Canc`es. Asymptotic behavior of two-phase flows in heterogeneous porous media for capillarity dependingonly of the space. II. Occurrence of non-classical shocks to model oil-trapping. submitted, arXiv:0902.1872,2009.[17] C. Canc`es. Immiscible two-phase flows within porous media: Modeling and numerical analysis. In
Analyticaland Numerical Aspects of Partial Differential Equations . de Gruyter, Emmrich & Wittbold edition, 2009. Toappear.[18] C. Canc`es and T. Gallou¨et. On the time continuity of entropy solutions. arXiv:0812.4765v1, 2008.[19] C. Canc`es, T. Gallou¨et, and A. Porretta. Two-phase flows involving capillary barriers in heterogeneous porousmedia. to appear in Interfaces Free Bound., 2009.[20] J. Carrillo. Entropy solutions for nonlinear degenerate problems.
Arch. Ration. Mech. Anal. , 147(4):269–361,1999.[21] C. Chainais-Hillairet. Finite volume schemes for a nonlinear hyperbolic equation. Convergence towards theentropy solution and error estimate.
M2AN Math. Model. Numer. Anal. , 33(1):129–156, 1999.[22] J. Droniou. A density result in Sobolev spaces.
J. Math. Pures Appl. (9) , 81(7):697–714, 2002.[23] G. Ench´ery, R. Eymard, and A. Michel. Numerical approximation of a two-phase flow in a porous mediumwith discontinuous capillary forces.
SIAM J. Numer. Anal. , 43(6):2402–2422, 2006.[24] R. Eymard, T. Gallou¨et, M. Ghilani, and R. Herbin. Error estimates for the approximate solutions of anonlinear hyperbolic equation given by finite volume schemes.
IMA J. Numer. Anal. , 18(4):563–594, 1998.[25] R. Eymard, T. Gallou¨et, and R. Herbin. Finite volume methods. Ciarlet, P. G. (ed.) et al., in Handbook ofnumerical analysis. North-Holland, Amsterdam, pp. 713–1020, 2000.[26] G. Gagneux and M. Madaune-Tort. Unicit´e des solutions faibles d’´equations de diffusion-convection.
C. R.Acad. Sci. Paris S´er. I Math. , 318(10):919–924, 1994.[27] J. Jimenez. Some scalar conservation laws with discontinuous flux.
Int. J. Evol. Equ. , 2(3):297–315, 2007.2528] K. H. Karlsen, N. H. Risebro, and J. D. Towers. On a nonlinear degenerate parabolic transport-diffusionequation with a discontinuous coefficient.
Electron. J. Differential Equations , pages No. 93, 23 pp. (electronic),2002.[29] K. H. Karlsen, N. H. Risebro, and J. D. Towers. Upwind difference approximations for degenerate parabolicconvection-diffusion equations with a discontinuous coefficient.
IMA J. Numer. Anal. , 22(4):623–664, 2002.[30] K. H. Karlsen, N. H. Risebro, and J. D. Towers. L stability for entropy solutions of nonlinear degenerateparabolic convection-diffusion equations with discontinuous coefficients. Skr. K. Nor. Vidensk. Selsk. , (3):1–49,2003.[31] C. Mascia, A. Porretta, and A. Terracina. Nonhomogeneous Dirichlet problems for degenerate parabolic-hyperbolic equations.
Arch. Ration. Mech. Anal. , 163(2):87–124, 2002.[32] A. Michel, C. Canc`es, T. Gallou¨et, and S. Pegaz. Numerical comparison of invasion percolation models and finitevolume methods for buoyancy driven migration of oil in discontinuous capillary pressure fields. In preparation,2009.[33] A. Michel and J. Vovelle. Entropy formulation for parabolic degenerate equations with general Dirichletboundary conditions and application to the convergence of FV methods.
SIAM J. Numer. Anal. , 41(6):2262–2293 (electronic), 2003.[34] F. Otto. L -contraction and uniqueness for quasilinear elliptic-parabolic equations. J. Differential Equations ,131:20–38, 1996.[35] N. Seguin and J. Vovelle. Analysis and approximation of a scalar conservation law with a flux function withdiscontinuous coefficients.
Math. Models Methods Appl. Sci. , 13(2):221–257, 2003.[36] J. D. Towers. Convergence of a difference scheme for conservation laws with a discontinuous flux.
SIAM J.Numer. Anal. , 38(2):681–698 (electronic), 2000.[37] J. D. Towers. A difference scheme for conservation laws with a discontinuous flux: the nonconvex case.
SIAMJ. Numer. Anal. , 39(4):1197–1218 (electronic), 2001.[38] C. J. van Duijn, J. Molenaar, and M. J. de Neef. The effect of capillary forces on immiscible two-phase flowsin heterogeneous porous media.
Transport in Porous Media , 21:71–93, 1995.[39] J. Vovelle. Convergence of finite volume monotone schemes for scalar conservation laws on bounded domains.