A residual-based a posteriori error estimator for a fully-mixed

formulation of the Stokes-Darcy coupled problem

Gabriel N. Gatica

Ricardo Oyarz ua

Francisco-Javier Sayas

Abstract
In this paper we develop an a posteriori error analysis of a new fully mixed nite element
method for the coupling of uid ow with porous media ow in 2D. Flows are governed by
the Stokes and Darcy equations, respectively, and the corresponding transmission conditions
are given by mass conservation, balance of normal forces, and the Beavers-Joseph-Saman
law. We consider dual-mixed formulations in both media, which yields the pseudostress and
the velocity in the uid, together with the velocity and the pressure in the porous medium,
and the traces of the porous media pressure and the uid velocity on the interface, as the
resulting unknowns. The set of feasible nite element subspaces includes Raviart-Thomas
elements of lowest order and piecewise constants for the velocities and pressures, respectively,
in both domains, together with continuous piecewise linear elements for the traces. We derive
a reliable and ecient residual-based a posteriori error estimator for the coupled problem.
The proof of reliability makes use of the global inf-sup condition, Helmholtz decompositions
in both media, and local approximation properties of the Clement interpolant and Raviart-
Thomas operator. On the other hand, inverse inequalities, the localization technique based
on triangle-bubble and edge-bubble functions, and known results from previous works, are
the main tools for proving the eciency of the estimator. Finally, some numerical results
conrming the theoretical properties of this estimator, and illustrating the capability of the
corresponding adaptive algorithm to localize the singularities of the solution, are reported.
Key words: a posteriori error analysis, eciency, reliability, Stokes, Darcy, fully-mixed
Mathematics Subject Classications (2000): 65N15, 65N30, 74F10, 74S05
1 Introduction
The derivation of new nite element methods for the Stokes-Darcy coupled problem, in which
the respective interface conditions are given by mass conservation, balance of normal forces,
and the Beavers-Joseph-Saman law, has become a very active research area lately (see, e.g.
[5], [10], [13], [14], [20], [22], [23], [24], [29], [30], [35], [38], [40], [41], [42], [43], [47] and the
references therein). The above list includes porous media with cracks, nonlinear problems, and
the incorporation of the Brinkman equation in the model (see [10], [23], and [47]). In addition,

CI
2
MA and Departamento de Ingeniera Matematica, Universidad de Concepcion, Casilla 160-C, Concepcion,
Chile, email: [email protected]

Departamento de Ingeniera Matematica, Universidad de Concepcion, Casilla 160-C, Concepcion, Chile, email:
[email protected]

Departamento de Matematica Aplicada, Centro Politecnico Superior, Universidad de Zaragoza, Mara de

Luna, 3 - 50018 Zaragoza, Spain, e-mail: [email protected] Present address: School of Mathematics, University
of Minnesota, 206 Church St. SE, Minneapolis, MN 55455, USA.
1
most of the formulations employed are based on appropriate combinations of stable elements
for the free uid ow and for the porous medium ow, and the rst theoretical results in this
direction go back to [22] and [35]. Indeed, an iterative subdomain method employing the primal
variational formulation and standard nite element subspaces in both domains is proposed in
[22], whereas the primal method in the uid and the dual-mixed method in the porous medium
are applied in [35]. In this way, the approach from [35] yields the velocity and the pressure in
both domains, together with the trace of the porous medium pressure on the interface, as the
main unknowns of the coupled problem. This trace unknown is motivated by the fact that one of
the transmission conditions becomes essential. Then, new mixed nite element discretizations
of the variational formulation from [35] have been introduced and analyzed in [29] and [30].
The stability of a specic Galerkin method is the main result in [29], and the resulting mixed
nite element method is the rst one that is conforming for the primal/dual-mixed formulation
proposed in [35]. The results from [29] are improved in [30] where it is shown that the use of
any pair of stable Stokes and Darcy elements implies the stability of the corresponding Stokes-
Darcy Galerkin scheme. The analysis in [30] hinges on the fact that the operator dening the
continuous variational formulation is given by a compact perturbation of an invertible mapping.
Further techniques utilized in the literature include mortar nite element methods, discontinuous
Galerkin (DG) schemes, and stabilized formulations (see, e.g. [5], [13], [14], [20], [21], [24], [38],
[40], [41], [42], [43]). In particular, the main motivation for employing stabilized formulations
either in both domains or in one of them, is the possibility of approximating the Stokes and
Darcy ows with the same nite element subpaces. Certainly, dierent nite element subspaces
in each ow region may lead to dierent approximation properties for each subproblem. On the
contrary, using the same spaces guarantees the same accurateness along the entire domain and
leads to simpler and more ecient computational codes.
Now, in the recent paper [31] we have developed a new variational approach for the 2D
Stokes-Darcy coupled problem, which allows, on one hand, the introduction of further unknowns
of physical interest, and on the other hand, the utilization of the same family of nite element
subspaces in both media, without requiring any stabilization term. More precisely, in [31]
we consider dual-mixed formulations in both domains, which yields the pseudostress and the
velocity in the uid, together with the velocity and the pressure in the porous medium, as the
main unknowns. The pressure and the gradient of the velocity in the uid can then be computed
as a very simple postprocess of the above unknowns, in which no numerical dierentiation is
applied, and hence no further sources of error arise. In addition, since the transmission conditions
become essential, we impose them weakly and introduce the traces of the porous media pressure
and the uid velocity, which are also variables of importance from a physical point of view, as
the corresponding Lagrange multipliers. Then, we apply the well known Fredholm and Babuska-
Brezzi theories to prove the unique solvability of the resulting continuous formulation and derive
sucient conditions on the nite element subspaces ensuring that the associated Galerkin scheme
becomes well posed. Among the several dierent ways in which the equations and unknowns
can be ordered, we choose the one yielding a doubly mixed structure for which the inf-sup
conditions of the o-diagonal bilinear forms follow straightforwardly. In this way, the arguments
of the continuous analysis can be easily adapted to the discrete case. In particular, a feasible
choice of subspaces is given by Raviart-Thomas elements of lowest order and piecewise constants
for the velocities and pressures, respectively, in both domains, together with continuous piecewise
linear elements for the Lagrange multipliers.
2
On the other hand, it is well known that in order to guarantee a good convergence behaviour
of most nite element solutions, specially under the eventual presence of singularities, one usually
needs to apply an adaptive algorithm based on a posteriori error estimates. These are represented
by global quantities that are expressed in terms of local indicators
T
dened on each element
T of a given triangulation T . The estimator is said to be ecient (resp. reliable) if there
exists C
eff
> 0 (resp. C
rel
> 0), independent of the meshsizes, such that
C
eff
+ h.o.t. |error| C
rel
+ h.o.t. ,
where h.o.t. is a generic expression denoting one or several terms of higher order. In particular,
the a posteriori error analysis of variational formulations with saddle-point structure has already
been widely investigated by many authors (see, e.g. [2], [3], [4], [11], [15], [17], [27], [33], [36], [37],
[39], [44], and the references therein). These contributions refer mainly to reliable and ecient
a posteriori error estimators based on local and global residuals, local problems, postprocessing,
and functional-type error estimates. In addition, the applications include Stokes and Oseen
equations, Poisson problem, linear elasticity, and general elliptic partial dierential equations
of second order. However, up to our knowledge, the rst a posteriori error analysis for the
Stokes-Darcy coupled problem has been provided recently in [8], where a reliable and ecient
residual-based a posteriori error estimator for the variational formulation analyzed in [29] is
derived. Partially following known approaches, the proof of reliability makes use of suitable
auxiliary problems, diverse continuous inf-sup conditions satised by the bilinear forms involved,
and local approximation properties of the Clement interpolant and Raviart-Thomas operator.
Similarly, Helmholtz decomposition, inverse inequalities, and the localization technique based
on triangle-bubble and edge-bubble functions, are the main tools for proving the eciency of
the estimator.
Motivated by the discussion in the above paragraphs, our purpose now is to additionally
contribute in the direction of [8] and provide the a posteriori error analysis of the fully-mixed
variational approach introduced in [31]. According to this, the rest of this work is organized as
follows. In Section 2 we recall from [31] the Stokes-Darcy coupled problem and its continuous
and discrete fully-mixed variational formulations. The kernel of the present work is given by
Section 3, where we develop the a posteriori error analysis. In Section 3.1 we employ the
global continuous inf-sup condition, Helmholtz decompositions in both domains, and the local
approximation properties of the Clement and Raviart-Thomas operators, to derive a reliable
residual-based a posteriori error estimator. An interesting feature of our proof of reliability
is the previous transformation of the global continuous inf-sup condition into an equivalent
estimate involving global inf-sup conditions for each one of the components of the product
space to which the vector of unknowns belongs. Then, in Section 3.2 we apply again Helmholtz
decompositions, inverse inequalities, and the localization technique based on triangle-bubble and
edge-bubble functions to prove the eciency of the estimator. This proof benets partially from
the fact that some components of the a posteriori error estimator coincide with those obtained in
[8] and the related work [15]. Finally, numerical results conrming the reliability and eciency
of the a posteriori error estimator and showing the good performance of the associated adaptive
algorithm, are presented in Section 4.
We end this section with some notations to be used below. In particular, in what follows we
utilize the standard terminology for Sobolev spaces. In addition, if O is a domain, is a closed
3
Lipschitz curve, and r R, we dene
H
r
(O) := [H
r
(O)]
2
, H
r
(O) := [H
r
(O)]
22
, and H
r
() := [H
r
()]
2
.
However, for r = 0 we usually write L
2
(O), L
2
(O), and L
2
() instead of H
0
(O), H
0
(O), and
H
0
(), respectively. The corresponding norms are denoted by | |
r,O
(for H
r
(O), H
r
(O), and
H
r
(O)) and | |
r,
(for H
r
() and H
r
()). Also, the Hilbert space
H(div ; O) :=

w L
2
(O) : div w L
2
(O)

,
is standard in the realm of mixed problems (see, e.g. [12] or [32]). The space of matrix valued
functions whose rows belong to H(div ; O) will be denoted H(div; O). The Hilbert norms of
H(div ; O) and H(div; O) are denoted by | |
div ;O
and | |
div;O
, respectively. On the other
hand, the symbol for the L
2
() and L
2
() inner products
', `

:=

, L
2
(), ', `

:=

, L
2
()
will also be employed for their respective extensions as the duality products H
1/2
() H
1/2
()
and H
1/2
() H
1/2
(). Finally, we employ 0 as a generic null vector, and use C and c, with
or without subscripts, bars, tildes or hats, to mean generic positive constants independent of
the discretization parameters, which may take dierent values at dierent places.
2 The Stokes-Darcy coupled problem
In this section we follow very closely the presentation from [31] to introduce the model problem
and the corresponding continuous and discrete mixed variational formulations.
2.1 The model problem
The Stokes-Darcy coupled problem consists of an incompressible viscous uid occupying a region

S
, which ows back and forth across the common interface into a porous medium living in
another region
D
and saturated with the same uid. Physically, we consider a simplied 2D
model where
D
is surrounded by a bounded region
S
(see Figure 2.1 below). Their common
interface is supposed to be a Lipschitz curve and we assume that
D
= . The remaining
part of the boundary of
S
is also assumed to be a Lipschitz curve
S
. For practical purposes,
we can assume that both
S
and are polygons. The unit normal vector eld on the boundaries
n is chosen pointing outwards from
S
(and therefore inwards to
D
when seen on ). On
we also consider a unit tangent vector eld t in any xed orientation of this closed curve.
The governing equations in
S
are those of the Stokes problem, which are written in the
following non-standard velocity-pressure-pseudostress formulation:

S
= p
S
I + u
S
in
S
, div
S
+ f
S
= 0 in
S
,
div u
S
= 0 in
S
, u
S
= 0 on
S
,
(2.1)
where > 0 is the viscosity of the uid, u
S
is the uid velocity, p
S
is the pressure,
S
is the
pseudostress tensor, I is the 2 2 identity matrix, and f
S
L
2
(
S
) are known source terms.
4

t
n
n

D

S
S

Figure 2.1: Geometry of the problem

Here, div is the usual divergence operator acting on vector elds, and div denotes the action
of div along the rows of each tensor. On the other hand, the ow equations in
D
are those of
the linearized Darcy model:
u
D
= Kp
D
in
D
, div u
D
= f
D
in
D
, (2.2)
where the unknowns are the pressure p
D
and the ow u
D
, and the source term, given by f
D

L
2
(
D
), satises

D
f
D
= 0. The matrix valued function K, describing permeability of
D
divided by the viscosity , is symmetric, has L

(
D
) components and is uniformly elliptic.
Finally, the transmission conditions on are given by
u
S
n = u
D
n on ,

S
n +
1
(u
S
t) t = p
D
n on ,
(2.3)
where :=

( Kt) t

is the friction coecient, and is a positive parameter to be determined

experimentally. The rst equation in (2.3) corresponds to mass conservation on , whereas the
normal and tangential components of the second one constitute the balance of normal forces
and the Beavers-Joseph-Saman law, respectively. Throughout the rest of the paper we assume,
without loss of generality, that is a positive constant.
We complete the description of our model problem by observing that the equations in the
Stokes domain (cf. (2.1)) can be rewritten equivalently as

d
S
= u
S
in
S
, div
S
+ f
S
= 0 in
S
,
p
S
=
1
2
tr
S
in
S
, u
S
= 0 on
S
,
(2.4)
where tr stands for the usual trace of tensors, that is tr :=
11
+
22
, and

d
:=
1
2
(tr ) I
is the deviatoric part of the tensor := (
ij
)
22
.
5
We end this section by remarking that, though the geometry described by Figure 2.1 was
choosen to simplify the presentation, the case of a uid owing only across a part of the boundary
of the porous medium does not yield further complications for the a posteriori error analysis
of the problem. We already discussed this issue in [31, Section 2.1], in connection with the
respective a priori error analysis, and further details can be found in [24].
2.2 The fully-mixed variational formulation
We rst dene the global unknows := (
S
, u
D
, , ) and u := (u
S
, p
D
), where and are
the traces := u
S
[

and := p
D
[

. Then we recall from [31, Lemma 3.5] that the coupled

problem given by (2.2), (2.3), and (2.4) has the one-dimensional kernel dened by
((
S
, u
D
, , ), (u
S
, p
D
)) :
S
= c I, u
D
= 0, = 0, = c, u
S
= 0, p
D
= c ; c R .
Hence, in order to solve this indetermination, we introduce
L
2
0
(
D
) :=

q L
2
(
D
) :

D
q = 0

,
and dene the product spaces
X := H(div;
S
) H(div ;
D
) H
1/2
() H
1/2
() , M := L
2
(
S
) L
2
0
(
D
) ,
endowed with the product norms
||
X
:= |
S
|
div,
S
+ |v
D
|
div ;
D
+ ||
1/2,
+ ||
1/2,
:= (
S
, v
D
, , ) X,
and
|v|
M
:= |v
S
|
0.
S
+ |q
D
|
0,
D
v := (v
S
, q
D
) M.
In this way, as explained in [31, Sections 2 and 3]), it suces to consider from now on the
following modied variational formulation of (2.2), (2.3), and (2.4): Find (, u) XM such
that
/(, ) + B(, u) = T() := (
S
, v
D
, , ) X,
B(, v) = ((v) v := (v
S
, q
D
) M,
(2.5)
where
T() := 0, ((v) = (((v
S
, q
D
)) := (f
S
, v
S
)
S
(f
D
, q
D
)
D
, (2.6)
and / and B are the bounded bilinear forms dened by
/(, ) := a((
S
, u
D
), (
S
, v
D
)) + b((
S
, v
D
), (, ))
+ b((
S
, u
D
), (, )) c((, ), (, )) ,
(2.7)
with
a((
S
, u
D
), (
S
, v
D
)) :=
1
(
d
S
,
d
S
)
S
+ (K
1
u
D
, v
D
)
D
,
b((
S
, v
D
), (, )) := '
S
n, `

'v
D
n, `

,
c((, ), (, )) :=
1
' t, t`

+ ' n, `

' n, `

,
6
and
B(, v) := (div
S
, v
S
)
S
(div v
D
, q
D
)
D
. (2.8)
Hereafter we utilize, for each S, D, the following notations
(u, v)

:=

uv, (u, v)

:=

u v, (, )

:=

: ,
for all u, v L
2
(

), u, v L
2
(

), and , L
2
(

), where : := tr(
t
).
We nd it important to remark that and constitute the Lagrange multipliers associated
with the transmission conditions (2.3). In addition, we notice that (2.5) is equivalent to the
variational formulation dened in [31, Section 3.2, eq. (3.2)], in which
S
is decomposed into

S
= + , with H
0
(div;
S
) and R, where
H
0
(div;
S
) :=

H(div;
S
) :

S
tr() = 0

.
The following result taken from [31] establishes, in particular, the well-posedness of (2.5).
Theorem 2.1 For each pair (T, () X

there exists a unique (, u) X M solution

to (2.5), and there exists a constant C > 0, independent of the solution, such that
|(, u)|
XM
C

|T|
X
+ |(|
M

. (2.9)
Proof. See [31, Theorem 3.9].
We end this section with the converse of the derivation of (2.5). More precisely, the following
theorem establishes that the unique solution of (2.5), with T and ( given by (2.6), solves the
original transmission problem described in Section 2.1. This result will be used later on in
Section 3.2 to prove the eciency of our a posteriori error estimator. We remark that no extra
regularity assumptions on the data, but only f
S
L
2
(
S
) and f
D
L
2
(
D
), are required here.
Theorem 2.2 Let (, u) H Q be the unique solution of the variational formulation (2.5)
with T and ( given by (2.6). Then div
S
= f
S
in
S
,
1

d
S
= u
S
in
S
, u
S
H
1
(
S
),
div u
D
= f
D
in
D
, u
D
= Kp
D
in
D
, p
D
H
1
(
D
), u
D
n + n = 0 on ,

S
n + n

( t) t = 0 on , = p
D
on , = u
S
on , and u
S
= 0 on
S
.
Proof. It basically follows by applying integration by parts backwardly in (2.5) and using suitable
test functions. We omit further details.
2.3 The Galerkin formulation
Although the analysis in [31] provides general hypotheses for the well-posedness of a Galerkin
scheme of (2.5), it suces to consider in what follows the particular case described in [31, Section
5]. Let T
S
h
and T
D
h
be respective triangulations of the domains
S
and
D
formed by shape-
regular triangles T of diameter h
T
, and assume that T
S
h
and T
D
h
match in , so that their union
is a triangulation of
S

D
. Then, for each T T
S
h
T
D
h
we let RT
0
(T) be the local
Raviart-Thomas space of order 0, that is
RT
0
(T) : = span

1
0

0
1

x
1
x
2

,
7
where x :=

x
1
x
2

is a generic vector of R
2
, and for each S, D we dene the global spaces
H
h
(

) :=

v
h
H(div ;

) : v
h
[
T
RT
0
(T) T T

, (2.10)
and
L
h
(

) :=

q
h
:

R : q
h
[
T
P
0
(T) T T

.
Hereafter, given a non-negative integer k and a subset S of R
2
, P
k
(S) stands for the space of
polynomials dened on S of degree k. Next, we let
h
be the partition of inherited from
T
S
h
(or T
D
h
), and assume, without loss of generality, that the number of edges of
h
is even. The
case of an odd number of edges is easily reduced to the even case (see [31]). Then, we let
2h
be the partition of arising by joining pairs of adjacent edges of
h
. Note that because
h
is
inherited from one of the interior triangulations, it is automatically of bounded variation (that
is, the ratio of lengths of adjacent edges is bounded) and, therefore, so is
2h
.
Employing the above notations, we now introduce
H
h
(
S
) := :
S
R
22
: c
t
H
h
(
S
) c R
2
,
L
h
(
S
) := L
h
(
S
) L
h
(
S
) ,
L
h,0
(
D
) := L
h
(
D
) L
2
0
(
D
) ,

h
() :=
h
C() :
h
[
e
P
1
(e) e edge of
2h
,

h
() :=
h
()
h
() ,
and the product spaces
X
h
:= H
h
(
S
) H
h
(
D
)
h
()
h
() and M
h
:= L
h
(
S
) L
h,0
(
D
) .
In this way, the Galerkin scheme of (2.5) becomes: Find (
h
, u
h
) X
h
M
h
such that
/(
h
, ) +B(, u
h
) = T() := (
S
, v
D
, , ) X
h
,
B(
h
, v) = ((v) v := (v
S
, q
D
) M
h
,
(2.11)
where
h
= (
S,h
, u
D,h
,
h
,
h
) and u
h
:= (u
S,h
, p
D,h
).
The following theorems, also taken from [31], provide the well-posedness of (2.11), the asso-
ciated Cea estimate, and the corresponding theoretical rate of convergence.
Theorem 2.3 Assume that T
S
h
and T
D
h
are quasiuniform in a neighborhood of . Then the
Galerkin scheme (2.11) has a unique solution (
h
, u
h
) X
h
M
h
. Moreover, there exist
C
1
, C
2
> 0, independent of h, such that
|(
h
, u
h
)|
XM
C
1

|T[
X
h
|
X

h
+ |([
M
h
|
M

,
and
|
h
|
X
+ |u u
h
|
M
C
2

inf

h
X
h
|
h
|
X
+ inf
v
h
M
h
|u v
h
|
M

.
8
Proof. See [31, Theorems 5.3 and 5.4].
Theorem 2.4 Assume the same hypotheses of Theorem 2.3, and let (, u) X M and
(
h
, u
h
) X
h
M
h
be the unique solutions of the continuous and discrete formulations (2.5)
and (2.11), respectively. Assume that there exists (0, 1] such that
S
H

(
S
), div
S

H

(
S
), u
D
H

(
D
), and div u
D
H

(
D
). Then, u
S
H
1+
(
S
), p
D
H
1+
(
D
),
H
1/2+
(), H
1/2+
(), and there exists C > 0, independent of h and the continuous
and discrete solutions, such that
|(, u) (
h
, u
h
)|
XM
C h

|
S
|
,
S
+ |div
S
|
,
S
+ |u
D
|
,
D
+ |div u
D
|
,
D
+ |u
S
|
1+,
S
+ |p
D
|
1+,
D

.
(2.12)
Proof. See [31, Theorem 5.5].
3 A residual-based a posteriori error estimator
We rst introduce some notations. For each T T
S
h
T
D
h
we let c(T) be the set of edges of
T, and we denote by c
h
the set of all edges of T
S
h
T
D
h
, that is
c
h
= c
h
(
S
) c
h
(
S
) c
h
(
D
) c
h
() ,
where c
h
(
S
) := e c
h
: e
S
, c
h
(

) := e c
h
: e

for each S, D,
and c
h
() := e c
h
: e . Note that c
h
() is the set of edges dening the partition

h
. Analogously, we let c
2h
() be the set of double edges dening the partition
2h
. In what
follows, h
e
stands for the diameter of a given edge e c
h
c
2h
(). Now, let D, S and
let q [L
2
(

)]
m
, with m 1, 2, such that q[
T
[C(T)]
m
for each T T

h
. Then, given
e c
h
(

), we denote by [q] the jump of q across e, that is [q] := (q[

T
)[
e
(q[
T
)[
e
, where T

and
T

are the triangles of T

h
having e as an edge. Also, we x a unit normal vector n
e
:= (n
1
, n
2
)
t
to the edge e, which points either inward T

or inward T

, and let t
e
:= (n
2
, n
1
)
t
be the
corresponding xed unit tangential vector along e. Hence, given v L
2
(

) and L
2
(

)
such that v[
T
[C(T)]
2
and [
T
[C(T)]
22
, respectively, for each T T

h
, we let [v t
e
] and
[ t
e
] be the tangential jumps of v and , across e, that is [v t
e
] := (v[
T
)[
e
(v[
T
)[
e
t
e
and [ t
e
] := ([
T
)[
e
([
T
)[
e
t
e
, respectively. From now on, when no confusion arises, we
simply write t and n instead of t
e
and n
e
, respectively. Finally, for suently smooth scalar,
vector and tensors elds q, v := (v
1
, v
2
)
t
and := (
ij
)
22
, respectively, we let
curl v :=

v
1
x
2

v
1
x
1
v
2
x
2

v
2
x
1

, curl q :=

q
x
2
,
q
x
1

t
,
rot v :=
v
2
x
1

v
1
x
2
, and rot :=

12
x
1

11
x
2
,

22
x
1

21
x
2

t
.
Next, let (, u) XM and (
h
, u
h
) := ((
S,h
, u
D,h
,
h
,
h
), (u
S,h
, p
D,h
)) X
h
M
h
be the
unique solutions of (2.5) and (2.11), respectively. Then, we introduce the global a posteriori
9
error estimator:
:=

TT
S
h

2
S,T
+

TT
D
h

2
D,T

1/2
, (3.1)
where, for each T T
S
h
:

2
S,T
:= |f
S
+div
S,h
|
2
0,T
+ h
2
T
|rot
d
S,h
|
2
0,T
+ h
2
T
|
d
S,h
|
2
0,T
+

eE(T)E
h
(
S
)
h
e
|[
d
S,h
t]|
2
0,e
+

eE(T)E
h
(
S
)
h
e
|
d
S,h
t|
2
0,e
+

eE(T)E
h
()
h
e
|u
S,h
+
h
|
2
0,e
+

eE(T)E
h
()

h
e

S,h
n +
h
n

(
h
t) t

2
0,e
+ h
e

d
S,h
t +
h
t

2
0,e

,
and for each T T
D
h
:

2
D,T
:= |f
D
div u
D,h
|
2
0,T
+ h
2
T
|rot (K
1
u
D,h
)|
2
0,T
+ h
2
T
|K
1
u
D,h
|
2
0,T
+

eE(T)E
h
(
D
)
h
e

[K
1
u
D,h
t]

2
0,e
+

eE(T)E
h
()
h
e

K
1
u
D,h
t +
d
h
dt

2
0,e
+

eE(T)E
h
()

h
e
|u
D,h
n +
h
n|
2
0,e
+ h
e
|p
D,h

h
|
2
0,e

.
3.1 Reliability of the a posteriori error estimator
The main result of this section is stated as follows.
Theorem 3.1 There exists C
rel
> 0, independent of h, such that
|
h
|
X
+ |u u
h
|
M
C
rel
. (3.2)
We begin the derivation of (3.2) by recalling that the continous dependence result given by
(2.9) is equivalent to the global inf-sup condition for the continuous formulation (2.5). Then,
applying this estimate to the error (
h
, u u
h
) X M, we obtain
|(
h
, u u
h
)|
XM
C sup
( ,v)XM
(,v)=0
[R(, v)[
|(, v)|
XM
, (3.3)
where R : X MR is the residual operator dened by
R(, v) := /(
h
, ) + B(, u u
h
) + B(
h
, v), (, v) X M.
More precisely, according to (2.5) and the denitions of / and B (cf. (2.7), (2.8)), we nd that
for any (, v) := ((
S
, v
D
, , ), (v
S
, q
D
)) X M there holds
R(, v) = R
1
(
S
) + R
2
(v
D
) + R
3
() + R
4
() + R
5
(v
S
) + R
6
(q
D
) ,
10
where
R
1
(
S
) :=
1

d
S,h
:
d
S

S
u
S,h
div
S
'
S
n,
h
`

,
R
2
(v
D
) :=

D
K
1
u
D,h
v
D
+

D
p
D,h
div v
D
+ 'v
D
n,
h
`

,
R
3
() := '
S,h
n, `

' n,
h
`

' t,
h
t`

,
R
4
() := 'u
D,h
n, `

+ '
h
n, `

,
R
5
(v
S
) :=

S
v
S
(f
S
+div
S,h
) ,
and
R
6
(q
D
) :=

D
q
D
(f
D
div u
D,h
) .
Hence, the supremum in (3.3) can be bounded in terms of R
i
, i 1, ..., 6, which yields
|(
h
, u u
h
)|
XM
C

sup

S
H(div;
S
)

S
=0
[R
1
(
S
)[
|
S
|
div;
S
+ sup
v
D
H(div;
D
)
v
D
=0
[R
2
(v
D
)[
|v
D
|
div ;
D
+ sup
H
1/2
()
=0
[R
3
()[
||
1/2,
+ sup
H
1/2
()
=0
[R
4
()[
||
1/2,
+ sup
v
S
L
2
(
S
)
v
S
=0
[R
5
(v
S
)[
|v
S
|
0,
S
+ sup
q
D
L
2
0
(
D
)
q
D
=0
[R
6
(q
D
)[
|q
D
|
0,
D

.
(3.4)
Throughout the rest of this section we provide suitable upper bounds for each one of the terms
on the right hand side of (3.4). The following lemma, whose proof follows from straightforward
applications of the Cauchy-Schwarz inequality, is stated rst.
Lemma 3.1 There hold
sup
v
S
L
2
(
S
)
v
S
=0
[R
5
(v
S
)[
|v
S
|
0,
S
|f
S
+div
S,h
|
0,
S
=

TT
S
h
|f
S
+div
S,h
|
2
0,T

1/2
, (3.5)
and
sup
q
D
L
2
0
(
D
)
q
D
=0
[R
6
(q
D
)[
|q
D
|
0,
D
|f
D
div u
D,h
|
0,
D
=

TT
D
h
|f
D
div u
D,h
|
2
0,T

1/2
. (3.6)
The next lemma estimates the suprema on the spaces dened in the interface .
11
Lemma 3.2 There exist C
3
, C
4
> 0, independent of h, such that
sup
H
1/2
()
=0
[R
3
()[
||
1/2,
C
3

eE
h
()
h
e

S,h
n +
h
n

(
h
t) t

2
0,e

1/2
, (3.7)
and
sup
H
1/2
()
=0
[R
4
()[
||
1/2,
C
4

eE
h
()
h
e
|u
D,h
n +
h
n|
2
0,e

1/2
. (3.8)
Proof. It is clear from the denition of R
3
that
R
3
() = '
S,h
n +
h
n

(
h
t) t, `

H
1/2
() ,
and hence
sup
H
1/2
()
=0
[R
3
()[
||
1/2,
=

S,h
n +
h
n

(
h
t) t

1/2,
. (3.9)
In order to estimate

S,h
n +
h
n

(
h
t) t

1/2,
in terms of local quantities we now apply
a technical result from [16]. In fact, taking
S
= 0, v
D
= 0 and = 0 in the rst equation of
(2.11), we have
'
S,h
n +
h
n

(
h
t) t, `

= 0
h
() ,
which says that
S,h
n +
h
n

(
h
t) t is L
2
()-orthogonal to
h
(). Hence, applying [16,
Theorem 2], and recalling that
h
and
2h
are of bounded variation, we deduce that

S,h
n +
h
n

(
h
t) t

2
1/2,
C

eE
2h
()
h
e

S,h
n +
h
n

(
h
t) t

2
0,e
C

eE
h
()
h
e

S,h
n +
h
n

(
h
t) t

2
0,e
,
which, together with (3.9), yields (3.7).
The proof of (3.8) proceeds analogously. In fact, it is easy to see that
sup
H
1/2
()
=0
[R
4
()[
||
1/2,
= |u
D,h
n +
h
n|
1/2,
,
and hence, noting also from the rst equation of (2.11) that u
D,h
n+
h
n is L
2
()-orthogonal
to
h
(), another straightforward application of [16, Theorem 2] yields the required estimate.
We omit further details here.
Our next goal is to bound the rst two suprema on the right hand side of (3.4), for which
we need several preliminary results. We begin with the following lemma showing the existence
of stable Helmholtz decompositions for H(div;
D
) and H(div;
S
).
12
Lemma 3.3
a) For each v
D
H(div;
D
) there exist w H
1
(
D
) and H
1
(
D
), with

S
= 0,
such that there hold v
D
= w + curl in
D
, and
|w|
1,
D
+ ||
1,
D
C
D
|v
D
|
div ;
D
,
where C
D
is a positive constant independent of v
D
.
b) For each
S
H(div;
S
) there exist H
1
(
S
) and H
1
(
S
) such that there hold

S
= + curl in
S
, and
||
1,
S
+ ||
1,
S
C
S
|
S
|
div;
S
,
where C
S
is a positive constant independent of
S
.
Proof. Given v
D
H(div;
D
), we let G be a smooth convex domain containing
D
, and let
z H
1
0
(G) H
2
(G) be the unique solution of
z =

div v
D
in
D
0 in G `

in G, z = 0 on G.
It follows that
|z|
2,G
C |div v
D
|
0,
D
C |v
D
|
div ;
D
,
and hence, dening w := z in
D
, we nd that
div w = div v
D
in
D
and |w|
1,
D
|z|
2,
D
|z|
2,G
C |v
D
|
div ;
D
.
In addition, since div (v
D
w) = 0 and
D
is connected, there exists H
1
(
D
), with

D
= 0, such that v
D
w = curl in
D
. In this way, using the generalized Poincare
inequality and the above estimate for w, we deduce that
||
1,
D
C [[
1,
D
= C |curl |
0,
D
= C |v
D
w|
0,
D
C|v
D
|
div ;
D
,
which completes the proof of a).
We now let
S
H(div;
S
). Since
S
is not necessarily connected, we rst perform a
suitable extension of
S
to the domain :=
S

D
, and then apply a) to each row of
the resulting tensor. More precisely, let
S,i
H(div ;
S
) be the i-th row of
S
, i 1, 2,
and let
i
H
1
(
D
) be the unique solution of the Neumann problem:

i
=
'
S,i
n, 1`

[
D
[
in
D
,

i
n
=
S,i
n on ,

i
= 0 .
Then we dene
ext
i
=

S,i
in
S

i
in
D
, and notice that
ext
i
H(div ; ) and
|
ext
i
|
div ;
|
S,i
|
div ;
S
+ |
i
|
div ;
D
|
S,i
|
div ;
S
+ C |
S,i
n|
1/2,
C |
S,i
|
div ;
S
.
13
Proceeding as in the proof of a), but now for
ext
i
H(div ; ), we deduce the existence of
w
i
H
1
() and
i
H
1
(), with

i
= 0, such that
ext
i
= w
i
+ curl
i
in , and
|w
i
|
1,
+ |
i
|
1,
C |
ext
i
|
div ;
C |
S,i
|
div ;
S
.
Hence, the proof of b) follows by dening i-th row of := w
i
[

S
and := (
1
[

S
,
2
[

S
).
The Raviart-Thomas interpolation operator

h
: H
1
(

) H
h
(

) (cf. (2.10)), S, D,
which, given v H
1
(

), is characterized by

h
(v) H
h
(

) and

h
(v) n =

e
v n edge e of T

h
, (3.10)
will also be needed in what follows. Note that as a consequence of (3.10), there holds
div (

h
(v)) = {

h
(div v) , (3.11)
where {

h
, S, D, is the L
2
(

)-orthogonal projector onto the piecewise constant functions

on

. A tensor version of

h
, say

h
: H
1
(

) H
h
(

), which is dened row-wise by

h
,
and a vector version of {

h
, say P

h
, which is the L
2
(

)-orthogonal projector onto the piecewise

constant vectors on

, might also be required. The local approximation properties of

h
(and
hence of

h
) are stated as follows.
Lemma 3.4 For each S, D there exist constants c
1
, c
2
> 0, independent of h, such that
for all v H
1
(

) there hold
|v

h
(v)|
0,T
c
1
h
T
|v|
1,T
T T

h
,
and
|v n

h
(v) n|
0,e
c
2
h
1/2
e
|v|
1,Te
edge e of T

h
,
where T
e
is a triangle of T

h
containing e on its boundary.
Proof. See [12].
We will also utilize the Clement interpolation operators I

h
: H
1
(

) X
,h
(cf. [19]),
where
X
,h
:= v C(

) : v[
T
P
1
(T) T T

h
for each S, D .
In addition, we will make use of a vector version of I

h
, say I

h
: H
1
(

) X
,h
:= X
,h
X
,h
,
which is dened componentwise by I

h
. The following lemma establishes the local approximation
properties of I

h
(and hence of I

h
).
Lemma 3.5 For each S, D there exist constants c
3
, c
4
> 0, independent of h, such that
for all v H
1
(

) there hold
|v I

h
(v)|
0,T
c
3
h
T
|v|
1,(T)
T T

h
,
and
|v I

h
(v)|
0,e
c
4
h
1/2
e
|v|
1,(e)
e c
h
,
where

(T) := T

h
: T

T = 0 and

(e) := T

h
: T

e = 0 .
14
Proof. See [19].
Finally, we require the technical results given by the following two lemmas.
Lemma 3.6 Let H
1
(
S
) and H
1
(
S
). Then there hold
[R
1
(
S
h
())[ c
1

1

TT
S
h
h
T
|
d
S,h
|
0,T
||
1,T
+ c
2

eE
h
()
h
1/2
e
|u
S,h
+
h
|
0,e
||
1,Te
,
and
[R
1
(curl ( I
S
h
()))[
c
3

1

TT
S
h
h
T
|rot (
d
S,h
)|
0,T
||
1,
S
(T)
+ c
4

1

eE
h
(
S
)
h
1/2
e
|[
d
S,h
t]|
0,e
||
1,
S
(e)
+c
4

1

eE
h
(
S
)
h
1/2
e
|
d
S,h
t|
0,e
||
1,
S
(e)
+ c
4

eE
h
()
h
1/2
e

d
S,h
t +
h
t

0,e
||
1,
S
(e)
.
Proof. We rst let :=
S
h
() and observe, according to (3.10) and (3.11), that

e
p n = 0 p [P
0
(e)]
2
, edge e of T
S
h
, and div = div P
S
h
(div ) .
Then, since
d
S,h
:
d
=
d
S,h
: and u
S,h
is a constant vector on each T T
S
h
, we deduce from
the denition of R
1
and the above identities that
R
1
() =
1

TT
S
h

d
S,h
:
d

TT
S
h

T
u
S,h
div

eE
h
()

h
n
=
1

TT
S
h

d
S,h
:

eE
h
()

h
n
=
1

TT
S
h

d
S,h
:

eE
h
()

e
(u
S,h
+
h
) n.
On the other hand, we now let := I
S
h
(). Then, using that div(curl ()) = 0, noting
that curl () n = t on , integrating by parts on each T T
S
h
and on , and observing
that
h
t L
2
(), we obtain
R
1
(curl ()) =
1

d
S,h
: curl () 'curl () n,
h
`

=
1

TT
S
h

T
rot (
d
S,h
) +

T

d
S,h
t

eE
h
()

e
(
h
t)
=

TT
S
h

T
rot (
d
S,h
) +

eE
h
(
S
)

e
[
d
S,h
t]
+

eE
h
(
S
)

e

d
S,h
t +

eE
h
()

d
S,h
+
h

t .
15
Hence, straighforward applications of the Cauchy-Schwarz inequality to the above equations,
together with the approximation properties provided by Lemmas 3.4 and 3.5, namely
|
S
h
()|
0,T
c
1
h
T
||
1,T
, | n
S
h
() n|
0,e
c
2
h
1/2
e
||
1,T
| I
S
h
()|
0,T
c
3
h
T
||
1,
S
(T)
, and | I
S
h
()|
0,e
c
4
h
1/2
e
||
1,
S
(e)
,
for each T T
S
h
and for each e c(T), imply the required estimates and nish the proof.
Lemma 3.7 Let w H
1
(
D
) and H
1
(
D
). Then there hold
[R
2
(w
D
h
(w))[ c
1

TT
D
h
h
T
|K
1
u
D,h
|
0,T
|w|
1,T
+ c
2

eE
h
()
h
1/2
e
|p
D,h

h
|
0,e
|w|
1,Te
,
and
[R
2
(curl ( I
D
h
()))[ c
3

TT
D
h
h
T
|rot (K
1
u
D,h
)|
0,T
||
1,
D
(T)
+ c
4

eE
h
(
D
)
h
1/2
e
|[K
1
u
D,h
t]|
0,e
||
1,
D
(e)
+ c
4

eE
h
()
h
1/2
e

K
1
u
D,h
t +
d
h
dt

0,e
||
1,
D
(e)
.
Proof. Since R
1
and R
2
have analogue structures, the proof proceeds similarly as for Lemma
3.6.
We are now in a position to bound the suprema depending on R
1
and R
2
.
Lemma 3.8 There exists C
1
> 0, independent of h, such that
sup

S
H(div;
S
)

S
=0
[R
1
(
S
)[
|
S
|
div;
S
C
1

TT
S
h

2
S,T

1/2
, (3.12)
where, for each T T
S
h
:

2
S,T
:= h
2
T
|rot
d
S,h
|
2
0,T
+ h
2
T
|
d
S,h
|
2
0,T
+

eE(T)E
h
(
S
)
h
e
|[
d
S,h
t]|
2
0,e
+

eE(T)E
h
(
S
)
h
e
|
d
S,h
t|
2
0,e
+

eE(T)E
h
()

h
e

d
S,h
t +
h
t

2
0,e
+ h
e
|u
S,h
+
h
|
2
0,e

Proof. Given
S
H(div;
S
) we know from Lemma 3.3 that there exist H
1
(
S
) and
H
1
(
S
) such that
S
= + curl in
S
and
||
1,
S
+ ||
1,
S
C |
S
|
div;
S
. (3.13)
16
Then, since R
1
(
S,h
) = 0
S,h
H
h
(
S
), which follows from the rst equation of the
Galerkin scheme (2.11) taking (v
D
, , ) = (0, 0, 0), and thanks to the fact that R
1
is linear,
we obtain
R
1
(
S
) = R
1
(
S

S,h
)
S,h
H
h
(
S
) . (3.14)
In particular, we let
S,h
:=
S
h
() + curl (I
S
h
()), which can be seen as a discrete Helmholtz
decomposition of
S,h
, and obtain
R
1
(
S
) = R
1
(
S
h
()) + R
1
(curl ( I
S
h
()) . (3.15)
Hence, applying Lemma 3.6 and then the discrete Cauchy-Schwarz inequality to the resulting
terms, noting that the numbers of triangles in
S
(T) and
S
(e) are bounded, and nally using
the estimate (3.13), we conclude the upper bound (3.12).
Lemma 3.9 There exists C
2
> 0, independent of h, such that
sup
v
D
H(div;
D
)
v
D
=0
[R
2
(v
D
)[
|v
D
|
div ;
D
C
2

TT
D
h

2
D,T

1/2
, (3.16)
where, for each T T
D
h
:

2
D,T
:= h
2
T
|rot (K
1
u
D,h
)|
2
0,T
+ h
2
T
|K
1
u
D,h
|
2
0,T
+

eE(T)E
h
(
D
)
h
e

[K
1
u
D,h
t]

2
0,e
+

eE(T)E
h
()

h
e

K
1
u
D,h
t +
d
h
dt

2
0,e
+ h
e
|p
D,h

h
|
2
0,e

.
Proof. It follows basically the same lines of the proof of Lemma 3.8. In fact, given v
D

H(div;
D
) we rst apply Lemma 3.3 to deduce the existence of w H
1
(
D
) and H
1
(
D
)
such that v
D
= w+ curl and
|w|
1,
D
+ ||
1,
D
C |v
D
|
div ;
D
. (3.17)
Then, since R
2
(v
D,h
) = 0 v
D,h
H
h
(
D
), which corresponds to the rst equation of the
Galerkin scheme (2.11) with (
S
, , ) = (0, 0, 0), and thanks to the fact that R
2
is linear, we
obtain
R
2
(v
D
) = R
2
(v
D
v
D,h
) v
D,h
H
h
(
D
) . (3.18)
Next, we choose v
D,h
=
D
h
(w) + curl

I
D
h
()

, notice that
R
2
(v
D
) = R
2
(w
D
h
(w)) + R
2

curl ( I
D
h
())

,
and apply Lemma 3.7. Thus, using again the discrete Cauchy-Schwarz inequality, noting that
the numbers of triangles in
D
(T) and
D
(e) are bounded, and employing now the upper bound
(3.17), we conclude (3.16).
We end this section by observing that the reliability estimate (3.2) (cf. Theorem 3.1) is a
direct consequence of Lemmas 3.1, 3.2, 3.8, and 3.9.
17
3.2 Eciency of the a posteriori error estimator
The main result of this section is stated as follows.
Theorem 3.2 There exists C
eff
> 0, independent of h, such that
C
eff
|
h
|
X
+ |u u
h
|
M
+ h.o.t. , (3.19)
where h.o.t. stands, eventually, for one or several terms of higher order.
We remark in advance that the proof of (3.19) makes frequent use of the identities provided
by Theorem 2.2. We begin with the estimates for the zero order terms appearing in the denition
of
2
S,T
and
2
D,T
.
Lemma 3.10 There hold
|f
S
+ div
S,h
|
0,T
|
S

S,h
|
div;T
T T
S,h
and
|f
D
div u
D,h
|
0,T
|u
D
u
D,h
|
div ;T
T T
D,h
.
Proof. It suces to recall, as established by Theorem 2.2, that f
S
= div
S
in
S
and
f
D
= div u
D
in
D
.
In order to derive the upper bounds for the remaining terms dening the global a posteriori
error estimator (cf. (3.1)), we proceed similarly as in [8], using results from [15], [17] and [25],
and apply Helmholtz decomposition, inverse inequalities, and the localization technique based
on triangle-bubble and edge-bubble functions. To this end, we now introduce further notations
and preliminary results. Given T T
S
h
, T
D
h
, T T , and e c(T), we let
T
and
e
be
the usual triangle-bubble and edge-bubble functions, respectively (see (1.5) and (1.6) in [46]).
In particular,
T
satises
T
P
3
(T), supp(
T
) T,
T
= 0 on T, and 0
T
1 in
T. Similarly,
e
[
T
P
2
(T), supp(
e
) w
e
:= T

T : e c(T

),
e
= 0 on T`e,
and 0
e
1 in w
e
. We also recall from [45] that, given k N 0, there exists an
extension operator L : C(e) C(T) that satises L(p) P
k
(T) and L(p)[
e
= p p P
k
(e).
A corresponding vector version of L, that is the componentwise application of L, is denoted by
L. Additional properties of
T
,
e
, and L are collected in the following lemma.
Lemma 3.11 Given k N 0, there exist positive constants c
1
, c
2
and c
3
, depending only
on k and the shape regularity of the triangulations (minimum angle condition), such that for
each triangle T and e c(T), there hold
|q|
2
0,T
c
1
|
1/2
T
q|
2
0,T
q P
k
(T), (3.20)
|q|
2
0,e
c
2
|
1/2
e
q|
2
0,e
q P
k
(e), (3.21)
and
|
1/2
e
L(q)|
2
0,T
c
3
h
e
|q|
2
0,e
q P
k
(e) . (3.22)
Proof. See Lemma 1.3 in [45].
The following inverse estimate for polynomials will also be used.
18
Lemma 3.12 Let k, l, m N 0 such that l m. Then, there exists c > 0, depending only
on k, l, m and the shape regularity of the triangulations, such that for each triangle T there
holds
[q[
m,T
c h
lm
T
[q[
l,T
, q P
k
(T) . (3.23)
Proof. See Theorem 3.2.6 in [18].
In addition, we need to recall a discrete trace inequality, which establishes the existence of a
positive constant c, depending only on the shape regularity of the triangulations, such that for
each T T
S
h
T
D
h
and e c(T), there holds
|v|
2
0,e
c

h
1
e
|v|
2
0,T
+ h
e
[v[
2
1,T

v H
1
(T) . (3.24)
For a proof of inequality (3.24) we refer to Theorem 3.10 in [1] (see also eq. (2.4) in [6]).
The following lemma summarizes known eciency estimates for ten terms dening
2
S,T
and

2
D,T
. In fact, their proofs, which apply the preliminary results described above, are already
available in the literature (see, e.g. [8], [9], [15], [25], [26], [28]). From now on we assume,
without loss of generality, that K
1
u
D,h
is polynomial on each T T
D
h
. Otherwise, additional
higher order terms, given by the errors arising from suitable polynomial approximations, should
appear in the corresponding bounds below, which explains the expression h.o.t. in (3.19).
Lemma 3.13 There exist positive constants C
i
, i 1, ..., 10, independent of h, such that
a) h
2
T
|rot (K
1
u
D,h
)|
2
0,T
C
1
|u
D
u
D,h
|
2
0,T
T T
D
h
,
b) h
2
T
|rot
d
S,h
|
2
0,T
C
2
|
S

S,h
|
2
0,T
T T
S
h
,
c) h
e
[[K
1
u
D,h
t]|
2
0,e
C
3
|u
D
u
D,h
|
2
0,we
e c
h
(
D
), where the set w
e
is given by
w
e
:=

T
D
h
: e c(T

,
d) h
e
|[
d
S,h
t]|
2
0,e
C
4
|
S

S,h
|
2
0,we
e c
h
(
S
), where the set w
e
is given by
w
e
:=

T
S
h
: e c(T

,
e) h
e
|
d
S,h
t|
2
0,e
C
5
|
S

S,h
|
2
0,T
e c
h
(
S
), where T is the triangle of T
S
h
having
e as an edge,
f) h
2
T
|K
1
u
D,h
|
2
0,T
C
6

|p
D
p
D,h
|
2
0,T
+ h
2
T
|u
D
u
D,h
|
2
0,T

T T
D
h
,
g) h
2
T
|
d
S,h
|
2
0,T
C
7

|u
S
u
S,h
|
2
0,T
+ h
2
T
|
S

S,h
|
2
0,T

T T
S
h
,
h) h
e
|p
D,h

h
|
2
0,e
C
8

|p
D
p
D,h
|
2
0,T
+ h
2
T
|u
D
u
D,h
|
2
0,T
+ h
e
|
h
|
2
0,e

e c
h
(),
where T is the triangle of T
D
h
having e as an edge,
i)

eE
h
()
h
e

K
1
u
D,h
t +
d
h
dt

2
0,e
C
9

eE
h
()
|u
D
u
D,h
|
2
0,Te
+ |
h
|
2
1/2,

,
where, given e c
h
(), T
e
is the triangle of T
D
h
having e as an edge, and
19
j)

eE
h
(
S
)
h
e

d
S,h
t +
h
t

2
0,e
C
10

eE
h
(
S
)
|
S

S,h
|
2
0,Te
+ |
h
|
2
1/2,

,
where, given e c
h
(
S
), T
e
is the triangle of T
S
h
having e as an edge.
Proof. For a) and b) we refer to [15, Lemma 6.1]. Alternatively, a) and b) follow from straight-
forward applications of the technical result provided in [9, Lemma 4.3] (see also [28, Lemma
4.9]). Similarly, for c), d), and e) we refer to [15, Lemma 6.2] or apply the technical result given
by [9, Lemma 4.4] (see also [28, Lemma 4.10]). Then, for f) and g) we refer to [15, Lemma
6.3] (see also [28, Lemma 4.13] or [25, Lemma 5.5]). On the other hand, the estimate given
by h) corresponds to [8, Lemma 4.12]. In particular, its proof makes use of the discrete trace
inequality (3.24). Finally, the proofs of i) and j) follow from very slight modications of the
proof of [25, Lemma 5.7]. Alternatively, an elasticity version of i) and j), which is provided in
[26, Lemma 20], can also be adapted to our case.
We nd it important to remark that the estimates i) and j) in the previous lemma provide the
only non-local bounds of the present eciency analysis. However, under additional regularity
assumptions on and , one is able to prove the following local bounds.
Lemma 3.14 Assume that [
e
H
1
(e) for each e c
h
(), and that [
e
H
1
(e) for each
e c
h
(
S
). Then there exist

C
9
,

C
10
> 0, such that
h
e

K
1
u
D,h
t +
d
h
dt

2
0,e

C
9

|u
D
u
D,h
|
2
0,Te
+ h
e

d
dt

2
0,e

e c
h
() ,
and
h
e

d
S,h
t +
h
t

2
0,e

C
10

|
S

S,h
|
2
0,Te
+ h
e

d
dt

2
0,e

e c
h
(
S
) .
Proof. Similarly as for i) and j) from Lemma 3.13, it follows by adapting the corresponding
elasticity version from [26]. We omit details here and refer to [26, Lemma 21].
It remains to provide the eciency estimates for three residual terms dened on the edges of
the interface . They have to do with the transmision conditions and with the trace equation
u
S
+ = 0 on . More precisely, we have the following lemmas.
Lemma 3.15 There exists C > 0, independent of h, such that for each e c
h
(), there holds
h
e
|u
D,h
n +
h
n|
2
0,e
C

|u
D
u
D,h
|
2
0,T
+ h
2
T
|div (u
D
u
D,h
)|
2
0,T
+ h
e
|
h
|
2
0,e

,
where T is the triangle of T
D
h
having e as an edge.
Proof. We proceed similarly as in [8, Lemma 4.7]. Given e c
h
(), we let T be the triangle of
T
D
h
having e as an edge, and dene v
e
:= u
D,h
n +
h
n on e. Then, applying (3.21), recalling
that
e
= 0 on T`e, extending
e
L(v
e
) by zero in
D
`T so that the resulting function belongs
to H
1
(
D
), and using that u
D
n + n = 0 on , we get
|v
e
|
2
0,e
c
2
|
1/2
e
v
e
|
2
0,e
= c
2

e
v
e
(u
D,h
n +
h
n) = c
2
'u
D,h
n +
h
n,
e
L(v
e
)`

= c
2
'u
D,h
n u
D
n,
e
L(v
e
)`

+ c
2
'
h
n n,
e
L(v
e
)`

,
(3.25)
20
where ', `

stands here for the duality pairing between H

1/2
() and H
1/2
(). Next, integra-
ting by parts in
D
, and noting that

h
n n

L
2
(), we nd, respectively, that
'u
D,h
n u
D
n,
e
L(v
e
)`

e
L(v
e
)

(u
D,h
u
D
) +

e
L(v
e
) div (u
D,h
u
D
) ,
and
'
h
n n,
e
L(v
e
)`

h
n n

e
v
e
.
Thus, replacing the above expressions back into (3.25), applying the Cauchy-Schwarz inequality
and the inverse estimate (3.23), and recalling that 0
e
1, we obtain
|v
e
|
2
0,e
C

h
1
T
|u
D
u
D,h
|
0,T
+ |div (u
D
u
D,h
)|
0,T

|
e
L(v
e
)|
0,T
+ c |v
e
|
0,e
|
h
|
0,e
.
But, using again that 0
e
1 and thanks to (3.22), we get
|
e
L(v
e
)|
0,T
|
1/2
e
L(v
e
)|
0,T
c
1/2
3
h
1/2
e
|v
e
|
0,e
, (3.26)
whence the previous inequality yields
|v
e
|
0,e
C h
1/2
e

h
1
T
|u
D
u
D,h
|
0,T
+ |div (u
D
u
D,h
)|
0,T

+ c |
h
|
0,e
.
Finally, it is easy to see that this estimate and the fact that h
e
h
T
imply the required upper
bound for h
e
|v
e
|
2
0,e
, which nishes the proof.
Lemma 3.16 There exists C > 0, independent of h, such that for each e c
h
(), there holds
h
e
|
S,h
n +
h
n

(
h
t) t|
2
0,e
C

|
S

S,h
|
2
0,T
+ h
2
T
|div(
S

S,h
)|
2
0,T
+ h
e
|
h
|
2
0,e
+ h
e
|
h
|
2
0,e

,
where T is the triangle of T
S
h
having e as an edge.
Proof. We proceed as in the previous lemma (see also [8, Lemma 4.6]). Indeed, given e c
h
(),
we let T be the triangle of T
S
h
having e as an edge, and dene v
e
:=
S,h
n +
h
n

(
h
t) t
on e. Then, applying (3.21), recalling that
e
= 0 on T`e, extending
e
L(v
e
) by zero in

S
`T so that the resulting function belongs to H
1
(
S
), using that
S
n + n

( t) t = 0
on , and then integrating by parts in
S
, we arrive at
|v
e
|
2
0,e
c
2
|
1/2
e
v
e
|
2
0,e
= c
2

e
v
e

S,h
n +
h
n

(
h
t) t

= c
2

T
(
e
L(v
e
)) : (
S,h

S
) + c
2

e
L(v
e
)) div(
S,h

S
)
+ c
2

e
v
e

(
h
) n

(
h
t t) t

.
Next, applying the Cauchy-Schwarz inequality and the inverse estimate (3.23), recalling that
0
e
1, and employing the vector version of (3.26), we deduce that
|v
e
|
0,e
C h
1/2
e

h
1
T
|
S

S,h
|
0.T
+ |div(
S

S,h
)|
0,T

+ C

|
h
|
0,e
+ |
h
|
0,e

,
which easily yields the required estimate, thus nishing the proof.
21
Lemma 3.17 There exists C > 0, independent of h, such that for each e c
h
(), there holds
h
e
|u
S,h
+
h
|
2
0,e
C

|u
S
u
S,h
|
2
0,T
+ h
2
T
|
S

S,h
|
2
0,T
+ h
e
|
h
|
2
0,e

,
where T is the triangle of T
S
h
having e as an edge.
Proof. Let e c
h
() and let T be the triangle of T
S
h
having e as an edge. We follow the proof
of [8, Lemma 4.12] and obtain rst an upper bound of h
2
T
[u
S
u
S,h
[
2
1,T
. Indeed, using that
u
S
=
1

d
S
in
S
(cf. Theorem 2.2) and that u
S,h
is constant in T, adding and substracting

d
S,h
, and then applying the estimate g) from Lemma 3.13, we deduce that
h
2
T
[u
S
u
S,h
[
2
1,T
=
h
2
T

2
|
d
S
|
2
0,T
C h
2
T

|
S

S,h
|
2
0,T
+ |
d
S,h
|
2
0,T

|u
S
u
S,h
|
2
0,T
+ h
2
T
|
S

S,h
|
2
0,T

.
(3.27)
Next, since = u
S
on (cf. Theorem 2.2), we nd that
h
e
|u
S,h
+
h
|
2
0,e
2 h
e

|u
S
u
S,h
|
2
0,e
+ |
h
|
2
0,e

,
which, employing the discrete trace inequality (3.24) and the estimate (3.27), yields
h
e
|u
S,h
+
h
|
2
0,e
C

|u
S
u
S,h
|
2
0,T
+ h
2
T
[u
S
u
S,h
[
2
1,T
+ h
e
|
h
|
2
0,e

|u
S
u
S,h
|
2
0,T
+ h
2
T
|
S

S,h
|
2
0,T
+ h
e
|
h
|
2
0,e

,
which completes the proof.
We end this section by observing that the eciency estimate (3.19) follows straightforwardly
from Lemmas 3.10, 3.13, 3.15, 3.16, and 3.17. In particular, the terms h
e
|
h
|
2
0,e
and
h
e
|
h
|
2
0,e
, which appear in Lemma 3.13 (item h)), 3.15, 3.16, and 3.17, are bounded as
follows:

eE
h
()
h
e
|
h
|
2
0,e
h|
h
|
2
0,
C h|
h
|
2
1/2,
,
and

eE
h
()
h
e
|
h
|
2
0,e
h|
h
|
2
0,
C h|
h
|
2
1/2,
.
4 Numerical results
In [31, Section 5] we presented several numerical results illustrating the performance of the
Galerkin scheme (2.11) with the subspaces X
h
:= H
h
(
S
) H
h
(
D
)
h
()
h
() and
M
h
:= L
h
(
S
) L
h,0
(
D
) dened in Section 2.3. We now provide three examples conrming
the reliability and eciency of the respective a posteriori error estimator derived in Section
3, and showing the behaviour of the associated adaptive algorithm.
22
In what follows, N stands for the number of degrees of freedom dening X
h
and M
h
, and
the individual and total errors are dened by:
e(
S
) := |
S

S,h
|
div,
S
, e(u
S
) := |u
S
u
S,h
|
div ;
S
,
e(u
D
) := |u
D
u
D,h
|
div ;
D
, e(p
D
) := |p
D
p
D,h
|
0,
D
,
e() := |
h
|
1/2,
, e() := |
h
|
1/2,
,
and
e(, u) :=

(e(
S
))
2
+ (e(u
S
))
2
+ (e(u
D
))
2
+ (e(p
D
))
2
+ (e())
2
+ (e())
2

1/2
,
whereas the eectivity index with respect to is given by
eff() := e(, u)/,
where
(, u) := ((
S
, u
D
, , ), (u
S
, p
D
)) X M
and
(
h
, u
h
) := ((
S,h
, u
D,h
,
h
,
h
), (u
S,h
, p
D,h
)) X
h
M
h
denote the unique solutions of (2.5) and (2.11), respectively.
Also, we let r(
S
), r(u
S
), r(u
D
), r(p
D
), r(), r(), and r(, u) be the individual and global
experimental rates of convergence given by
r(%) :=
log(e(%)/e

(%))
log(h/h

)
for each %

S
, u
S
, u
D
, p
D
, ,

,
and
r(, u) :=
log(e(, u)/e

(, u))
log(h/h

)
,
where h and h

denote two consecutive meshsizes with errors e and e

. However, when the

adaptive algorithm is applied (see details below), the expression log(h/h

) appearing in the
computation of the above rates is replaced by
1
2
log(N/N

), where N and N

denote the
corresponding degrees of freedom of each triangulation.
The examples to be considered in this section are described next. In all of them we choose
for simplicity = 1, K = I, the identity matrix of R
22
, and = 1. Example 1 is employed
to conrm the reliability and eciency of the a posteriori error estimator . Then, Examples 2
and 3 are utilized to illustrate the behaviour of the associated adaptive algorithm, which applies
the following procedure from [46]:
1) Start with a coarse mesh T
h
:= T
D
h
T
S
h
.
2) Solve the discrete problem (2.11) for the actual mesh T
h
.
3) Compute
,T
for each triangle T T

h
, D, S.
4) Evaluate stopping criterion and decide to nish or go to next step.
23
5) Use blue-green procedure to rene each T

h
, D, S, whose indicator
,T
satises

T

1
2
max
i{D,S}

max

i,T
: T T
i
h

.
6) Dene resulting meshes as actual meshes T
D
h
and T
S
h
, and go to step 2.
In Example 1 we consider the regions
D
:=] 0.5, 0.5[
2
and
S
:=] 1, 1[
2
`

D
, which
yields a porous medium completely surrounded by a uid, and choose the data f
S
and f
D
so
that the exact solution is given by the regular functions
u
S
(x) =

2 sin
2
(x
1
) sin(x
2
) cos(x
2
)
2 sin(x
1
) sin
2
(x
2
) cos(x
1
)

x := (x
1
, x
2
)
S
,
p
S
(x) = x
3
1
e
x
2
x := (x
1
, x
2
)
S
,
and
p
D
(x) = x
3
1
sin(x
2
) x := (x
1
, x
2
)
D
.
In Example 2 we consider
D
:=] 1, 0[
2
and let
S
be the L-shaped domain given by
] 1, 1[
2
`

D
, which yields a porous medium partially surrounded by a uid. Then we choose
the data f
S
and f
D
so that the exact solution is given by
u
S
(x) = curl

0.1

x
2
2
1

2
sin
2
(x
1
)

x := (x
1
, x
2
)
S
,
p
S
(x) =
1
100 (x
2
1
+ x
2
2
) + 0.1
x := (x
1
, x
2
)
S
,
and
p
D
(x) =

x
1
+ 1
10

2
sin
3
(2 (x
2
+ 0.5)) x := (x
1
, x
2
)
D
.
Note that the uid pressure p
S
has high gradients around the origin.
Finally, in Example 3 we take
D
:=] 1, 1[ ] 2, 1[ and
S
:=] 1, 1[
2
` [0, 1]
2
, which
yields a porous medium below a uid, and choose the data f
S
and f
D
so that the exact solution
is given by
u
S
(r, ) = curl

0.1 r
5/3
(r
2
cos
2
() 1)
2
(r sin() 1)
2
sin
2

2
3

(r, )
S
,
p
S
(x) = 0.1 x
1
sin(x
2
) x := (x
1
, x
2
)
S
,
and
p
D
(x) = 0.1 (x
2
+ 2)
2
sin
3
(x
1
) x := (x
1
, x
2
)
D
.
Note that u
S
is dened in polar coordinates and that its derivatives are singular at the origin.
The numerical results shown below were obtained using a MATLAB code. In Table 4.1 we
summarize the convergence history of the mixed nite element method (2.11), as applied to
Example 1, for a sequence of quasi-uniform triangulations of the domain. We observe there,
24
Table 4.1: Example 1, quasi-uniform scheme
N h e(
S
) r(
S
) e(u
S
) r(u
S
) e(u
D
) r(u
D
) e(p
D
) r(p
D
)
321 0.5000 35.4015 0.6875 0.1996 0.0117
1201 0.2500 20.0107 0.8647 0.4266 0.7234 0.1121 0.8743 0.0057 1.0798
4641 0.1250 10.0700 1.0160 0.1615 1.4370 0.0531 1.1046 0.0023 1.3213
18241 0.0625 5.0492 1.0087 0.0801 1.0238 0.0259 1.0490 0.0011 1.0967
72321 0.0312 2.5268 1.0052 0.0401 1.0064 0.0129 1.0178 0.0005 1.0234
288001 0.0156 1.2637 1.0029 0.0200 1.0031 0.0064 1.0062 0.0003 1.0062
N h e() r() e() r() e(, u) r(, u) eff()
321 0.5000 4.2653 0.0981 35.6649 39.0015 0.9144
1201 0.2500 4.3919 0.0973 0.0124 20.4920 0.8399 22.6847 0.9033
4641 0.1250 1.7410 1.3690 0.0537 0.8781 10.2209 1.0292 11.1965 0.9129
18241 0.0625 0.8088 1.1202 0.0259 1.0670 5.1144 1.0117 5.5954 0.9140
72321 0.0312 0.3949 1.0408 0.0126 1.0516 2.5579 1.0060 2.7969 0.9145
288001 0.0156 0.1962 1.0123 0.0062 1.0266 1.2791 1.0031 1.3982 0.9148
looking at the corresponding experimental rates of convergence, that the O(h) predicted by
Theorem 2.4 (when = 1) is attained in all the unknowns. In addition, we notice that the
eectivity index eff() remains always in a neighborhood of 0.91, which illustrates the reliability
and eciency of in the case of a regular solution.
Next, in Tables 4.2 - 4.5 we provide the convergence history of the quasi-uniform and adaptive
schemes, as applied to Examples 2 and 3. We observe that the errors of the adaptive procedures
decrease faster than those obtained by the quasi-uniform ones, which is conrmed by the global
experimental rates of convergence provided there. This fact is also illustrated in Figures 4.1 and
4.3 where we display the total errors e(, u) vs. the degrees of freedom N for both renements.
As shown by the values of r(, u), the adaptive method is able to keep the quasi-optimal rate of
convergence O(h) for the total error. Furthermore, the eectivity indexes remain bounded from
above and below, which conrms the reliability and eciency of in these cases of non-smooth
solutions. Intermediate meshes obtained with the adaptive renements are displayed in Figures
4.2 and 4.4. Note that the method is able to recognize the region with high gradients in Example
2, and the singularity of the solution in Example 3.
References
[1] S. Agmon, Lectures on Elliptic Boundary Value Problems. Van Nostrand, Princeton, New
Jersey, 1965.
[2] M. Ainsworth and J.T. Oden, A unied approach to a posteriori error estimation based
on element residual methods. Numerische Mathematik, vol. 65, pp. 23-50, (1993).
[3] M. Ainsworth and J.T. Oden, A posteriori error estimators for the Stokes and Oseen
equations. SIAM Journal on Numerical Analysisi, vol. 34, 1, pp. 228-245, (1997).
25
Table 4.2: Example 2, quasi-uniform scheme
N h e(
S
) e(u
S
) e(u
D
) e(p
D
) e() e()
608 0.3536 4.5187 0.1198 0.2649 0.0184 0.5760 0.1120
2332 0.1768 4.9963 0.0529 0.1520 0.0035 0.2653 0.0347
9140 0.0884 6.7481 0.0253 0.0778 0.0005 0.1485 0.0096
36196 0.0442 4.2857 0.0125 0.0392 0.0002 0.0771 0.0042
144068 0.0221 2.4834 0.0062 0.0196 0.0001 0.0348 0.0022
N h e(, u) r(, u) eff()
608 0.3536 4.5660 5.4033 0.8450
2332 0.1768 5.0060 5.2805 0.9480
9140 0.0884 6.7503 6.8230 0.9894
36196 0.0442 4.2866 0.6599 4.3158 0.9932
144068 0.0221 2.4837 0.7901 2.4958 0.9952
1
10
100 1000 10000 100000
N
quasi-uniform renement
3
3
3
3
3
3
adaptive renement
+
+
+
+
+
+
+
+
+
+
+
+
+
+
Figure 4.1: Example 2, e(, u) vs. N for quasi-uniform/adaptive schemes
26
Table 4.3: Example 2, adaptive scheme
N e(
S
) e(u
S
) e(u
D
) e(p
D
) e() e()
608 4.5188 0.1199 0.2649 0.0184 0.5760 0.1121
1118 5.3792 0.0709 0.2262 0.0091 0.3185 0.0402
1391 7.2290 0.0661 0.2098 0.0082 0.2846 0.0215
1636 5.1151 0.0657 0.2094 0.0110 0.2591 0.0236
1884 3.9177 0.0657 0.2093 0.0108 0.2577 0.0229
3558 2.6519 0.0491 0.2020 0.0037 0.1626 0.0128
7164 1.8814 0.0320 0.1751 0.0067 0.1160 0.0171
13073 1.3945 0.0237 0.1591 0.0034 0.0742 0.0109
26227 0.9771 0.0165 0.1222 0.0030 0.0730 0.0103
35611 0.8163 0.0140 0.1089 0.0018 0.0384 0.0075
55318 0.6608 0.0114 0.0808 0.0005 0.0375 0.0039
70434 0.5825 0.0099 0.0747 0.0005 0.0357 0.0038
149402 0.4052 0.0070 0.0548 0.0003 0.0208 0.0023
N e(, u) r(, u) eff()
608 4.5660 5.4033 0.8450
1118 5.3940 5.7977 0.9304
1391 7.2379 7.4956 0.9656
1636 5.1264 4.2524 5.4334 0.9435
1884 3.9324 3.7572 4.3145 0.9114
3558 2.6650 1.2238 2.9662 0.8985
7164 1.8934 0.9768 2.0913 0.9054
13073 1.4057 0.9902 1.5394 0.9132
26227 0.9876 1.0142 1.0951 0.9018
35611 0.8246 1.1796 0.9191 0.8972
55318 0.6669 0.9637 0.7388 0.9026
70434 0.5885 1.0359 0.6505 0.9046
149402 0.4095 0.9644 0.4550 0.8999
27
Figure 4.2: Example 2, adapted meshes with 1884, 7164, 26227, and 55318 degrees of freedom
Table 4.4: Example 3, quasi-uniform scheme
N h e(
S
) e(u
S
) e(u
D
) e(p
D
) e() e()
344 0.5000 16.8563 0.4452 0.7130 0.0674 1.8109 0.1615
1324 0.2500 11.3317 0.3329 0.3846 0.0130 2.5160 0.0826
5204 0.1250 7.0011 0.0849 0.1980 0.0038 0.8665 0.0458
20644 0.0625 4.4530 0.0412 0.0992 0.0018 0.3859 0.0203
82244 0.0312 2.8037 0.0206 0.0496 0.0009 0.1877 0.0097
N h e(, u) r(, u) eff()
344 0.5000 16.9751 18.8901 0.8986
1324 0.2500 11.6191 0.5626 13.1132 0.8861
5204 0.1250 7.0579 0.7284 7.8041 0.9044
20644 0.0625 4.4711 0.6626 5.0014 0.8940
82244 0.0312 2.8105 0.6717 3.1653 0.8879
28
Table 4.5: Example 3, adaptive scheme
N e(
S
) e(u
S
) e(u
D
) e(p
D
) e() e()
344 16.8564 0.4453 0.7131 0.0675 1.8109 0.1616
684 11.8048 0.3406 0.5828 0.0177 2.5165 0.0863
1367 10.6242 0.1330 0.4426 0.0099 0.8682 0.0530
1625 10.4486 0.1314 0.4426 0.0099 0.8682 0.0530
1863 10.3440 0.1278 0.4426 0.0098 0.8678 0.0530
2291 9.2480 0.1173 0.4427 0.0097 0.8672 0.0526
3109 7.5456 0.1013 0.4425 0.0098 0.8670 0.0522
11719 3.9053 0.0530 0.3296 0.0072 0.3875 0.0271
34611 2.2713 0.0202 0.2614 0.0058 0.1901 0.0092
60159 1.7281 0.0153 0.1723 0.0034 0.1759 0.0083
79482 1.5031 0.0111 0.1644 0.0032 0.1154 0.0072
115241 1.2620 0.0167 0.1498 0.0019 0.1101 0.0055
182014 0.9954 0.0130 0.1226 0.0012 0.0900 0.0027
N e(, u) r(, u) eff()
344 16.9751 18.8901 0.8986
684 12.0893 0.9877 13.6112 0.8882
1367 10.6698 0.3608 11.3264 0.9420
1625 10.4949 0.1912 11.1221 0.9436
1863 10.3907 0.1460 10.8244 0.9599
2291 9.3000 1.0724 9.9113 0.9383
3109 7.6090 1.3146 8.2092 0.9269
11719 3.9388 1.0924 4.2413 0.9362
34611 2.2943 0.9981 2.4691 0.9292
60159 1.7456 0.9889 1.8902 0.9235
79482 1.5165 1.0102 1.5941 0.9513
115241 1.2757 0.9309 1.3418 0.9507
182014 1.0070 1.0350 1.0817 0.9309
29
1
10
100 1000 10000 100000
N
quasi-uniform renement
3
3
3
3
3
3
adaptive renement
+
+
+++
+
+
+
+
+
+
+
+
+
Figure 4.3: Example 3, e(, u) vs. N for quasi-uniform/adaptive schemes
[4] A. Alonso, Error estimators for a mixed method. Numerische Mathematik, vol. 74, pp.
385-395, (1996).
[5] T. Arbogast and D.S. Brunson, A computational method for approximating a Darcy-
Stokes system governing a vuggy porous medium. Computational Geosciences, vol. 11, 3,
pp. 207-218, (2007).
[6] D.N. Arnold, An interior penalty nite element method with discontinuous elements.
SIAM Journal on Numerical Analysis, vol. 19, 4, pp. 742-760, (1982).
[7] I. Babuska and G.N. Gatica, On the mixed nite element method with Lagrange multi-
pliers. Numerical Methods for Partial Dierential Equations, vol. 19, 2, pp. 192-210, (2003).
[8] I. Babuska and G.N. Gatica, A residual-based a posteriori error estimator for the
Stokes-Darcy coupled problem. SIAM Journal on Numerical Analysis, vol. 48, 2, pp. 498-523,
(2010).
[9] T.P. Barrios, G.N. Gatica, M. Gonz alez, and N. Heuer, A residual based a pos-
teriori error estimator for an augmented mixed nite element method in linear elasticity.
ESAIM: Mathematical Modelling and Numerical Analysis, vol. 40, 5, pp. 843-869, (2006).
[10] C. Bernardi, F. Hecht, and O. Pironneau, Coupling Darcy and Stokes equations for
porous media with cracks. ESAIM: Mathematical Modelling and Numerical Analysis, vol.
39, 1, pp. 7-35, (2005).
[11] D. Braess and R. Verf urth, A posteriori error estimators for the Raviart-Thomas
element. SIAM Journal on Numerical Analysis, vol. 33, pp. 24312444, (1996).
30
Figure 4.4: Example 3, adapted meshes with 1863, 3109, 11719, and 60159 degrees of freedom
31
[12] F. Brezzi and M. Fortin, Mixed and Hybrid Finite Element Methods. Springer Verlag,
1991.
[13] E. Burman and P. Hansbo, Stabilized Crouzeix-Raviart elements for the Darcy-Stokes
problem. Numerical Methods for Partial Dierential Equations, vol. 21, 5, pp. 986-997,
(2005).
[14] E. Burman and P. Hansbo, A unied stabilized method for Stokes and Darcys equations.
Journal of Computational and Applied Mathematics, vol. 198, 1, pp. 35-51, (2007).
[15] C. Carstensen, A posteriori error estimate for the mixed nite element method. Mathe-
matics of Computation, vol. 66, 218, pp. 465-476, (1997).
[16] C. Carstensen, An a posteriori error estimate for a rst kind integral equation. Mathe-
matics of Computation, vol. 66, 217, pp. 139-155, (1997).
[17] C. Carstensen and G. Dolzmann, A posteriori error estimates for mixed FEM in
elasticity. Numerische Mathematique, vol. 81, pp. 187-209, (1998).
[18] P. G. Ciarlet, The nite Element Method for Elliptic Problems. North-Holland, Amster-
dam, New York, Oxorfd, 1978.
[19] P. Clement, Approximation by nite element functions using local regularisation. RAIRO
Modelisation Mathematique et Analyse Numerique, vol. 9, pp. 77-84, (1975).
[20] M.R. Correa, Stabilized Finite Element Methods for Darcy and Coupled Stokes-Darcy
Flows. D.Sc. Thesis, LNCC, Petropolis, Rio de Janeiro, Brasil (in portuguese), (2006).
[21] M.R. Correa and A.F.D. Loula, A unied mixed formulation naturally coupling Stokes
and Darcy ows. Computer Methods in Applied Mechanics and Engineering, vol. 198, 33-36,
pp. 2710-2722, (2009).
[22] M. Discacciati, E. Miglio, and A. Quarteroni, Mathematical and numerical models
for coupling surface and groundwater ows. Applied Numerical Mathematics, vol. 43, pp.
57-74, (2002).
[23] V.J. Ervin, E.W. Jenkins, and S. Sun, Coupled generalized nonlinear Stokes ow with
ow through a porous medium. SIAM Journal on Numerical Analysis, vol. 47, 2, pp. 929-952,
(2009).
[24] J. Galvis and M. Sarkis, Non-matching mortar discretization analysis for the coupling
Stokes-Darcy equations. Electronic Transactions on Numerical Analysis, vol. 26, pp. 350-
384, (2007).
[25] G.N. Gatica, A note on the eciency of residual-based a-posteriori error estimators for
some mixed nite element methods. Electronic Transactions on Numerical Analysos, vol 17,
pp. 218-233, (2004).
[26] G.N. Gatica, G.C. Hsiao, and S. Meddahi, A residual-based a posteriori error esti-
mator for a two-dimensional uid-solid interaction problem. Numerische Mathemaik, vol.
114, 1, pp. 63-106, (2009).
32
[27] G.N. Gatica and M. Maischak, A posteriori error estimates for the mixed nite element
method with Lagrange multipliers. Numerical Methods for Partial Dierential Equations,
vol. 21, 3, pp. 421-450, (2005).
[28] G.N. Gatica, A. M arquez, and M.A. S anchez, Analysis of a velocity-pressure-
pseudostress formulation for the stationary Stokes equations. Computer Methods in Applied
Mechanics and Engineering, vol. 199, 17-20, pp. 1064-1079, (2010).
[29] G.N. Gatica, S. Meddahi, and R. Oyarz ua, A conforming mixed nite-element method
for the coupling of uid ow with porous media ow. IMA Journal of Numerical Analysis,
vol. 29, 1, pp. 86-108, (2009).
[30] G.N. Gatica, R. Oyarz ua and F-J Sayas, Convergence of a family of Galerkin dis-
cretizations for the Stokes-Darcy coupled problem. Numerical Methods for Partial Dieren-
tial Equations DOI 10.1002/num, to appear.
[31] G.N. Gatica, R. Oyarz ua, and F.-J. Sayas, Analysis of fully-mixed nite element
methods for the Stokes-Darcy coupled problem. Preprint 2009-08, Departamento de Inge-
niera Matematica, Universidad de Concepcion, Chile, (2009).
[32] V. Girault and P. A. Raviart, Finite Element Methods for Navier-Stokes Equations.
Theory and Algorithms. Springer Series in Computational Mathematics, vol. 5, Springer-
Verlag, 1986.
[33] R.H.W. Hoppe and B.I. Wohlmuth, A comparison of a posteriori error estimators for
mixed nite element discretizations by Raviart-Thomas elements. Mathematics of Compu-
tation, vol. 68, 228, pp. 1347-1378, (1999).
[34] R. Kress, Linear Integral Equations. Springer Verlag, Berlin, 1989.
[35] W.J. Layton, F. Schieweck, and I. Yotov, Coupling uid ow with porous media
ow. SIAM Journal on Numerical Analysis, vol. 40, 6, pp. 2195-2218, (2003).
[36] M. Lonsing and R. Verf urth, A posteriori error estimators for mixed nite element
methods in linear elasticity. Numerische Mathematik, vol. 97, 4, pp. 757-778, (2004).
[37] C. Lovadina and R. Stenberg, Energy norm a posteriori error estimates for mixed nite
element methods. Mathematics of Computation, vol. 75, 256, pp. 1659-1674, (2006).
[38] A. Masud, A stabilized mixed nite element method for Darcy-Stokes ow. International
Journal for Numerical Methods in Fluids, vol. 54, 6-8, pp. 665-681, (2008).
[39] S. Repin, S. Sauter, and A. Smolianski, Two-sided a posteriori error estimates for
mixed formulations of elliptic problems. SIAM Journal on Numerical Analysis, vol. 45, 3,
pp. 928-945, (2007).
[40] B. Riviere, Analysis of a discontinuous nite element method for coupled Stokes and Darcy
problems. Journal of Scientic Computing, vol. 22-23, pp. 479-500, (2005).
[41] B. Riviere and I. Yotov, Locally conservative coupling of Stokes and Darcy ows. SIAM
Journal on Numerical Analysis, vol. 42, 5, pp. 1959-1977, (2005).
33
[42] H. Rui and R. Zhang, A unied stabilized mixed nite element method for coupling Stokes
and Darcy ows. Computer Methods in Applied Mechanics and Engineering, vol. 198, 33-36,
pp. 2692-2699, (2009).
[43] J.M. Urquiza, D. NDri, A. Garon, and M.C. Delfour, Coupling Stokes and Darcy
equations. Applied Numerical Mathematics, vol. 58, 5, pp. 525-538, (2008).
[44] R. Verf urth, A posteriori error estimators for the Stokes problem. Numerische Mathe-
matik, vol. 55, pp. 309-325, (1989).
[45] R. Verf urth, A posteriori error estimation and adaptive mesh-renement thecniques.
Journal of Computational and Applied Mathematics, vol 50, pp. 67-83, (1994).
[46] R. Verf urth, A Review of A posteriori Error Estimation and Adaptive Mesh-Renement
Thecniques. Wiley-Teubner (Chichester), 1996.
[47] X. Xie, J. Xu, and G. Xue, Uniformly stable nite element methods for Darcy-Stokes-
Brinkman models. Journal of Computational Mathematics, vol. 26, 3, pp. 437-455, (2008).
34