CONVERGENCE OF FINITE VOLUME APPROXIMATIONS FOR A ...

Report 2 Downloads 175 Views
c 2004 Society for Industrial and Applied Mathematics 

SIAM J. NUMER. ANAL. Vol. 42, No. 1, pp. 228–251

CONVERGENCE OF FINITE VOLUME APPROXIMATIONS FOR A NONLINEAR ELLIPTIC-PARABOLIC PROBLEM: A “CONTINUOUS” APPROACH∗ ¨ GUTNIC‡ , AND PETRA WITTBOLD‡ BORIS A. ANDREIANOV† , MICHAEL Abstract. We study the approximation by finite volume methods of the model parabolicelliptic problem b(v)t = div (|Dv|p−2 Dv) on (0, T ) × Ω ⊂ R × Rd with an initial condition and the homogeneous Dirichlet boundary condition. Because of the nonlinearity in the elliptic term, a careful choice of the gradient approximation is needed. We prove the convergence of discrete solutions to the solution of the continuous problem as the discretization step h tends to 0, under the main hypotheses that the approximation of the operator div (|Dv|p−2 Dv) provided by the finite volume scheme is still monotone and coercive, and that the gradient approximation is exact on the affine functions of x ∈ Ω. An example of such a scheme is given for a class of two-dimensional meshes dual to triangular meshes, in particular for structured rectangular and hexagonal meshes. The proof uses the rewriting of the discrete problem under a “continuous” form. This permits us to directly apply the Alt–Luckhaus variational techniques which are known for the continuous case. Key words. doubly nonlinear elliptic-parabolic equations, finite volume methods, convergence of approximate solutions, continuous approach AMS subject classifications. 35J60, 35K55, 35K65, 35M10, 65M12, 76M12 DOI. 10.1137/S00361429014000062

1. Introduction. Let Ω be an open bounded polygonal domain in Rd , d ≥ 1, and T > 0. We consider the initial boundary value problem for a system of nonlinear elliptic-parabolic equations:   b(v)t = div ap (Dv) on Q = (0, T ) × Ω, v=0 on Σ = (0, T ) × ∂Ω, (1.1)  b(v)(0, ·) = u0 on Ω, where 1 < p < ∞ and div ap (Dv) = div (|Dv|p−2 Dv) is the N -dimensional pLaplacian, i.e., d N

ap : ξ = (ξ1 , . . . , ξN ) ∈ (R )

p−2

→ |ξ|

ξ=

 i,j

|ξij |2

p/2−1

(ξ1 , . . . , ξN ) ∈ (Rd )N .

We assume that  b : RN → RN is continuous cyclically monotone; i.e., (1.2) there exists a convex differentiable function Φ : RN → R s.t. b = ∇Φ, normalized by b(0) = 0 and Φ(0) = 0. Moreover, we assume (1.3)

u0 ∈ L1 (Ω)N

with

Ψ(u0 ) ∈ L1 (Ω),

∗ Received by the editors December 21, 2001; accepted for publication (in revised form) May 30, 2003; published electronically January 28, 2004. http://www.siam.org/journals/sinum/42-1/40006.html † LATP, CMI, Universit´ e de Provence, Technopole de Chˆ ateau-Gombert, 39, rue Fr´ed´ eric JoliotCurie, 13453 Marseille Cedex 13, France ([email protected]). ‡ IRMA, Universit´ e Louis Pasteur, 7, rue Ren´e Descartes, 67084 Strasbourg Cedex, France ([email protected], [email protected]).

228

FINITE VOLUME APPROXIMATIONS: A “CONTINUOUS” APPROACH

229

where Ψ is the Legendre transform of Φ given by Ψ : z ∈ RN → sup



σ∈RN

1

0

(z − b(sσ))σ ds = sup (σz − Φ(σ)). σ∈RN

Equations of elliptic-parabolic type (1.1) arise as models of the flow of fluids through porous media (cf., e.g., [6, 12]). They have already been studied extensively in the literature in the last decade from a theoretical point of view (cf., e.g., [1, 21, 22, 12, 7, 26, 8, 10, 2]). Existence of weak solutions of general systems of ellipticparabolic equations has been proved in [1], using Galerkin approximations and timediscretization. Similar results have been obtained later by other authors using different methods (e.g., using a semigroup approach as in [7, 8] in the case N = 1). In particular, it is known that in the case of the system (1.1), for any u0 satisfying (1.3), there exists a weak solution of (1.1), where the weak solution is defined as follows. Denote by E the Banach space Lp (0, T ; W01,p (Ω))N and by E  its dual; E  =   Lp (0, T ; W −1,p (Ω))N , where p = p/(p − 1) is the conjugate exponent of p. Denote by ·, · E  ,E the duality pairing between E  and E. Definition 1.1. A function v ∈ E is a weak solution of the problem (1.1) if b(v) ∈ L∞ (0, T ; L1 (Ω))N and b(v)t ∈ D (Q)N can be extended to a functional χ on E satisfying  (1.4) ap (Dv) · Dφ = 0 for all φ ∈ E, χ, φ E  ,E + Q

(1.5) χ, ξ

E  ,E

  =− b(v) ξt − u0 (·) ξ(0, ·) Q



for all ξ ∈ E with ξt ∈ L∞ (Q)N , ξ(T, ·) = 0.

Note that if v is a weak solution of (1.1), then, by the “chain rule” lemma of [1], one has (1.6)

B(v) ∈ L∞ (0, T ; L1 (Ω))N , B : z ∈ RN

where  1  b(z)z − Φ(z) ≡ → (b(z) − b(sz))zds ≡ Ψ(b(z)) ∈ R. 0

From the results of [26, 10] it also follows that, in the scalar case N = 1, there is uniqueness of a weak solution of (1.1). To our knowledge, the question of uniqueness is open in the case N ≥ 2. In this paper we study the convergence of time-implicit approximations by finite volume numerical schemes for the model nonlinear elliptic-parabolic problem (1.1). Finite volume methods are well suited for numerical simulation of processes where extensive quantities are conserved, and they are popular methods among engineers in hydrology where equations of this type arise. Therefore justification of convergence of this numerical approximation process is of particular interest. In [17] the finite volume method has been studied and convergence of this approximation procedure has been proved for problem (1.1) in the particular case p = 2, N = 1. The same method has also been studied for this equation (i.e., p = 2, N = 1) in the presence of an additional convection term (cf. [18, 14]), and for a nonlinear diffusion problem in [16]. To our knowledge, in the case p = 2, only the convergence of finite element methods has been studied (cf. [19, 11, 5, 20] and their references).

230

B. A. ANDREIANOV, M. GUTNIC, AND P. WITTBOLD

Let us emphasize that our main object is not only to prove the convergence of some finite volume methods for (1.1), but also to develop a “continuous” approach for this proof. The main idea of this adaptation is to rewrite the discrete finite volume scheme under an equivalent continuous form and to apply known stability techniques for the continuous equation (cf. [1] and [2, Chap. V] for the version we use) in order to get convergence of discrete solutions to a solution of the continuous problem. The “continuous” approach and the convergence result have already been presented in [3]. In section 2, we describe the finite volume schemes and in particular the admissible flux approximations we use. We show the existence and uniqueness of the solution of a finite volume scheme and give some a priori estimates on discrete solutions. Then we state the convergence result. In section 3, we show in Proposition 3.3 that the solution of a finite volume scheme, originally satisfying a discrete system of algebraic equations, also verifies a “continuous” formulation similar to (1.4), (1.5). This representation makes clear in which sense finite volume schemes approximate the elliptic operator in (1.1); we prove that this approximation is consistent. In section 4 we prove the convergence theorem, passing to the limit in the “continuous” formulation of Proposition 3.3. In section 5, we analyze the two admissibility conditions imposed in section 2. For d = 2, we propose a scheme on meshes dual to triangular meshes that enters into our framework; in particular, we have the convergence result on structured rectangular and hexagonal meshes. We consider the p-Laplacian as a prototype of a class of the so-called Leray–Lionstype operators; in [4], we discuss the extension of the techniques presented above to a particular case of the p-Laplacian operator with convection, studied in [12]. In order to simplify the notation, we restrict the exposition to the scalar equation (N = 1). The proofs of the auxiliary results used in section 4 can be found in [4]. 2. The numerical method. In order to construct approximate solutions to the problem (1.1), we will use the implicit discretization in time and a finite volume scheme in space. 2.1. Finite volume meshes, discrete gradients and finite volume schemes for the problem (1.1). Let Ω be an open bounded polygonal subset of Rd . A finite volume mesh T of Ω is given by a family of open polygonal convex subsets of Ω with positive measure, called “control volumes,” a family of subsets of Ω contained in hyperplanes of Rd , with positive (d−1)-measure (these are the interfaces between control volumes), and a family of points of Ω, one per control volume (these are the “centers” of the volumes). For a volume K with center xK ∈ K , the interfaces contained in ∂Ω are considered as additional “boundary” volumes, unless xK ∈ ∂Ω. For the sake of simplicity, we shall denote by T the family (K )K ∈T of control volumes; (xK )K ∈T denotes the family of their centers. The set of all volumes K such that xK ∈ ∂Ω is denoted by Text , and the set of all volumes K with xK ∈ Ω is denoted by Tint . The set of interfaces K|L such that K or L or both belong to Tint is denoted by E , and K|L denotes the interface between two neighbors K , L ∈ T . For all K|L, KL denotes the “diamond” over K|L, i.e., the smallest convex set of Rd containing K|L, xK and xL . Whenever we use K , K|L, or n to index objects and make summations, we mean that K ∈ T , K|L ∈ E, and n ∈ {1, . . . , [T /k]+1}, where k is the time step of the scheme. Following [15], we give the following definition. Definition 2.1. We say that T is a finite volume mesh of Ω if the following hold: (2.1 i) The closure of the union of all control volumes is Ω.

FINITE VOLUME APPROXIMATIONS: A “CONTINUOUS” APPROACH

231

(2.1 ii) For any (K , L) ∈ (T )2 with K = L, either the (d−1)-dimensional measure of K ∩ L is 0 or K ∩ L = σ for some σ ∈ E (in which case we denote σ = K|L = L|K ). (2.1 iii) For any K ∈ T , there exists a subset EK of E such that ∂ K = ∪σ∈EK σ. Furthermore, E = ∪K∈T EK . We will denote by NK the set of volumes adjacent to K ; i.e., NK = {L ∈ T , K|L ∈ EK }. (2.1 iv) The family of points (xK )K is such that xK ∈ K for all K ∈ T , and it is assumed that the straight line joining xK and xL is orthogonal to K|L whenever L ∈ NK . Denote by m(K) and d(K) the d-dimensional measure and the diameter of K ∈ T , respectively; and denote by m(K|L) the (d−1)-dimensional measure of K|L ∈ E. A mesh T is characterized, in particular, by the following numbers: dist (xK , σ) , d(K) minK minσ∈EK dist (xK , σ) . M (T ) = max card(EK ), ζ∗ (T ) = K size(T )

A finite volume method for (1.1) requires a family (T h , k h ) h of meshes and corresponding time steps k h > 0 such that both the size of the mesh and the time step go to zero. We will assume in our notation that the family is parametrized with h in some subset of (0, 1) whose closure contains zero, and size(T h ) + k h ≤ h. A grid. couple (T h , k h ) will be called

a space-time In relation to a family (T h , k h ) h , we define the numbers ζ ∗ (T ) = min min

size(T ) = max d(K), K

K

σ∈EK

(2.1) M = sup M (T h ) ∈ N, ζ ∗ = inf ζ ∗ (T h ) ∈ R+ , and ζ∗ = inf ζ∗ (T h ) ∈ R+ . h

h

h

Definition 2.2. We say that the family of meshes (T h )h is weakly proportional if M < ∞ and ζ ∗ > 0. We say that the family of meshes (T h )h is strongly proportional if, in addition, ζ∗ > 0. Weak proportionality is standard (cf. [18]). Strong proportionality is a technical assumption which ensures that (T h )h has the interpolation property (cf. sections 2.5 and 5.2). Given a grid (T h , k h ), to each time-space volume QnK = I n × K , I n = (k h (n − n 1), k h n) one associates an unknown value vK ∈ RN . In order to obtain a finite volume n . scheme for (1.1), one “integrates” the equation in (1.1) over each grid volume QK The time derivative in the left-hand side is approximated by the corresponding finite difference. On the right-hand side, one uses the Green formula and then needs to replace the flux on the lateral boundary of QnK by some function of the unknowns n (vK )K ,n . For problem (1.1), this amounts to finding a substitution for Dv in the expression I n×K|L ap (Dv) · νK ,L (where νK ,L is the unit normal vector to K|L pointing from K into L). We will assume that this substitution is in Lp on each interface n I × K|L, typically constant in time and piecewise constant in space. We therefore consider “discrete gradient” operators Dh of the form

h n n )K ,n → (DK| D : (vK L )K|L,n , (2.2) n p n DK| for all K|L, n. L ∈ L (I × K|L) It seems natural, though not necessary, to require that Dh be a linear operator. A finite volume scheme for (1.1) is defined by a grid (T h , k h ) and a discrete gradient Dh associated with the grid. Finally, a finite volume method for (1.1) is

232

B. A. ANDREIANOV, M. GUTNIC, AND P. WITTBOLD



given by a family (T h , k h , Dh ) h of grids and associated discrete gradient operators Dh . In sections 2.3 and 2.5 we state the admissibility conditions for such methods. Now we are able to write the equations for a scheme (T h , k h , Dh ): n n−1   b(vK ) − b(vK ) n h for all K ∈ Tint , = ap (DK| L (x))dx · νK ,L (2.3) m(K) kh K|L h L∈NK for all n ∈ {1, . . . , [T /k ]+1}.

The homogeneous Dirichlet boundary condition is taken into account by assigning (2.4)

n =0 vK

for all

K

h , ∈ Text

for all n ∈ {1, . . . , [T /k h ]+1}.

0 The initial condition is given by any values vK ∈ b−1 (u0K ), where  1 0 h uK = (2.5) u0 for all K ∈ Tint . m(K) K  0 We denote by uh0 the piecewise constant initial function K uK 1lK , where 1lK is the characteristic function of the set K . Other choices of u0K are possible, provided one has uh0 → u0 a.e. on Ω and Ψ(uh0 ) → Ψ(u0 ) in L1 (Ω) as h → 0, where Ψ is defined in the introduction. These properties hold for u0K given by (2.5), due to the convexity of Ψ. We denote by (S h ) the system (2.3), (2.4), (2.5) corresponding to a given finite volume scheme (T h , k h , Dh ).

2.2. Memento on notation. In this section we collect the most used notation related to the finite volume schemes. T: Text , Tint : E: K, L: K|L : EK : NK : xK : d K ,L : dK ,K|L : νK ,L : : KL d(K), m(K): size(T ): m(K|L): |R|: n I : QnK : ΣnK : Υς (K ):

a finite volume mesh; the set of exterior,interior control volumes; the set of interfaces between control volumes; control volumes of T ; the interface between the two neighbors K and L; the set of all interfaces surrounding K ; the set of all neighbors of K ; the “center” of K ; the distance between xK and xL , dK ,L = |xK − xL |; the distance between xK and K|L; one has dK ,K|L + dL,K|L = dK ,L ; the unit normal vector to K|L pointing from K to L; the smallest convex set of Rd containing K|L, xK , and xL ; the diameter and the d-dimensional measure of K , respectively; the size of the mesh T , size(T ) = maxK d(K); the (d−1)-dimensional measure of K|L; the (d+1)-dimensional measure of a set R ⊂ R+ × Rd ; the time interval, I n = ((n−1)k, nk); n the time-space grid element, QK = In ×K; the lateral boundary of QnK , ΣnK = I n ×∂ K ; the union of all control volumes of (T ) that are separated from K by at most (ς − 1) other control volumes;

233

FINITE VOLUME APPROXIMATIONS: A “CONTINUOUS” APPROACH

1lA : h

the characteristic function of a set A; h

h

(T , k , D ):

a finite volume scheme (mesh, time step, discrete gradient);

(S h ): h: M, ζ ∗ :

the corresponding system of equations (2.3),(2.4),(2.5); the discretization parameter, h ≥ size(T h ) + k h ; the weak proportionality bounds for (T h )h ,

ζ∗ :

M = suph maxK card(NK ), and ζ ∗ = inf h minK the strong proportionality bound for (T h )h ,

n vK : h v : uh0 : n DK| L: h D :

K ,L∈NK ; ζ∗ = inf h size(T ) the unknown of the scheme (S h ) corresponding to the volume QnK ;  n 1lQnK ; a discrete solution for the scheme (T h , k h , Dh ), v h = K ,n vK  h 0 the discrete initial data, u0 = K uK 1lK ; n p n the discrete gradient values on I n × K|L, DK| L ∈ L (I × K|L); h n n the discrete gradient operator, D : (vK )K ,n → (DK| L )K|L,n ;

n : D⊥, K|L h D⊥ :

min

minL∈N

K

dK,K|L

d(K)

;

dK,K|L

v n −v n

n = LdK,LK , featuring in the “discrete Lp (0, T ; W 1,p (Ω)) the value D⊥, K|L norm” of v h ; h n n the corresponding operator, D⊥ : (vK )K ,n → (D⊥, ) . K|L K|L,n

h It is convenient to extend Dh (as wellas D⊥ ) to an operator acting from E into Lp (Q). h Let P be the operator from Ω to K|L which projects x ∈ K on ∂ K along the ray joining xK to x. We define the appropriate lifting operator Lh and averaging operator Mh by    n n h DK| Lh (DK| L )K|L,n (t, x) = L (P (x)) 1lI n×KL (t, x), K|L,n  1 h 1 h n N n M : η ∈ L (Q) → M [η] = (ηK )K ,n ⊂ R , ηK = n η. |QK | QnK

We will abusively write Dh for the operators Dh , Lh ◦ Dh , and Lh ◦ Dh ◦ Mh ; and the h . same for D⊥ The following notations, specific to the “continuous” approach, are introduced in sections 2.5 and 3.1. uh : vh : Gh: A: Ah :

the continuous in t interpolation of b(v h ), affine on each time interval I n ; an interpolated solution in E for v h (cf. Definition 2.8); the interpolated gradient operator produced by (T h , k h , Dh ) (cf. Definition 3.2); the elliptic operator in (1.1), A : η ∈ E → −div ap (Dη) ∈ E  ; the finite volume approximation of A produced by the scheme (T h , k h , Dh ), given by Ah : η ∈ E → −div ap (G h [η]) ∈ E  .

2.3. Admissible flux approximations. For simplicity, we consider only the h h gradients that yield fully implicit schemes; in this case D , D⊥ act independently on n each set vK K , and the dependence on n does not matter for their definition. h Let us introduce the operator D⊥ , which appears naturally in the a priori estimates of section 2.4: (2.6)

h : (vK )K → (D⊥,K|L )K|L , D⊥

D⊥,K|L =

vL − vK ∈ R. d K ,L

234

B. A. ANDREIANOV, M. GUTNIC, AND P. WITTBOLD

For ς ∈ N, denote by Υς (K ) the union of all control volumes of T that are  separated from K by at most (ς − 1) other control volumes; for instance, Υ1 (K ) = L∈NK L. The choice of ς corresponds to the choice of control volumes that are really involved in the construction of Dh on ∂ K . Now we can make precise the assumptions on discrete gradient operators of the form (2.2). Definition 2.3. Let (T h , Dh )h be a family of finite volume meshes and corresponding discrete gradient operators. The gradient approximation provided by Dh is admissible if the following hold. (2.3 i) Dh is linear and injective; vK )K ⊂ (Rd )N (2.3 ii) Dh provides a strictly monotone scheme; i.e., for all (vK )K , ( that do not coincide,    1   K|L (x)) dx · νK ,L > 0, (vL−vK ) − ( ap (DK|L (x)) − ap (D vL− vK ) d K|L K|L      K|L )K|L = Dh ( vK )K ; where (DK|L )K|L = Dh (vK )K , (D (2.3 iii) Dh provides a scheme coercive at zero; i.e., there exists a constant C∗> 0, in dependent of h, such that for all (vK )K ⊂ (Rd )N and (DK|L )K|L = Dh (vK )K , one has  1  vL − vK d K|L

 K|L

 p  h  ap (DK|L (x))dx · νK ,L ≥ C∗ D⊥ [(vK )K ]

Lp (Ω)

;

and there exists ς ∈ N, independent of h, such that the following hold. (2.3 iv) For each h, Dh is consistent with affine functions. More exactly, assume that, for K ∈ T h given, there exists an affine function w on Ω such that vL = 1 m(L) L w whenever L ⊂ Υς (K ). Then DK|L (x) = Dw = const for all x ∈ K|L for all L ∈ NK .  ∈ T h and (2.3 v) There exists a constant C ∗ , independent of h, such that, for all K N all sets of values (vK )K of R ,     p p  h  h   D [(vK )K ] ≤ C ∗ D⊥ [(vK )K ] . )  K Υς (K Conditions (2.3 ii) and (2.3 iv) imply strong restrictions on the gradient approximation. We provide some examples of methods with admissible gradient approximation in section 5.1. n 2.4. Discrete solutions. Recall that we consider as unknowns the values vK n on K ∈ Tint , assigning vK to be zero in K ∈ Text . We will repeatedly use the following “summation by parts” formula (cf., e.g., [15]). Remark 2.4. Let T be a finite volume mesh of Ω in the sense of Definition 2.1. Let (vK )K ∈T ⊂ RN , (FK ,L )(K ,L)∈T 2 ⊂ RN . Assume vK = 0 for all K ∈ Text and FK ,L = −FL,K for all K|L ∈ E. Then    vK FK , L = (vK − vL ) FK ,L . K

L∈NK

K|L

FINITE VOLUME APPROXIMATIONS: A “CONTINUOUS” APPROACH

235

 n n n is the corre1lQK )K ,n verifies (S h ), we say that the function v h = K ,n vK If (vK p sponding discrete solution. We prove the discrete version of the L (0, T ; W01,p (Ω))- a h h priori estimate on v h (which is exactly the estimate on D⊥ [v ] in Lp ), and the discrete version of (1.6).

Proposition 2.5. Let (T h , k h ) h be a family of finite volume grids and let

h D h be a family of corresponding discrete gradient operators satisfying property (2.3 iii) of Definition 2.3. Then, for any solution v h of the discrete problem (S h ), there exists a constant C which depends only on p, d, Ω, T , on C∗ in (2.3 iii), and on Ψ(u0 ) L1 (Ω) such that  v n − v n p   1    h h p K  = m(K|L)dK ,L  L (i)  ≤ C; D⊥ [v ] p d K|L,n d K ,L L (Q)      n (ii) ) ≤ C. = sup m(K) B(vK B(v h ) ∞ 1 L

(0,T ;L (Ω))

n∈{1,...,[T /kh ]+1}

K

i Proof. Take i ∈ {1, . . . , [T /k h ]+1} and multiply each term in (2.3) by vK . By (2.3 iii), using Remark 2.4 and (2.4), one gets    h h p i i−1 i h D⊥ [v ] ≤ 0. m(K)(b(vK ) − b(vK )) vK + C∗ k d Ω

K

i i−1 i i−1 i )) vK ≥ B(vK ) − B(vK ). Summing ) − b(vK By the convexity of Φ, one has (b(vK h over i from 1 to n ∈ {1, . . . , [T /k ]+1} and taking into account the convexity of Ψ, we infer  nkh   h h p n D⊥ [v ] m(K)B(vK ) + C∗ d K



 K

0

0



m(K)Ψ(uK ) =

 K



m(K)Ψ

1 m(K)



K

u0









Ψ(u0 ).

Next, let us prove the discrete version of the Poincar´e inequality and of the compact embedding of W 1,p (Ω) in L1 (Ω). Note that we do not need any proportionality assumptions on the mesh. Lemma 2.6. Let Ω ⊂ Rd be a polygonal domain of diameter d(Ω), and let T be a  finite volume mesh of Ω. Let v h = K vK 1lK such that (vK )K ∈T ⊂ R and vK = 0 for all K ∈ Text . Then there exists a constant C which depends only on p and d such that    h h v h Lp (Ω) ≤ C d(Ω)D⊥ [v ] ; (i) Lp (Ω)     h h (ii) for all ∆ > 0, sup |v h (x+∆x) − v h (x)| dx ≤ ∆ × D⊥ [v ] . L1 (Ω) |∆x|≤∆ Rd Proof. (i) For x ∈ Ω, set ψK|L (x) = 1 in the case that the orthogonal projection of on the hyperplane {x1 = 0} contains (0, x2 , . . . , xd ), and set ψK|L (x) = 0 otherwise. One has   1   |v h (x)|p ≤ ψK|L (x) |vL |p − |vK |p  2 K|L   |vL − vK |  ≤C |vK |p−1 + |vL |p−1 . ψK|L (x) dK ,L d K ,L K|L

K|L

236 Since

B. A. ANDREIANOV, M. GUTNIC, AND P. WITTBOLD



older inequality ψK|L (x) dx ≤ m(K|L) d(Ω), one has by the H¨



|v h (x)|p dx Ω 1   p−1  p  v − v p p  1 1

 L K  m(K|L) dK ,L  m(K|L) dK ,L |vK |p + |vL |p . ≤ Cd(Ω)  d d K ,L d K|L K|L

Denote h = size(T ). Assertion (i) will follow by the Young inequality if we show that 1 (2.7)

K|L

d



m(K|L) dK ,L |vK |p + |vL |p ≤ (1 + 2p )



m(K) |vK |p + 2(2h)p

 v − v p  L K  m(K|L) dK ,L   , d K ,L d

1 K|L

K

since h ≤ d(Ω). Denote by R the left-hand side of (2.7). We have dK ,L= dK ,K|L+dL,K|L ; thus R=



m(K)|vK |p +

K

Note that

K|L

p

|vK | dL,K|L ≤

1 d

m(K|L)(|vK |p dL,K|L + |vL |p dK ,K|L ).

2p |vL |p dL,K|L   −vK p (2h)p  vLdK,L  d K ,L

if |vK | ≤ 2|vL |, otherwise.

Indeed, if |vK | > 2|vL |, one has |vL − vK | > 12 |vK | so that  v − v p  L K  |vK |p dL,K|L ≤ |vK |p dK ,L ≤ 2p |vL − vK |p dK ,L ≤ 2p hp   d K ,L . d K ,L Using the same argument for |vL |p dK ,K|L , we obtain the desired estimate (2.7). (ii) Now for x ∈ Rd , set ψ K|L (x) = 1 in the case where the segment [x, x + ∆x] crosses K|L, and set ψ K|L (x) = 0 otherwise. Note that Rd ψ K|L (x) dx ≤ m(K|L) ∆; hence (ii) follows, since    |v h (x) − v h (x+∆x)| dx ≤ ψ K|L (x)|vL − vK | dx Rd

Rd

≤∆

K|L

 K|L

v − v   L K  m(K|L)dK ,L  . d K ,L

Now we can state the result for existence and uniqueness of a discrete solution. Theorem 2.7. Let T h be a finite volume mesh of Ω, k h > 0, and let Dh be a discrete gradient associated to T h . If Dh satisfies (2.3 iii), there exists a solution n (vK )K ,n to the discrete problem (S h ). If Dh satisfies (2.3 ii), the solution is unique. Proof of Theorem 2.7. Using Remark 2.4 and the coercivity of the scheme, we apply the Brouwer fixed point theorem and get existence. Uniqueness follows from the monotonicity of b(·) and the strict monotonicity of the scheme. See [3] for more detailed proofs.

FINITE VOLUME APPROXIMATIONS: A “CONTINUOUS” APPROACH

237

2.5. Interpolation property and main result. Consider a family of finite

n )K ,n be a solution to volume schemes (T h , k h , Dh ) h such that h tends to 0. Let (vK the scheme (S h ) and v h the corresponding discrete solution. We require the existence of what will be called “interpolated solutions” for v h , denoted by v h , such that v h ∈ E; these should be close to v h (asymptotically as h → 0) and satisfy the a priori estimate in E analogous to the estimate of Proposition 2.5(i) n on v h . Moreover, the values (vK )K ,n should be recoverable from v h . To this end, we n require Mh [v h ] = (vK )K ,n .

Definition 2.8. A family of grids T h , k h h has the interpolation property in E n if, for any family (v h )h of functions such that v h |QnK = vK ≡ const for each K ∈ T h ,  h h h n h n ∈ {1, . . . , [T /k ]+1}, with vK = 0 for K ∈ Text and with D⊥ [v ]Lp (Q) ≤ C for all h, there exists a family (v h )h ⊂ E such that vh − vh

(2.8)

→0

as h → 0,

n Mh [v h ] = (vK )K ,n ,

(2.9) (2.10)

Lp (Q)

vh

E

≤ I(C)

with some function I : R+ → R+ independent of h.

If v h is a solution to a finite volume scheme, we say that v h is an interpolated solution for v h . The interpolation property is the main technical assumption required by the “continuous” approach. In section 5.2 we give two conditions ensuring this property. Now let us state the main result

of this paper. Theorem 2.9. Let (T hm , k hm , Dhm ) m∈N be a sequence of finite volume schemes, where k hm + size(T hm ) ≤ hm → 0 as m → ∞. Assume that the family of meshes is weakly proportional, the gradient approximation is admissible, and the interpolation property holds (cf. Definitions 2.2, 2.3, and 2.8). For m ∈ N, let v hm be a discrete solution of (S hm ). Then there exists a subsequence (hml )l∈N such that v hml 0 v in Lp (Q) as l → ∞, where v is a weak solution of the problem (1.1). Note that it suffices to strengthen slightly assumption (2.3 ii) of Definition 2.3 in order to get the strong convergence of v hml to v in Lp (Q) (cf. [4, Corollary 1]). Moreover, in the case when N = 1 the whole sequence converges to the unique solution of (1.1). In this case error estimates can be proved (cf., e.g., [15] for the linear case), but this is not the purpose of the present paper. In what follows, we write k instead of k h and omit subscripts in sequences (hm ) and(hml ), simply writing that h tends to zero.  n n 3. The “continuous” approach. Take the discrete solution v h = K ,n vK 1lQK h h produced by the finite volume scheme (S ). Let v ∈ E be a corresponding interpolated solution. We will show that there exist functions uh ∈ L1 (Q) and Gh ∈ Lp (Q) such that uh (0, ·) = uh0 (·) and uh t = div ap (Gh ) in the weak sense of Definition 1.1, and the functions uh , Gh can be recovered from the interpolated solution. More exactly, we prove in Proposition 3.3 below that uh t ∈ D can be extended to χh ∈ E  and  h χ , φ E  ,E + (3.1) ap (G h [v h ]) · Dφ = 0 for all φ ∈ E, Q

238

B. A. ANDREIANOV, M. GUTNIC, AND P. WITTBOLD

  uh ξt − uh0 (·) ξ(0, ·) (3.2) χh , ξ E  ,E = − Q



for all ξ ∈ E with ξt ∈ L∞ (Q)N , ξ(T, ·) = 0,

with an operator G h : E → Lp (Q) to be defined. The analogy of (3.1), (3.2) with (1.4), (1.5) in Definition 1.1 plays the key role in the proof of the convergence result of Theorem 2.9. 3.1. Interpolated gradient and the “continuous” form of the scheme. First define uh as the piecewise affine in t interpolation of b(v h ):    t − kn n n−1 h n u = b(vK ) + (3.3) (b(vK ) − b(vK )) 1lQnK . k K ,n Then (3.2) holds, since uh (0, ·) = uh0 (·) and the piecewise constant function uh t extends to χh ∈ E  by  (3.4) uh t φ for all φ ∈ E. χh , φ E  ,E = Q

Next, note that in (S h ) the numerical flux is prescribed on the boundary of h ,n and a each control volume; we will extend it to Q as follows. For given K ∈ Tint function FKn : ∂ K → R, consider the following Neumann problem in the factor space W(K ) = W 1,p (K )/R:     div ap (Dw) = 1 FKn on K , m(K) (3.5) K|L K ∈NK  ap (Dw) · νK |∂ K = FKn , h with m(K) > 0, we where νK is the exterior unit normal vector to ∂ K . For K ∈ Text drop in (3.5) the Neumann boundary condition on ∂ K ∩ ∂Ω and seek w ∈ W 1,p (K ) with w|∂ K ∩∂Ω = 0.   h ). Then (3.5) Lemma 3.1. Let FKn ∈ Lp (∂ K ) (FKn ∈ Lp (∂ K \ ∂Ω), if K ∈ Text admits a unique solution. The proof is standard, using the coercivity and monotonicity argument [25, Chap. 2, Th. 2.1]. Now we can introduce the interpolated gradient operator. Definition 3.2. The interpolated gradient operator G h : E → Lp (Q) maps η ∈ E into G h [η] given by   h n n G [η] = DηK 1lQnK , where ηK ∈ W(K ) solves   K,n           1  n n n  − a (Dη ) · Dϕ + ϕ a (D ) · ν = ϕ ap (DK|  p p K K K|L L ) · νK m (K) K K|L K L∈N K|L L∈NK K  1,p 1,p h   for all ϕ ∈ W (K ) (for all ϕ ∈ W (K ) with ϕ|∂K∩∂Ω = 0, in case K ∈ Text )     n n h   and the values DK| L (x) are given by (DK|L )K|L,n = D [η].

If v h solves (S h ), we set Gh = G h [v h ]. We remark that Gh = G h [v h ] by property (2.9) of interpolated solutions v h . We show that (3.1) follows from (3.4) and the conservation of fluxes. n Proposition 3.3. Assume that (vK )K ,n is a solution of (S h ). Let v h be a correh sponding interpolated solution, let u and χh be defined by (3.3) and (3.4), respectively,

FINITE VOLUME APPROXIMATIONS: A “CONTINUOUS” APPROACH

239

and let G h be the interpolated gradient operator of Definition 3.2. Then (3.1), (3.2) hold. Proof. It remains to check (3.1). By (3.3), for all K ∈ T h and n, we have  n n−1 1  ) − b(vK ) b(vK h n h − u t − div ap (G ) = ap (DK| L (x))dx · νK ,L = 0 k m(K) K|L L∈NK

n everywhere on QK because of (2.3). Therefore, using (3.4) and integrating by parts n in each QK , we have   uh t φ + ap (Gh ) · Dφ χh , φ E  ,E + ap (G h [v h ]) · Dφ Q  Q     n φ ap (DK| (uh t − div ap (Gh ))φ + = L ) · νK , L n n×K|L I Q K K,n L ∈N K,n K   n =0+ φ ap (DK| L ) · (νK ,L + νL,K ) = 0. K|L,n

I n×K|L

3.2. Properties of the interpolated gradient and consistency. In view of (3.1) and (1.4), it is natural to compare the elliptic operator in (1.1), A : η ∈ E → −div ap (Dη) ∈ E  ,

(3.6) with the operators (3.7)

Ah : η ∈ E → −div ap (G h [η]) ∈ E  .

Indeed, Ah can be considered as the finite volume approximation of A, whence the following definition.

Definition 3.4. Let (T h , k h , Dh ) h be a family of finite volume schemes for the problem (1.1), with size(T h ) + k h ≤ h → 0. We say that the approximation of (1.1) by these schemes is consistent if, for all η ∈ E, one has Ah [η] → A[η] in E  as h → 0. In this section we prove

the following result. Theorem 3.5. Let (T h , k h , Dh ) h be a family of finite volume schemes with a weakly proportional family of meshes and an admissible gradient approximation (cf. Definitions 2.2 and 2.3). Then it provides a consistent approximation of (1.1), in the sense of Definition 3.4. The proof of Theorem 3.5 is based upon the following properties of the interpolated gradient operator G h .

Proposition 3.6. Let (T h , k h , Dh ) h be a family of finite volume schemes with admissible gradient approximation and weakly proportional family of meshes. (i) Thereexists a constant C such that for all η ∈ E and H ⊂ Q such that m ni H = i=1 QK , i      h p |Dη|p , G [η] ≤ C H

(3.8)

where Υς+1 (H) = bounded on E and

m

i=1 I

ni

Υς+1 (H)

× Υς+1 (K i ). In particular, G h h are uniformly

 h  G [η] p ≤C η L (Q)

E.

240

B. A. ANDREIANOV, M. GUTNIC, AND P. WITTBOLD

(ii) The operators (G h )h are locally equicontinuous on E. More exactly, there exists a constant C(R) such that, whenever η E ≤ R and µ E ≤ R,   h min{p−1, p−11 }

G [η] − G h [µ] p (3.9) ≤ C(R) η − µ E , L (Q)   min{(p−1)2 , p−11 }

ap (G h [η]) − ap (G h [µ]) p (3.10) ≤ C(R) η − µ E . L (Q) In the statement above and in the rest of this section, C denotes a generic constant that depends only on p, d, Ω, on M, ζ ∗ of (2.1), and on C∗ , C ∗ , ς of Definition 2.3, unless the additional dependence on R is specified. The proof uses the standard properties of the function ap (·) (cf. [19, 12]): for all y1 , y2 ∈ Rd ,   1 < p ≤ 2;  |ap (y1 ) − ap (y2 )|p ≤ C |y1 − y2 |p , p−2   (3.11)  |a (y ) − a (y )|p ≤ C |y − y |p |y |p + |y |p p−1 , p ≥ 2; 1 2 p 2 1 2 p 1

(3.12)

  

2−p  p2  2 |y1 − y2 |p ≤ C (ap (y1 ) − ap (y2 )) · (y1 − y2 ) |y1 |p + |y2 |p , 1 < p ≤ 2; |y1 − y2 |p ≤ C (ap (y1 ) − ap (y2 )) · (y1 − y2 ),

p ≥ 2.

Before turning to the proofs of Proposition 3.6 and Theorem 3.5, note the following three lemmas. Lemma 3.7. Let K ⊂ Rd be a bounded convex domain of Rd of diameter d(K) and d-dimensional measure m(K). Assume that K contains a ball of radius ζ ∗ d(K) > 0. 1 Then there exists a constant C such that, assigning w = m(K) w, one has K   |w − w|p ≤ C (d(K))p−1 |Dw|p ∂K

K

1,p

for all w ∈ W (K ), where w|∂ K is understood in the sense of traces. Proof. Applying, e.g., the proofs of [13, Theorems 59, 60, and 76] with p = 2 replaced by a general p ∈ (1, +∞), we obtain the claim of the lemma with C depending on p, d, and the Lipschitz continuity of ∂ K . Due to the convexity of K , C actually depends only on p, d, and ζ ∗ . h )h Lemma 3.8. Let (T h )h be a weakly proportional family of meshes, and let (D⊥ be the operators defined by (2.6). Then there exists a constant C such that for all K , n for all η ∈ E,      h p [η] ≤ C |Dη|p . D⊥  Qn K

I n×Υ1 (K)

n 1 n = Mh [η] and ηK| Proof. Let ηK L = km(K|L) I n×K|L η in the sense of traces. K ,n By definition,  n     n p  ηL − ηK  h p  1   k m(K|L) dK ,K|L  D⊥ [η] =  d d K ,L Qn K|L K      n p n n ηK| ηK|L − ηLn p 1 L − ηK + ≤C k m(K|L) dK ,K|L d (dK ,K|L )p dK ,K|L (dL,K|L )p−1 K|L    n  p n n 1 1  ηK|L − ηK   ηK|L − ηLn p     . k m(K|L) dK ,K|L  k m(K|L) dL,K|L  ≤C  +C  d d d d K ,K|L L,K|L K|L K|L

241

FINITE VOLUME APPROXIMATIONS: A “CONTINUOUS” APPROACH

By the convexity of the function z → |z|p and Lemma 3.7,   n n p n p p−1 k m(K|L) |ηK| − η | ≤ |η − η | ≤ C d (K) L K K I n×K|L

Qn K

|Dη|p ,

and the same holds if K and L are exchanged. Hence by (2.1) we have           h p p p ≤C |Dη| |Dη|p . |Dη| + D⊥ [η] ≤ C Qn K

L∈NK

Qn K

I n×Υ1 (K)

Qn L



Lemma 3.9. Let (T h , k h , Dh ) h be a family of finite volume schemes with a weakly proportional family of meshes and an admissible gradient approximation. Then the following hold: (i) For all R > 0 there exists a constant C(R) such that, whenever η E ≤ R and µ E ≤ R,   p min{p,p }     . d(K) ap (Dh [η]) − ap (Dh [µ]) ≤ C(R) η − µ E Σn K

K ,n

(ii) There exists a constant C such that for all η ∈ E and H, Υς+1 (H) as in Proposition 3.6 one has   m     h p d(K i ) [η] ≤ C |Dη|p . D n i=1

ΣKi

Υς+1 (H)

 p Proof. First take K and consider ϕK = ap (Dh [η])−ap (Dh [µ]) ∈ L1 (∂ K ). Recall that the values of Dh have been extended from ∂ K inside K by means of the projection operator P h (cf. section 2.2). Hence by (2.1) we have    K d(K) |ϕ | = d(K) |ϕK | ∂K K|L L∈NK  (3.13)  d(K)  d K h = d |ϕ ◦ P | ≤ ∗ |ϕK ◦ P h |. dK ,K|L KL∩K ζ K L∈NK If 1 < p ≤ 2, (3.13)   d(K) K ,n

i

and (3.11) yield  p     ap (Dh [η]) − ap (Dh [µ]) ≤ C

Σn K

K ,n

Qn K

 p  h  D [η] − Dh [µ] .

In turn, (2.3 i), (2.3 iv), Lemma 3.8, and (2.1) imply that  p p       h   h h [η] − D [µ] ≤ C D D⊥ [η − µ]  n n QK I ×Υς (K) K ,n   p (3.14) K ,n      |Dη − Dµ|p = C η − µ ≤C D(η − µ) ≤ C K ,n

I n×Υς+1 (K)

Q

p E

,

which was the claim of (i) for 1 < p ≤ 2. Furthermore, we remark that (3.14) also holds for p ≥ 2, in particular with η = 0 or µ = 0. Therefore for p ≥ 2, using (3.13), and then (3.11) and the H¨ older inequality, we get  K ,n

 d(K)

Σn K

  p     h h ap (D [η])−ap (D [µ]) ≤ C K ,n

Qn K

 p  h  D [η]−Dh [µ]

 pp

  p−2 p−1 × Rp .

242

B. A. ANDREIANOV, M. GUTNIC, AND P. WITTBOLD

Thus, in the case p ≥ 2, (i) also follows from (3.14).  The proof of (ii) is similar, using the identity |ap (y)|p = |y|p instead of inequalities (3.11). Proof of Proposition 3.6. Recalling Definition 3.2, for all K ∈ Tint and n, denote n (respectively, by µnK ) a function in W 1,p (K ) that solves (3.5) with FKn (x) = by ηK h n n n ap (DK|L (x))·νK for x ∈ K|L ∈ EK , where (DK| L )K|L,n = D [η] (respectively, (DK|L )K|L,n = h n n D [µ]). In other words, each of ηK , µK verifies the integral identity corresponding to n (3.5) with all test functions in W 1,p (K ). Taking for the test function (ηK − µnK ), n subtracting the two identities, and integrating in t ∈ I , we obtain      n n ap (DηK − DµnK ) − ap (DµnK ) · DηK Qn K    (3.15)   n n −µn ap (Dh [η]) − ap (Dh [µ]) · νK , ηK −µnK − ηK = K In ∂K

n n 1 n n −µn = where ηK older K m(K) KηK −µK for a.a. t ∈ I . Summing over K , n, using the H¨ inequality and Lemmas 3.7 and 3.9(i), we have from (3.15)      ap (G h [η]) − ap (G h [µ]) · G h [η] − G h [µ] Q   −1 n n − µn | d(K) p | ηK − µnK − ηK ≤ K n I ∂ K K ,n   1   × d(K) p ap (Dh [η]) − ap (Dh [µ])   p1   n n p 1−p n n ≤ d(K) | ηK − µK − ηK − µK | (3.16) In ∂K K ,n 1    p p    × d(K) ap (Dh [η]) − ap (Dh [µ]) K ,n



  

Σn K

 p1 n

min{p/p ,1}

n p

|DηK − DµK | η−µ E Qn K ,n K     min{p/p ,1} = G h [η] − G h [µ] η−µ E . p L (Q)

In the same manner, taking µ = 0 and using Lemma 3.9(ii), we get         1   h p  h p p G [η] C G [η] H

H

 1 Υς+1 (H)

|Dη|p

p

,

which proves (i). Now if 1 < p ≤ 2, (3.12), (3.16), and the H¨ older inequality yield  p  h  G [η] − G h [µ] p L (Q)    p  p  p  2  h p  h    p h ≤ C G [η] − G [µ] η−µ E + G h [µ] p G [η] p Lp (Q)

L (Q)

 2−p 2

L (Q)

Using (3.8), we obtain (3.9). Now (3.10) follows by (3.11). If p ≥ 2, (3.12) and (3.16) readily yield (3.9); hence (3.10) follows by (3.11).

.

FINITE VOLUME APPROXIMATIONS: A “CONTINUOUS” APPROACH

243

  Proof of Theorem 3.5. We have to prove that ap (Dη) − ap (G h [η])Lp (Q) → 0 as h → 0. Let us first prove the theorem for the case of η ∈ E that is piecewise constant in t and piecewise affine in x. Let J ⊂ Q be the set of discontinuities of Dη. Clearly, J is of finite d-dimensional Hausdorff measure Hd (J).  For ς given in Definition 2.3, let us introduce H h = {K,n | I n×Υς (K)∩J=Ø} QnK . Note that |H h | ≤ (ς + 1)h Hd (J) → 0 as h → 0; likewise, |Υς+1 (H h )| → 0 as h → 0. Therefore, by Proposition 3.6(i), we have      p   p p h →0 |Dη| + |Dη| ap (Dη) − ap (G [η]) ≤ C Hh

Hh

Υς+1 (H h )

as h → 0. Moreover, for all QnK such that QnK ∩ H h = Ø, we have G h [η] ≡ Dη on QnK . Indeed, we have D[η] ≡ const on Υς+1 (QnK ). Therefore Dh [η]|QnK ≡ Dη = const by property (2.3 iv) of admissible gradient approximations. Hence Dw = Dη satisfies the boundary condition in (3.5); the equation is also satisfied, since div ap (Dη) ≡ 0 1 on K and m(K) a (Dh [η]) · νK = ap (Dη) · ∂K νK = 0. ∂K p  It follows that ap (Dη) − ap (G h [η]) p → 0 as h → 0, which was our claim. L (Q)

Now let us approximate an arbitrary function η in E by functions ηl that are piecewise constant in t and piecewise affine in x. Note that we can always choose this sequence ηl in E such that ηl → η in E and a.e. on Q as l → ∞, and |Dηl |p are dominated by an L1 (Q) function independent of l. We have         ≤ ap (Dη) − ap (Dηl )  ap (Dη) − ap (G h [η]) p L (Q) Lp (Q)     (3.17)     +ap (Dηl ) − ap (G h [ηl ])  + ap (G h [ηl ]) − ap (G h [η])  . Lp (Q)

Lp (Q)

As l → 0, the first term in the right-hand side of (3.17) converges to zero by the Lebesgue dominated convergence theorem, independently of h. The second one converges to zero as h → 0 for all l fixed. Finally, by Proposition 3.6(ii), the third one converges to zero as l → ∞ uniformly in h. Hence the result follows. 4. Proof of Theorem 2.9. In the context of continuous dependence upon the data of weak solutions to “general” elliptic-parabolic problems (cf. [2, Chap.V]), the proof of convergence of weak solutions of approximating problems is based upon the three essential arguments (A), (B), and (C) below. (A) A priori estimates, using (1.2) and the Alt–Luckhaus chain rule lemma (cf. [1, 26, 10]). (B) Strong compactness in the parabolic term, using a variant of the Kruzhkov lemma (cf. [23]): Lemma 4.1 (cf. [4], [2, Chap. V]). Let Ω be an open domain in Rd , Q = h 1 h (0, T )×Ω, and let the families of functions (u )h , (Fα )h,α be bounded in L (Q) ∂ h and satisfy ∂t u = |α|≤m Dα Fαh in D (Q). Assume that uh can be extended by zero outside Q, and one has  (4.1) sup |uh (t, x + ∆x) − uh (t, x)| dxdt ≤ ω(∆), with lim ω(∆) = 0, |∆x|≤∆

Rd+1

∆→0

where ω(·) does not depend on h. Then (uh )h is relatively compact in L1 (Q).

244

B. A. ANDREIANOV, M. GUTNIC, AND P. WITTBOLD

(C) Convergence in the elliptic term, using a variant of the Minty–Browder argument (cf., e.g., [25]). Lemma 4.2 (cf. [4], [2, Chap. V]). Let E be a Banach space, E  its dual and ·, · E  ,E denote the duality product of elements of E  and E. Let (v h )h ⊂ E and v h 0 v as h → 0. Let Ah be a sequence of monotone operators from E ∗ to E  such that Ah [v h ] 0 −χ for some χ ∈ E  . Assume that Ah converge pointwise to some operator A, and A is hemicontinuous (i.e., continuous in the weak-∗ topology of E  along each direction). Assume that (4.2)

lim inf Ah [v h ], v h E  ,E ≤ −χ, v E  ,E . h→0

Then χ + A[v] = 0, and (4.2) necessarily holds with equality. Taking advantage of the “continuous” form (3.1), (3.2) of the discrete problem (S h ), we can prove the convergence of finite volume approximate solutions in the same way, using the discrete a priori estimates shown in Propositions 2.5 and 3.6(i), using next Lemma 4.1, and then using finally Lemma 4.2 together with the essential consistency result of Theorem 3.5. Proof of Theorem 2.9. Let v h be the solution of (S h ). Let v h be a corresponding interpolated solution, and let Ah be the finite volume approximate of the operator A in (1.1) (cf. (3.6), (3.7)). Note that all the convergences we state below take place up to extraction of a subsequence.   h h  [v ] Lp (Q) ≤ const uniformly in h so that the (A) By Proposition 2.5(i), D⊥

family (v h )h is bounded in E, by (2.10). Hence there exists a function v ∈ E such that v h 0 v in E as h → 0. By (2.8), one also has v h 0 v in Lp (Q). (B) We claim that the family (uh )h given by (3.3) is relatively compact in L1 (Q). Indeed, let us check the assumptions of Lemma 4.1. We have uh t = div ap (G h [v h ]) in

 D (Q) by (3.1), (3.2), and the family ap (G h [v h ]) h is bounded in Lp (Q) by Propo h h sition 2.5(i), equation (2.10), and Proposition 3.6(i) (note that A [v ] h is thus bounded in E  ). Furthermore, (3.3) yields   uh L1 (Q) ≤ 2 |b(v h )| + k h m(K)|u0K |, Q

K

and one has |b(z)| ≤ δB(z) + sup|ζ|≤1/δ |b(ζ)| for all δ > 0 (cf., e.g., [1]). By Propo sition 2.5(ii) and since uh0 = K m(K)|u0K | → u0 in L1 (Ω) as h → 0, it follows that (uh )h is bounded in L1 (Q). Finally, by Proposition 2.5(i) and Lemma 2.6(ii), we obtain (4.1) with uh replaced by v h . Hence the estimate (4.1) for uh follows by (3.3), as in the continuous case (cf. [1]); see [4] for the detailed proof. Thus the claim of (B) follows, and there exists a function u ∈ L1 (Q) such that h u → u in L1 (Q) and a.e. on Q. In addition, we claim that u = b(v), where v is the weak limit of v h in E. It suffices to show that v h 0 v in L1 (Q) and b(v h ) → u in L1 (Q), and then apply the monotonicity argument of [9]; see [4] for the detailed proof. (C) By (A), we have v h 0 v in E. We claim that v is a weak solution of (1.1). h h  ] = is By Proposition 3.3, χh + A

[v 0 in E and the  initial condition (3.2) h h verified for all h. The family A [v ] h is bounded in E (cf. (B)), thus (χh )h is weak-∗ relatively compact in E  . By (3.4), (B), and Definition 1.1, we also have ∗ χh = uh t → b(v)t = χ in D (Q). Hence Ah [v h ] = −χh 0 −χ in E  .

FINITE VOLUME APPROXIMATIONS: A “CONTINUOUS” APPROACH

245

Moreover, passing to the limit in (3.2), using (B) and the convergence of uh0 to u0 in L1 (Ω), we get (1.5). Consequently, by the chain rule argument [1, Lemma 1.5] we have   (4.3) Ψ(u0 ). Ψ(b(v(T, ·))) + −χ, v E  ,E = − Ω



On the other hand, by (3.4), (3.3), (2.9), and the monotonicity of b(·), we have −χh , v h E  ,E = − = −



 1 n n−1 (b(vK ) − b(vK )) vh n k QK K,n

n−1 n n m(K)(b(vK ) − b(vK )) vK

K ,n

≤ −

 [T /kh ]+1 )) + m(K)Ψ(u0K ) m(K)Ψ(b(vK K  h = − Ψ(b(v (T, ·))) + Ψ(uh0 ).

 K





Recall that Ψ(uh0 ) → Ψ(u0 ) in L1 (Ω). Without loss of generality, we can assume that v h (T, ·) → v(T, ·) a.e. on Ω; hence by the Fatou lemma and (4.3) we get (4.2). Next, the operators Ah are monotone. Indeed, take ϕ ∈ E and (ϕnK )K ,n = Mh [ϕ]. Arguing as in the proof of Proposition 3.3, integrating by parts in QnK , and cancelling the boundary terms, we get h



(4.4)

 

ap (G [η]) · Dϕ = − ϕ div ap (G h [v h ]) n QK K ,n     1 n = − ap (DK| ϕ× L (x))dx · νK ,L n m (K) K|L QK K ,n L∈NK    n = −k ϕnK ap (DK| L (x))dx · νK ,L .

A [η], ϕ

E  ,E

=

h

Q

K ,n

L∈NK

K|L

Substituting (4.4) and applying Remark 2.4, we infer by property (2.3 ii) of Definition 2.3 that η ], η − η E  ,E Ah [η] − Ah [    1  n n n n n (ηL−ηK ) − ( ap (DK| ηLn− ηK ) = k L (x)) − ap (DK|L (x)) dx · νK ,L ≥ 0, d K|L,n K|L n h n )K ,n = Mh [η], (DK| . where (ηK L )K|L,n = D [η], and the same for η h Finally, by Theorem 3.5, A converge pointwise to the hemicontinuous operator A · = −div ap (D ·). By Lemma 4.2 we conclude that χ + A[v] = 0 in E  . Thus (1.4) holds and v is a weak solution of (1.1).

5. Examples of admissible methods. By an admissible method, we mean a method which provides an admissible gradient approximation and weakly proportional meshes satisfying the interpolation property. Recall that in this case, we have the convergence of the finite volume approximation (cf. Theorem 2.9). In this section, we prove that such admissible methods exist.

246

B. A. ANDREIANOV, M. GUTNIC, AND P. WITTBOLD

5.1. On discrete gradients. In this section we construct an admissible gradient for a family (T h )h of finite volume meshes of the Vorono¨ı kind dual to a family (T h )h of triangular meshes.  to denote a triangle of the mesh T h ; for Let us introduce some notation. We use O h h  = ∆xK xL xM (the triangle with the  ∈ T , there exist K , L, M ∈ T such that O all O corners xK , xL , xM ). The three interfaces K|L, L|M , M|K intersect at point xO, which is . We require it to be inside O . Let us the center of the circumscribed circle of triangle O denote by SO, SK ,L , SL,M , and SM ,K the surfaces of ∆xK xL xM , ∆xOxK xL , ∆xOxL xM , and ∆xOxM xK , respectively. One has SO = SK ,L +SL,M +SM ,K . −−−→ −→ −−−→ x− Recall that νK ,L = − K xL /dK ,L , νL,M = xL xM /dL,M , νM ,K = xM xK /dM ,K . Note the following elementary lemma.  = ∆xK xL xM be a triangle in R2 , let xO be the center of its Lemma 5.1. Let O circumscribed circle, and let xO ∈ ∆xK xL xM . With the above notation, for all r in R2 , we have r=

 2  SK ,L (r · νK ,L )νK ,L + SL,M (r · νL,M )νL,M + SM ,K (r · νM ,K )νM ,K . SO

This property can be generalized to any polygon in R2 which admits the circumscribed circle.  = ∆xK xL xM , let v h,0 O : R2 → RN be the  ∈ T h such that O Furthermore, for O affine function that takes the values vK , vL , vM at the points xK , xL , xM , respectively. The discrete gradient operator Dh,0 = Lh ◦ Dh,0 is defined by  Dh,0 : (vK )K → Dv h,0 O(x)1lO(x). O ∈Th In the case of structured hexagonal meshes, as well as that of structured rectangular ones, the family (Dh,0 )h is admissible (this will be proved in Proposition 5.2, as a particular case). In general, this construction does not work. Indeed, if the points xK are not the barycenters of K ∈ T h , property (2.3 iv) fails. This can be overcome, for instance, in the following way. For all K ∈ T h , let yK  ∈ T h such that O  = ∆xK xL xM , be the barycenter of K and set σK = xK − yK . For O h 2 N let v : R → R be the affine function that takes the values vK , vL , vM at the points  O yK , yL , yM , respectively. The discrete gradient operator Dh = Lh ◦ Dh is defined by (5.1)

Dh : (vK )K →

  

O ∈T h

DvOh(x)1lO(x);

i.e., the affine interpolation over the triangle ∆yK yL yM is actually used in the triangle ∆xK xL xM . We will take advantage of considering Dh as a perturbation of Dh,0 . For all  ∈ T h , let us define the correction operators O

    σL − σK σM − σL 2 2 SK ,L r · νK ,L + SL,M r · νL,M RO : r ∈ R → SO d K ,L dL,M    σK − σM νM ,K , + SM ,K r · d M ,K

FINITE VOLUME APPROXIMATIONS: A “CONTINUOUS” APPROACH

247

with the notation introduced above. We need to guarantee that the Euclidiean norm  ∈ T h. of RO is less than min{p − 1, 1/(p − 1)} for all O Proposition 5.2. Assume that (T h )h is a family of meshes dual to a family of  ∈ T h are triangles with angles less than or equal to π/2. meshes (T h )h such that all O  ∈ T h , Assume that for all h, for all O 

2 |σM − σL | |σK − σM | |σL − σK | < min{p − 1, 1/(p − 1)}, SK ,L + SL,M + SM ,K d K ,L dL,M d M ,K SO where σK is the difference between the “center” xK of the volume K and its barycenter, etc. Then the family of discrete gradient operators (Dh )h on (T h )h defined by (5.1) is admissible in the sense of Definition 2.3. 1 Proof. Since, for any affine function w on K , one has m(K) w(x)dx = w(yK ), K h where yK is the barycenter of K , property (2.3 iv) holds for D (with ς = 1, by construction). Next, (2.3 i) is clear. Let us establish the relation between Dh,0 and Dh . Denote by D0 , DO the values O   of Dh,0 [(vK )K ] and Dh [(vK )K ], respectively. Let us show that for all O  ∈ T h , on O (5.2)

DO0 = (I − RO)DO.

 = ∆xK xL xM , one has Indeed, if O v h (xL ) − v h (xK ) (vL + DO · σL ) − (vK + DO · σK ) O  = DO · νK ,L = O d K ,L d K ,L vL − vK σL − σK σL − σK = + DO · = DO0 · νK ,L + DO · . d K ,L d K ,L d K ,L Writing the same relation for L, M and RODO, whence (5.2) follows.

M , K,

from Lemma 5.1, we get DO − D0 = O 

h h  ∈ T h such that By Lemma 5.1 and the definition of D⊥ = Lh ◦ D⊥ for all O  = ∆xK xL xM we have O    p p  h,0    D [(vK )K ] = SODh,0 [(vK )K ]  O

  v − v p  v − v p  v − v p    L K  M L K M  ∗ (5.3) SK ,L  ≤ C  + SL,M   + SM ,K   d dL,M d M ,K   K ,L p  h  = C∗ D⊥ [(vK )K ]  O

 ∈ T h and O  ∈ T h we with a constant C ∗ that depends only on p. Since for given K  ∩O  = Ø if and only if O  ∈ Υ1 ( have K K ), it follows that property (2.3 v) holds for the discrete gradient Dh,0 , with ς = 1. Now set θO = RO . We have θO < 1. One has 1 |D0 |; therefore (2.3 v) also holds for Dh . |DO| ≤ (I − RO)−1 |D0 | ≤ 1−θ O  O  O Next, each term in the sum in (2.3 iii) splits into two terms corresponding to 1 , O 2 ∈ T h . Let us the two parts of the interface K|L included in different triangles O h  ∈T , O  = ∆xK xL xM , write down all the terms corresponding to the same triangle O

248

B. A. ANDREIANOV, M. GUTNIC, AND P. WITTBOLD

combine them using Lemma 5.1, and estimate using (5.2): SK ,L (ap (DO) · νK ,L )(DO0 · νK ,L ) + SL,M (ap (DO) · νL,M )(DO0 · νL,M )

SO ap (DO) · DO0 2 1 − θO 1 − θO S S |D0 |p . = O ap (DO) · (I − RO)DO ≥ SO|DO|p ≥ 2(1 + θO)p O O 2 2

+ SM ,K (ap (DO) · νM ,K )(DO0 · νM ,K ) =

−vK | h = D⊥, Property (2.3 iii) for Dh follows, because one has |D0 | ≥ |D0 ·νK ,L | = |vLdK,L K|L O  so that  p  p     h  ≥ D⊥ [(vK )] . SO|DO0|p = Dh,0 [(vK )K ] Lp (Ω) Lp (Ω) O ∈Th

 of Dh,0 [( 0 , D vK )K ] and The proof of (2.3 ii) is similar. Denoting the values D O  O , one can rewrite the sum in (2.3 ii) as D [( vK )K ], respectively, on O      )) · νK ,L (D0 − D  0 ) · νK , L SK ,L (ap (DO) − ap (D O  O  O  ∈Th O     )) · νL,M (D0 − D  0 ) · νL,M + SL,M (ap (DO) − ap (D O  O  O      0 ) · νM , K  )) · νM ,K (D0 − D + SM ,K (ap (DO) − ap (D O  O  O  1  0  )) · (D − D  0 ). SO(ap (DO) − ap (D = O  O  O  2 ∈Th O h

Using (5.2) and denoting by H the Hessian matrix of the function x ∈ R2 → p1 |x|p , we get  )) · (D0 − D 0 ) (ap (DO)−ap (D O  O  O 



 )t = (DO − D O 

0

1

!    ). H(DO +τ (DO − DO)) dτ (I −RO) (DO − D O 

2

For all x ∈ R , x = 0, H(x) is a symmetric matrix with positive eigenvalues λ1 , λ2 such that λ1 /λ2 = p − 1. Thus the condition RO < min{p − 1, 1/(p − 1)} ensures that, for all τ ∈ [0, 1],      + τ (D − D  + τ (D − D  )) (I − R ) r ≥ a rt H(D  )) r > 0 rt H(D O O O O O O O        for all r ∈ R2 , r = 0, with some constant a > 0. Now (2.3 ii) follows. 5.2. On interpolated solutions. First note that it is sufficient to prove the interpolation property in W01,p (Ω) if we require, in addition to the time-independent analogs of (2.8), (2.9) (referred to as (2.8 ), (2.9 )), that    h h (2.10  ) v h W 1,p (Ω) ≤ c × D⊥ [v ] p 0

L (Ω)

with a constant c independent of h. We obtain the interpolation property in E with the function I : C → c × C by taking v h constant on each I n and summing in n ∈ {1, . . . , [T /k h ]+1}.

FINITE VOLUME APPROXIMATIONS: A “CONTINUOUS” APPROACH

249

Lemma 5.3. Let (T h )h be a strongly proportional family of finite volume meshes of Ω ⊂ Rd . Then it has the interpolation property in W01,p (Ω). In order to prove the lemma, we first show that the strong proportionality allows us to majorate the Lp norm of the translates of the discrete solutions v h in Lemma 2.6(ii) by const∆ (h + ∆)p−1 . Then we convolute v h with the appropriate mollifier; finally, we restore the average over each mesh volume as in Lemma 5.4 below. The complete proof is given in [4]. Note that the interpolation property can fail on weakly proportional meshes, at least for p > 2. Indeed, consider Ω = (0, 1)2 . For s ≥ 2, let T s be the finite volume mesh of Ω 1 1 s ), , 4s = {K s , Ls }, where K s = {(x, y) ∈ Ω | x + y < 1/s} with xK s = ( 4s such that Tint 1 1 s s s and L is the interior of the complementary of K with xLs = ( 2 , 2 ). Take v such that  h s p v s ≡ s1/p on K s and v s ≡ 0 on Ls . Then Ω D⊥ [v ] ≤ const uniformly in s. If there 1,p s s exist v ∈ W0 (Ω) interpolated solutions for v , we have v s W 1,p (Ω) ≤ const. Hence 0 by the standard embedding theorem, v s are uniformly bounded. This contradicts the 1 s 1/p → +∞ as s → +∞. fact that m(K s) Ks v = s Nevertheless, we have the following result in the situation close to that of Proposition 5.2. Lemma 5.4. Assume that (T h )h is a weakly proportional family of meshes of  ∈ T h are triangles with angles Ω ⊂ R2 dual to a family of meshes (T h )h such that all O h less than or equal to π/2. Then (T )h has the interpolation property in W01,p (Ω).  h h  Take discrete solutions v = K vK 1lK on each of T such that, for all ∗h  hProof. D [v h ] p ≤ C. Denote by c the generic constant that depends only on p and ζ . ⊥ L (Q) Let v h,0 be the continuous piecewise affine function on Ω that interpolates the values  for all O  ∈ T h (we use the construction and vK , vL , vM at the points xK , xL , xM over O notation of section 5.1). We have v h,0 ∈ W01,p (Ω) and Dv h,0 ≡ Dh,0 [v h ] so that (5.3) yields Dv h,0 Lp (Ω) ≤ c × C. Note that for all x ∈ K ∈ T h , |v h,0 (x) − v h (x)| = |v h,0 (x) − v h,0 (xK )| ≤ d(K) |Dv h,0 (x)|.

(5.4)

Hence v h,0 − v h Lp (Ω) ≤ h Dv h,0 Lp (Ω) → 0 as h → 0. Now take a continuously 2 supp π = differentiable function π : R2 → R+ such that   {x ∈ R | |x| ≤ 1} and x−xK h (for boundary volumes K π = 1. For all K ∈ Tint set ϕK = (ζ ∗m(K) d(K))2 π ζ ∗ d(K) R2 of nonzero measure; i.e., if xK ∈ ∂Ω, an easy modification is needed in order to keep the trace on ∂Ω equal to zero). Set   1 v h = v h,0 + αK ϕK , with αK = vK − v h,0 . m (K) K K Since supp ϕK ⊂ K , by the choice of αK , the family (v h )h verifies (2.9 ). Moreover, m(K) h , for all h, by the H¨ older inequality we get since d(K) 2 ≤ c for all K ∈ T 



h

|v − v

h,0 p

| =



p

|αK |



K

|ϕK |p

 p     m(K) p  h h,0  ≤ c v  m(K)  v − p  m(K)  K d(K)2 K K    1 h h,0 p p /p m (K) |v − v | m (K) = c |v h − v h,0 |p → 0 ≤ c m (K)p K Ω K 

K

1

250

B. A. ANDREIANOV, M. GUTNIC, AND P. WITTBOLD

as h → 0. Thus (v h )h satisfies (2.8 ). In the same manner, using (5.4) we have    h,0 p p h |Dv − Dv | = |αK | |DϕK |p K Ω K   m(K) p   1 p /p h h,0 p |v − v | × m (K) ≤ c m (K) m(K)p d(K)3 K K K     1  1 h h,0 p h,0 p p |v − v | ≤ c d (K) |Dv | = c |Dv h,0 |p . ≤c p p d (K) d (K) K K Ω K K Hence v h

W01,p (Ω)

≤ c v h,0

W01,p (Ω)

≤ c × C, so (v h )h satisfies (2.10 ). Thus (v h )h

can be chosen as interpolated solutions for (v h )h . REFERENCES [1] H. W. Alt and S. Luckhaus, Quasilinear elliptic-parabolic differential equations, Math. Z., 183 (1983), pp. 311–341. [2] B. Andreianov, Quelques probl` emes de la th´ eorie de syst` emes paraboliques d´ eg´ en´ er´ es nonlin´ eaires et de lois de conservation, Ph.D. thesis, Laboratoire de Math´ematiques, Universit´ e de Franche-Comt´e, Besan¸con, France, 2000. [3] B. Andreianov, M. Gutnic, and P. Wittbold, L’approche “continue” pour une m´ ethode de volumes finis, C. R. Acad. Sci. Paris S´er. I Math., 332 (2001), pp. 477–482. [4] B. Andreianov, M. Gutnic, and P. Wittbold, Convergence of Finite Volume Approximations for a Nonlinear Elliptic-Parabolic Problem: A “Continuous” Approach, prepublication 03036, Institut de Recherche Math´ematique Avanc´ee, Universit´ e Louis Pasteur, Strasbourg, France, 2003; also available online from http://www-irma. u-strasb.fr/irma/publications/2003/03036.ps.gz. [5] J. W. Barrett and W. B. Liu, Finite element approximation of the parabolic p-Laplacian, SIAM J. Numer. Anal., 31 (1994), pp. 413–428. [6] J. Bear, Dynamics of Fluids in Porous Media, Elsevier, New York, 1972. ´nilan and P. Wittbold, On mild and weak solutions of elliptic-parabolic problems, [7] Ph. Be Adv. Differential Equations, 1 (1996), pp. 1053–1073. [8] F. Bouhsiss, Etude d’un probl` eme parabolique par les semi-groupes non lin´ eaires, in Analyse Non Lin´ eaire, Publ. Math. UFR Sci. Tech. Besan¸con 15, Univ. Franche-Comt´e, Besan¸con, France, 1995/1997, pp. 133–141. ´zis and W. A. Strauss, Semi-linear second-order elliptic equations in L1 , J. Math. [9] H. Bre Soc. Japan, 25 (1973), pp. 565–590. [10] J. Carrillo and P. Wittbold, Uniqueness of renormalized solutions of degenerate ellipticparabolic problems, J. Differential Equations, 156 (1999), pp. 93–121. [11] S.-S. Chow, Finite element error estimates for nonlinear elliptic equations of monotone type, Numer. Math., 54 (1989), pp. 373–393. ´lin, On a nonlinear parabolic problem arising in some models related [12] J. I. Diaz and F. de The to turbulent flows, SIAM J. Math. Anal., 25 (1994), pp. 1085–1111. [13] Yu. V. Egorov and V. A. Kondratiev, On Spectral Theory of Elliptic Operators, Oper. Theory Adv. Appl. 89, Birkh¨ auser, Basel, 1996. [14] R. Eymard, T. Gallouet, M. Gutnic, R. Herbin, and D. Hilhorst, Approximation by the finite volume method of an elliptic-parabolic equation arising in environmental studies, Math. Models Methods Appl. Sci., 11 (2001), pp. 1505–1528. ¨t, and R. Herbin, Finite volumes methods, in Handbook of Nu[15] R. Eymard, T. Galloue merical Analysis, P. G. Ciarlet and J.-L. Lions, eds., North–Holland, Amsterdam, 2000, pp. 715–1022. ¨t, D. Hilhorst, and Y. Na¨ıt Slimane, Finite volumes and non[16] R. Eymard, T. Galloue linear diffusion equations, RAIRO Math. Model. Numer. Anal., 32 (1998), pp. 747–761. [17] R. Eymard, M. Gutnic, and D. Hilhorst, The finite volume method for an elliptic-parabolic equation, Acta Math. Univ. Comenian. (N.S.), 67 (1998), pp. 181–195. [18] R. Eymard, M. Gutnic, and D. Hilhorst, The finite volume method for Richards equation, Comput. Geosci., 3 (1999), pp. 259–294 (2000). [19] R. Glowinski and A. Marrocco, Sur l’approximation par ´ el´ ements finis d’ordre un, et la r´ esolution, par p´ enalisation-dualit´ e, d’une classe de probl` emes de Dirichlet non lin´ eaires, RAIRO Anal. Num´ er., 9 (1975), pp. 41–76.

FINITE VOLUME APPROXIMATIONS: A “CONTINUOUS” APPROACH

251

[20] N. Ju, Numerical analysis of the parabolic p-Laplacian: Approximation of trajectories, SIAM J. Numer. Anal., 37 (2000), pp. 1861–1884. ˇur, On a solution of degenerate elliptic-parabolic problems in Orlicz-Sobolev spaces. I, [21] J. Kac Math. Z., 203 (1990), pp. 153–171. ˇur, On a solution of degenerate elliptic-parabolic problems in Orlicz-Sobolev spaces. [22] J. Kac II, Math. Z., 203 (1990), pp. 569–579. [23] S. N. Kruzhkov, Results on the nature of the continuity of solutions of parabolic equations and some of their applications, Mat. Zametki, 6 (1969), pp. 97–108 (in Russian); English translation in Math. Notes, 6 (1969), pp. 517–523. [24] J. Leray and J.-L. Lions, Quelques r´ esultats de Viˇsik sur les probl` emes elliptiques non lin´ eaires par les m´ ethodes de Minty-Browder, Bull. Soc. Math. France, 93 (1965), pp. 97–107. [25] J.-L. Lions, Quelques m´ ethodes de r´ esolution des probl` emes aux limites non lin´ eaires, Dunod, Gauthier-Villars, Paris, 1969. [26] F. Otto, L1 -Contraction and uniqueness for quasilinear elliptic-parabolic problems, J. Differential Equations, 131 (1996), pp. 20–38.