Polynomia l -degree-robust a po steriori estimates in a uni fied setting
for conforming, nonconforming, dis co nt inuous Galerkin, and mixed
discretizations
∗
Alexandre Ern
†
Martin Vohral´ık
‡
July 30, 2014
Abstract
We present equilibrated fl ux a posteriori error estimates in a un ifi ed setting for conforming, noncon-
forming, discont inuous Galerkin, and mixed finite element discretizations of the two-dimensional Poisson
problem. Relying on the equilibration by mixed finite element solution of patchwise Neumann problems,
the estimates are guaranteed, locally comput ab le, locally efficient, and robust with respect to polynomial
degree. Maximal local overestimation is guaranteed as well. Nu merical experiments suggest asymptotic
exactness for th e incomplete interior penalty discontinuous Galerkin scheme.
Key words: a posteriori error estimate, equilibrated flux, unified framework, robustness, p olynomial degree,
conforming finite element method, nonconforming finite element method, discontinuous Galer kin method,
mixed finite element method
1 Introdu c tion
A posteriori error estimates in the conforming finite element setting have already received a large attention.
In pa rticular, following the concept of Prager and Synge [
64], cf. also Synge [72], Aubin and Burchard [13],
and Hlav´aˇcek et al. [
50], and invoking fluxes in the H(div, Ω) spa ce, guaranteed upper bounds on the error
can be obtained. A g eneral functional framework delivering guaranteed upper bounds, independent of the
numerical method, has been der ived by Repin [
66, 67, 68]. It does not rely on Galerkin orthogonality
neither on local equilibration and accommodates an arbitra ry flux reco nstruction. The idea of us ing a local
residual equilibration procedure for the normal fa c e fluxes rec onstruction has been pr oposed by Ladev`eze [55],
Ladev`eze and Le guillon [
56], Kelly [51], Ainsworth and Oden [7, 8], and Par´es et al. [61, 62]. In this context,
guaranteed upper bounds typically require solving infinite-dimensional element problems , which, in practice,
are approximated. On the other hand, an essential proper ty achieved by means of local equilibration
procedures is local efficiency, meaning that the derived estimators also represent local lower bounds of the
error, up to a generic constant. This appears to be crucial in view of local mesh refinement, as well as
in order to obtain robustness in singularly perturbed problems. Cheap local flux equilibrations leading to
a fully computable guaranteed upper bound have been obtained by Destuynder and M´etivet [
38]. Later,
mixed finite element solutions of local Neumann problems pose d over patches of (sub)elements, where one
minimizes locally the estimator contributions, were prop osed, see Luce and Wohlmuth [
57], Braess and
Sch¨oberl [19], and [77, 30, 79]. As a matter of fact, lifting the normal face fluxes of the equilibrated residual
method to H(div, Ω) immediately yields equilibra ted fluxes , cf. Nicaise et al. [59]. T he n, both a guaranteed
bound and local efficiency are obtained. For computational comparisons o f some of these approaches in the
lowest-order case, s e e Carstensen and Merdon [
28].
The theory in the nonc onforming setting, where the discrete solution (potential) is not in the energy
space H
1
(Ω), appears to be less developed. First contributions are those of Agouzal [
1] and Dari et al.
∗
This work was s upported by the ERC-CZ project MORE “MOdelling REvisited + MOdel REduction” LL1202.
†
Universit´e Paris-Est, CERMICS (ENPC), 77455 Marne-la-Vall´ee cedex 2, France (
ern@cermics.enpc.fr).
‡
INRIA Paris-Rocquencourt, B.P. 105, 78153 Le Chesnay, France (
martin.vohralik@inria.fr).
1
[35], whereas a guaranteed error upper bound in the lowest-order Crouzeix–Raviart c ase can be obtained
along the lines of Destuynder and M´etivet [37], se e Ainsworth [2], Kim [52, 53], or [76]. Different flux
equilibrations exist and tight links hold between them, see [46]. Higher-order methods have been treated
by Ainsworth and Rankin [
9], and a survey and a computatio nal comparison in the lowest-order case can
be found in Carstensen and Merdon [29]. For the discontinuous Galerkin method, fir st guaranteed upper
bounds by loc ally equilibrated fluxe s and H
1
(Ω)-conforming reconstructed potentials have been obtained by
Ainsworth [
3], Kim [52, 53], Cochez-Dhondt and Nicaise [3 2], Ainsworth and Rankin [10], and [41, 43, 42],
see also the refer ences therein. Similar results for mixed finite elements can be found in Kim [52, 54],
Ainsworth [
4], Ainsworth and Ma [6], and [76, 78]. All the cited references typically prove local efficiency
as well.
When the flux equilibration is achieved by mixed finite element solutions of local Neumann problems,
the local efficiency result, in the conforming finite element setting, can be sharp ened by showing that the
efficiency constant is independent of the underlying polynomial degree. This important result was recently
proven by Braes s et al. [
18], and we refer to it as polynomial-degree robustness. This robustness prop erty
stands in contrast to the class of usua l residual-based estimators, cf. Verf¨urth [74], which yield local efficiency,
but where such a robustness does not hold, see Melenk and Wohlmuth [58]. The first key ingredient for
the proof in [
18] are continuous-level problems on patches of elements around vertices fea turing the hat
functions, similar to those considered already in Carstensen and Funken [24]. This replaces the usua l
bubble function technique. The second key ingredient is the polynomial-degree-robust stability of mixed
finite elements of [
18, Theorem 7], hinging on the poly nomial-degree robust elementwise construction of a
right inve rse of the divergence operator in polynomial spaces by Costabel and McIntosh [34] and on the
polynomial extension operators by Demkowicz et al. [
36].
We finally mention that unified frameworks for different discretization methods have been conce ived
recently, see Carstense n et al. [22, 27, 26, 23], Ainsworth [5], and, using the eq uilibrated fluxes, in [44, 4 9, 45].
In the present paper, we unify the potentials and equilibrated fluxes approach for most standard dis-
cretization schemes, including conforming, nonconfo rming (where the potential interface jumps satisfy some
orthogonality conditions ), discontinuous Galerkin, and mixed finite elements. The construction of the es-
timators becomes method independent, being close to that of Destuynder and M´etivet [
38] and coinciding
with that of Braes s and Sch¨oberl [19] for fluxes in the confor ming case, while being closely related to that
of Carstensen and Merdon [
29] for potentials in the nonconforming case. In the discontinuous Galerkin and
mixed finite element cases, s uch an approach appears to be new. The potentials and fluxes a re actually
constructed by the same patchwise problems with different right-hand sides in the present two-dimensional
setting. Most importantly, we prove the polynomial-degree robustness in this unified setting comprising all
the discussed discretization schemes. Moreover, we can also guarantee a maximal overestimation factor, a
feature which can be important in optimal convergence proofs. Additionally, numerical experiments for the
incomplete interior penalty discontinuous Gale rkin scheme suggest asymptotic exactness.
The paper is or ganized as follows: The setting is des crib e d in Section
2. The main results together with
their proofs are collected in Section
3. Applications to most standard numerical methods are showcased
in Section 4, and a numerical illustration is presented in Section 5. Concluding remarks in Section 6 close
the paper.
2 Setting
We sta rt by intro ducing the c ontinuous and discrete se ttings.
2.1 Sobolev spaces
Let Ω ⊂ R
2
be a polygonal domain (open, bounded, and connected set). We denote by H
1
(Ω) the Sobolev
space of L
2
(Ω) functions with weak gradients in [L
2
(Ω)]
2
and by H
1
0
(Ω) its zero-trace s ubspace. H(div, Ω)
stands for the space of [L
2
(Ω)]
2
functions with weak divergences in L
2
(Ω). The no tations ∇ and ∇· are
used respectively for the weak gradient and divergence. Le t R
π
2
:=
0 −1
1 0
be the matr ix of rotation by
π
2
; then R
π
2
∇ sta nds for the weak curl, i.e., the rotated gradient: for v ∈ H
1
(Ω), R
π
2
∇v = (−∂
y
v, ∂
x
v). For
a subdoma in ω of Ω, we denote by (·, ·)
ω
the L
2
(ω)-inner product, by k·k
ω
the asso ciated norm (we omit the
2
index when ω = Ω), and by |ω| the Leb esgue measure of ω. For ω ⊂ R
1
, h·, ·i
ω
stands for the 1-dimensional
L
2
(ω)-inner product or for the appropriate duality pairing on ω.
2.2 Meshes
We consider partitions T
h
of Ω which consist either of closed triangles or of closed rectangles K such that
Ω =
S
K∈T
h
K. We suppose that T
h
is matching , i.e., such that for two distinct elements, their intersection
is e ither an empty set or a co mmon edge or a co mmon vertex. For any K ∈ T
h
, n
K
stands for the outward
unit normal vector to K and h
K
denotes the diameter of K. The edges of the mesh form the set E
h
divided
into interior edges E
int
h
and boundary edges E
ext
h
. A g eneric edge is denoted by e and its diameter by h
e
.
For any e ∈ E
h
, n
e
stands for the unit normal vector to e; the orientation is arbitrary but fixed fo r e ∈ E
int
h
and points outwards of Ω for e ∈ E
ext
h
. The set of vertices is denoted by V
h
; it is decomposed into interior
vertices V
int
h
and boundary vertices V
ext
h
. For a ∈ V
h
, T
a
denotes the patch of the elements of T
h
which share
a and ω
a
the corresponding o pen subdomain of diameter h
ω
a
. For K ∈ T
h
, V
K
denotes the set of vertices
of K. From Section
3.2 o nwards, we will need the shape-regularity assumption requesting the existence of
a constant κ
T
> 0 such that max
K∈T
h
h
K
/
K
≤ κ
T
for all tr iangulations T
h
, with
K
being the diameter
of the largest ball inscribed in K. We will also invoke the average operator {{·}} yielding the mean value of
the traces from adjacent mesh elements on inner edges a nd the actual trace on boundary edges; similarly,
the jump operator [[·]] y ields the difference evaluated along n
e
on e ∈ E
int
h
and the actual trace on e ∈ E
ext
h
.
2.3 Broken spaces
At some places, we will use the mesh-related broken Sobolev spaces H
1
(T
h
) := {v ∈ L
2
(Ω); v|
K
∈ H
1
(K)
for all K ∈ T
h
} as well a s H(div, T
h
) := {v ∈ [L
2
(Ω)]
2
; v|
K
∈ H(div, K) for all K ∈ T
h
}. Then, ∇ stands
for the broken (elementwis e) weak gradient, ∇· for the broken (elementwise) weak divergence, and R
π
2
∇ for
the broken (elementwise) weak curl.
2.4 Finite element spaces
We use P
p
(K) (respectively, Q
p
(K)), p ≥ 0, to denote polynomials in K ∈ T
h
of to tal degree at most
p (respectively, at most p in each variable), and P
p
(T
h
) and Q
p
(T
h
) to denote the corresponding broken
spaces. For a vertex a ∈ V
h
, let ψ
a
stand for the “hat” function fro m P
p
(T
h
) ∩ H
1
(Ω) or Q
p
(T
h
) ∩ H
1
(Ω)
which takes value 1 at the vertex a and zero at the other vertices. Following Brezzi and Fortin [
21] or
Roberts and Thomas [
69], let RT
p
:= {v
h
∈ H(div, Ω); v
h
|
K
∈ RT
p
(K)}, p ≥ 0, with the loc al spaces
RT
p
(K) := [P
p
(K)]
2
+ P
p
(K)x on triangles and RT
p
(K) := Q
p+1,p
(K) × Q
p,p+1
(K) on rec tangles, where
Q
·,·
(K) sets the maximal polynomial degree separately for each variable. We will employ this Raviart–
Thomas (RT) family, with P
p
(T
h
) or Q
p
(T
h
) for the corresponding L
2
(Ω) approximations, a nd we use the
abstract notatio n V
h
:= RT
p
, Q
h
:= P
p
(T
h
) or Q
p
(T
h
), V
h
(K) := RT
p
(K), and Q
h
(K) := P
p
(K) or
Q
p
(K); this allows us to discuss other families, like the Brez z i–Douglas–Marini one in Rema rk
3.21.
2.5 The model problem
We study in this paper the Poisson problem for the Laplace equation: for f ∈ L
2
(Ω), find u such that
−∆u = f in Ω, (2.1a)
u = 0 on ∂Ω. (2.1b)
The weak formulation cons ists in finding u ∈ H
1
0
(Ω) such that
(∇u, ∇v) = (f, v) ∀v ∈ H
1
0
(Ω). (2.2)
Existence and uniqueness of the solution u to (
2.2) follow from the Riesz representation theorem. We term
the s c alar-valued function u the potential and the vector-valued function σ := −∇u the flux. Extensions
to inhomogeneous Dirichlet and Neumann boundary conditions, more general meshes, meshes with hanging
nodes, and approximations with varying polynomial degree ar e possible modulo necessary technicalities.
3
3 Main results
We present in this section our main results. The gua ranteed error upper bound is presented in Section
3.1 and
a lower bound robust with respect to the polynomial degree is stated in Section 3.2. Maximal overestimation
is investigated in Section
3.3.
3.1 Guaranteed reliability
Let u
h
denote the given approximate solution to problem (2.2). In this section, we only need u
h
∈ H
1
(T
h
).
3.1.1 Equilibrated fl ux and potential reconstructions
Discrete so lutions are typically such that u
h
6∈ H
1
0
(Ω), −∇u
h
6∈ H(div, Ω), or ∇·(−∇u
h
) 6= f, while the weak
solution satisfies u ∈ H
1
0
(Ω), σ ∈ H(div, Ω), and ∇·σ = f with σ := −∇u. We begin by restoring /mimick ing
these three properties of the weak solution:
Definition 3.1 (Equilibrated flux reconstruction). We call an equilibrated flux reconstruction any function
σ
h
constructed from u
h
which satisfies
σ
h
∈ H(div, Ω), (3.1a)
(∇·σ
h
, 1)
K
= (f, 1)
K
∀K ∈ T
h
. (3.1b)
Definition 3. 2 (Potential r econstruction). We call a potential reco nstruction any fun ction s
h
constructed
from u
h
which satisfies
s
h
∈ H
1
0
(Ω). (3.2)
3.1.2 Guaranteed reliability
The error upper bound is straightforward:
Theorem 3.3 (A guara nteed a poster iori error es tima te). Let u be the weak solution of (
2.2) and let
u
h
∈ H
1
(T
h
) be arbitrary. Let σ
h
and s
h
be respectively the equilibrated flux and potential reconstructions
of D efinitions
3.1 and 3.2. Then
k∇(u − u
h
)k
2
≤
X
K∈T
h
k∇u
h
+ σ
h
k
K
+
h
K
π
kf − ∇·σ
h
k
K
2
+
X
K∈T
h
k∇(u
h
− s
h
)k
2
K
. (3.3)
Proof. The proof is straightforward along [
64, 55, 35, 60, 66, 57, 2, 52, 76, 41, 68, 19, 43, 28, 29]. We ske tch
it for self-completeness. As in [
60, 52], define s ∈ H
1
0
(Ω) by
(∇s, ∇v) = (∇u
h
, ∇v) ∀v ∈ H
1
0
(Ω). (3.4)
Its existence and uniqueness follow from the Riesz representation theorem. From this projection-type
construction results the Pythagore an equality
k∇(u − u
h
)k
2
= k∇(u − s)k
2
+ k∇(s − u
h
)k
2
(3.5)
and the minimization property
k∇(s − u
h
)k
2
= min
v∈ H
1
0
(Ω)
k∇(v − u
h
)k
2
. (3.6)
For the first term in (
3.5), using that u − s ∈ H
1
0
(Ω), (3.4) yields
k∇(u − s)k = sup
v∈ H
1
0
(Ω); k∇vk=1
(∇(u − s), ∇v) = sup
v∈ H
1
0
(Ω); k∇vk=1
(∇(u − u
h
), ∇v).
Let v ∈ H
1
0
(Ω) be fixed. Using (
2.2) and adding and subtra cting (σ
h
, ∇v),
(∇(u − u
h
), ∇v) = (f − ∇·σ
h
, v) − (∇u
h
+ σ
h
, ∇v),
4
where we have also employed the Green theorem (σ
h
, ∇v) = −(∇·σ
h
, v). The C auchy–Schwarz inequality
yields
−(∇u
h
+ σ
h
, ∇v) ≤
X
K∈T
h
k∇u
h
+ σ
h
k
K
k∇vk
K
,
whereas the approximate equilibrium property (
3.1b), the Poincar´e inequality
kw − Π
0
K
wk
K
≤ C
P,K
h
K
k∇wk
K
∀w ∈ H
1
(K), (3.7)
with Π
0
K
w the mean value of w on K and C
P,K
= 1/π thanks to the convexity of the mesh elements K (see
Payne and Weinberger [63] and Bebendorf [17]), and the Cauchy–Schwarz inequality yield
(f − ∇·σ
h
, v) =
X
K∈T
h
(f − ∇·σ
h
, v − Π
0
K
v)
K
≤
X
K∈T
h
h
K
π
kf − ∇·σ
h
k
K
k∇vk
K
. (3.8)
Combining these results with the Cauchy–Schwarz ineq uality and since any s
h
∈ H
1
0
(Ω) bounds (
3.6), we
infer the assertion.
3.1.3 Mixed finite element solution of Neumann problems on patches using the partition of
unity
This sectio n describes a pr actical way to obta in the equilibrated flux and potential reconstructions introduced
in Definitions
3.1 and 3.2. Fo r the flux reconstruction, we rewrite equivalently the technique of [19], see
also [
38], pr oc e eding as in [45]. The potential reconstruction is close to that of [29, Section 6.3]. In both cases,
the equilibration goes over patches of elements ω
a
sharing a generic vertex a ∈ V
h
, with V
h
(ω
a
) × Q
h
(ω
a
)
denoting restrictio ns to ω
a
of the mixed finite element spaces discussed in Section 2.4. We still only assume
u
h
∈ H
1
(T
h
).
Construction 3.4 (Flux σ
h
). Let u
h
satisfy the hat-function orthogonality
(∇u
h
, ∇ψ
a
)
ω
a
= (f, ψ
a
)
ω
a
∀a ∈ V
int
h
. (3.9)
For each a ∈ V
h
, prescribe ς
a
h
∈ V
a
h
and ¯r
a
h
∈ Q
a
h
by solving
(ς
a
h
, v
h
)
ω
a
− (¯r
a
h
, ∇·v
h
)
ω
a
= −(τ
a
h
, v
h
)
ω
a
∀v
h
∈ V
a
h
, (3.10a)
(∇·ς
a
h
, q
h
)
ω
a
= (g
a
, q
h
)
ω
a
∀q
h
∈ Q
a
h
, (3.10b)
with the spaces
V
a
h
:= {v
h
∈ V
h
(ω
a
); v
h
·n
ω
a
= 0 on ∂ω
a
},
Q
a
h
:= {q
h
∈ Q
h
(ω
a
); (q
h
, 1)
ω
a
= 0},
a ∈ V
int
h
, (3.11a)
V
a
h
:= {v
h
∈ V
h
(ω
a
); v
h
·n
ω
a
= 0 on ∂ω
a
\ ∂Ω},
Q
a
h
:= Q
h
(ω
a
),
a ∈ V
ext
h
, (3.11b)
and the right-hand sides
τ
a
h
:= ψ
a
∇u
h
, (3.12a)
g
a
:= ψ
a
f − ∇ψ
a
·∇u
h
. (3.12b)
Then, set
σ
h
:=
X
a∈V
h
ς
a
h
. (3.13)
In (
3.11), a homogeneous Neumann (no-flux) bounda ry condition on the whole boundary of the patch
ω
a
together with mean value zero is imposed for interior vertices, whereas the no-flux condition is only
imposed in the interior of Ω for boundary vertices. Also note that by (
3.9), (g
a
, 1)
ω
a
= 0 for interior vertices
a, which is the Neumann compatibility condition. Existence and uniqueness of the solution to (
3.10) are
standard, see [
21, 69, 78]. We now verify the requir ements of Definition 3.1:
5
