Adaptive Multilevel Techniques for Mixed Finite Element Discretizations of Elliptic Boundary Value Problems.
We consider mixed finite element discretizations of linear second-order elliptic boundary value problems with respect to an adaptively generated hierarchy of possibly highly nonuniform simplicial triangulations. By a well-known postprocessing technique the discrete problem is equivalent to a modified nonconforming discretization which is solved by preconditioned CG iterations using a multilevel preconditioner in the spirit of Bramble, Pasciak, and Xu designed for standard nonconforming approximations. Local refinement of the triangulations is based on an a posteriori error estimator which can be easily derived from superconvergence results. The performance of the preconditioner and the error estimator is illustrated by several numerical examples.
Introduction
.
In this work, we are concerned with adaptive multilevel techniques for the efficient
solution of mixed finite element discretizations of linear second order
elliptic boundary value problems. In recent years, mixed finite element methods
have been increasingly used in applications, in particular for such problems
where instead of the primal variable its gradient is of major interest. As examples
we mention the flux in stationary flow problems or neutron diffusion
and the current in semiconductor device simulation (cf. e.g. [4], [13], [14], [22],
[27], [36], [42] and [44]). An excellent treatment of mixed methods and further
references can be found in the monography of Brezzi and Fortin [12].
Mixed discretization give rise to linear systems associated with saddle point
problems whose characteristic feature is a symmetric but indefinite coefficient
matrix. Since the systems typically become large for discretized partial differential
equations, there is a need for fast iterative solvers. We note that preconditioned
iterative methods for saddle point problems have been considered by
Bank, Welfert and Yserentant [8] based on a modification on Uzawa's method
leading to an outer/inner iterative scheme and by Rusten and Winther [43] relying
on the minimum residual method. Moreover, there are several approaches
using domain decomposition techniques and related multilevel Schwarz iterations
(cf. e.g. Cowsar [15], Ewing and Wang [23, 24, 25], Mathew [32, 33]
and Vassilevski and Wang [46]). A further important aspect is to increase efficiency
by using adaptively generated triangulations. In contrast to the existing
concepts for standard conforming finite element discretizations as realized for
example in the finite element codes PLTMG [5] and KASKADE [19, 20], not
much work has been done concerning local refinement of the triangulations
in mixed discretizations. There is some work by Ewing et al. [21] in case
of quadrilateral mixed elements but the emphasis is more on the appropriate
treatment of the slave nodes then on efficient and reliable indicators for local
refinement. It is the purpose of this paper to develop a fully adaptive algorithm
for mixed discretizations based on the lowest order Raviart-Thomas elements
featuring a multilevel iterative solver and an a posteriori error estimator as
indicator for local refinement. The paper is organized as follows:
In section 2 we will present the mixed discretization and a postprocessing
technique due to Fraeijs de Veubeke [26] and Arnold and Brezzi [1]. This technique
is based on the elimination of the continuity constraints for the normal
components of the flux on the interelement boundaries from the conforming
Raviart-Thomas ansatz space. Instead, the continuity constraints are taken
care of by appropriate Lagrangian multipliers resulting in an extended saddle
point problem. Static condensation of the flux leads to a linear system which
is equivalent to a modified nonconforming approach involving the lowest order
Crouzeix-Raviart elements augmented by cubic bubble functions. Section 3 is
devoted to the numerical solution of that nonconforming discretization by a
multilevel preconditioned cg-iteration using a BPX-type preconditioner. This
preconditioner has been designed by the authors [30, 49] for standard non-conforming
approaches and is closely related to that of Oswald [39]. By an
application of Nepomnyaschikh's fictitious domain lemma [34, 35] it can be
verified that the spectral condition number of the preconditioned stiffness matrix
behaves like O(1). In section 4 we present an a posteriori error estimator
in terms of the L 2 -norm which can be easily derived from a superconvergence
result for mixed discretizations due to Arnold and Brezzi [1]. It will be shown
that the error estimator is equivalent to a weighted sum of the squares of
the jumps of the approximation of the primal variable across the interelement
boundaries. Finally, in section 5 some numerical results are given illustrating
both the performance of the preconditioner and the error estimator.
Mixed discretization and postprocessing.
We consider linear, second order elliptic boundary value problems of the form
\Gammadiv(a
@\Omega
stands for a bounded, polygonal domain in the Euclidean space IR 2
with boundary \Gamma and f is a given function in L
2(\Omega\Gamma4 We further assume that
i;j=1 is a symmetric 2 \Theta 2 matrix-valued function with a ij 2 L
and b is a function in L
satisfying
for almost all x 2 \Omega\Gamma We note that only for simplicity we have chosen homogeneous
Dirichlet boundary conditions in (2.1). Other boundary conditions of
Neumann type or mixed boundary conditions can be treated as well. Introducing
the Hilbert space
ae
oe
and the flux
as an additional unknown, the standard mixed formulation of (2.1) is given as
follows:
Find (j; u) 2 H(div; \Omega\Gamma \Theta L
2(\Omega\Gamma such that
where the bilinear forms a : H(div; \Omega\Gamma \Theta
are given by
R
\Omega
R
\Omega
divq
R
\Omega
bu
and (\Delta; \Delta) 0 stands for the usual L 2 -inner product. Note that under the above
assumption on the data of the problem the existence and uniqueness of a
solution to (2.3) is well established (cf. e.g. [12]). For the mixed discretization
of (2.3) we suppose that a regular simplicial triangulation T h
of\Omega is given. In
particular, for an element K 2 T h we refer to e i its edges and
we denote by E h the set of edges of T h and by E 0
the
subsets of interior and boundary edges, respectively. Further, for D
'\Omega we
refer to jDj as the measure of D and we denote by P k (D), k - 0, the linear
space of polynomials of degree - k on D. Then, a conforming approximation
of the flux space
H(div;\Omega\Gamma is given by V h := RT
RT
and RT 0 (K) stands for the lowest order Raviart-Thomas element
Note that any q h 2 RT 0 (K) is uniquely determined by its normal components
on the edges e i denotes the outer normal
vector of K. In particular, the conformity of the approximation is guaranteed
by specifying the basis in such a way that continuity of the normal components
is satisfied across interelement boundaries. Consequently, we have dimV
the standard mixed discretization of (2.3) is given by:
Find (j h ; u h h such that
For D
'\Omega we denote by (\Delta; \Delta) k;D , k - 0, the standard inner products and
by k \Delta k k;D the associated norms on the Sobolev spaces H k (D) and
respectively. For simplicity, the lower index D will be omitted if
it is well known that assuming u
the following a
priori error estimates hold true
where h stands as usual for the maximum diameter of the elements of T h and
C is a positive constant independent of h, u and j (cf. e.g. [1]; Thm. 1.1).
We further observe that the algebraic formulation of (2.5) gives rise to a linear
system with coefficient matrix
which is symmetric but indefinite. There exist several efficient iterative solvers
for such systems, for example those proposed by Bank et al. [8], Cowsar [15],
Ewing and Wang [23, 24, 25], Mathew [32], Rusten and Winther [43] and
Vassilevski and Wang [46]. However, we will follow an idea suggested by
Fraeijs de Veubeke [26] and further analyzed by Arnold and Brezzi in [1] (cf.
also [12]). Eliminating the continuity constraints (2.4) from V h results in the
nonconforming Raviart-Thomas space -
ae
oe
Since there are now two basic vector fields associated with each e
h , we
have -
h . Instead, the continuity constraints are taken
care of by Lagrangian multipliers living in M h := M 0
and
Then the nonconforming mixed discretization of (2.3) is to find (j h ; u
h \Theta W h \Theta M h such that
are given
by
R
R
R
As shown in [1] the above multiplier technique has two significant advantages.
The first one is some sort of a superconvergence result concerning the approximation
of the solution u in (2.1) in the L 2 -norm while the second one is
related to the specific structure of (2.7) and has an important impact on the
efficiency of the solution process. To begin with the first one we denote by \Pi h
the L 2 -projection onto M h . Then it is easy to see that there exists a unique
(cf. [1] Lemma 2.1). The function -
represents a nonconforming interpolation
of - h which can be shown to provide a more accurate approximation of u in
the L 2 -norm. In particular, if u
1(\Omega\Gamma then there exists a
constant c ? 0 independent of h, u and j such that
(cf. [12] Theorem 3.1, Chap. 5). The preceding result will be used for the
construction of a local a posteriori error estimator to be developed in Section
4.
As far as the efficient solution of (2.7) is concerned we note that the algebraic
formulation leads to a linear system with a coefficient matrix of the formB @
In particular, -
A stands for a block-diagonal matrix, each block being a 3 \Theta 3
matrix corresponding to an element K 2 T h . Hence, -
A is easily invertible
which suggests block elimination of the unknown flux (also known as static
resulting in a 2 \Theta 2 block system with a symmetric, positive definite
coefficient matrix. This linear system is equivalent to a modified nonconforming
approximation involving the lowest order Crouzeix-Raviart elements
augmented by cubic bubble functions. Denoting by m e the midpoint of an
and we set
Note, that dimCR
Further, we denote
by P h and -
P c the L 2 -projections onto W h and -
the latter with respect to
the weighted L 2 -inner product (\Delta; \Delta) As shown in [1], (Lemma 2.3
and Lemma 2.4) there exists a unique \Psi h 2 N h such that
Originally, Lemma 2.4 is only proved for b j 0 but the result can be easily
generalized for functions b - 0. Moreover, \Psi h is the unique solution of the
variational problem
where the bilinear form aN h
is given by
Z
We will solve (2.11) numerically by preconditioned cg-iterations using a multilevel
preconditioner of BPX-type. The construction of that preconditioner
will be dealt with in the following section.
3 Iterative solution by multilevel preconditioned
cg-iterations.
We assume a hierarchy (T k
k=0 of possibly highly nonuniform triangulations
of\Omega obtained by the refinement process due to Bank et al. [6] based on
regular refinements (partition into four congruent subtriangles) and irregular
refinements (bisection). For a detailed description including the refinement
rules we refer to [5] and [17]. We remark that the refinement rules are such
that each K 2 geometrically similar either to an element of
or to an irregular refinement of a triangle in T 0 . Consequently, there exist
constants depending only on the local geometry of T 0 such that
for all K 2 its edges e ae @K
Moreover, the refinement rules imply the property of local quasiuniformity,
i.e., there exists a constant depending only on the local geometry of T 0
such that for all K;K
where hK := diamK.
We consider the modified nonconforming approximation (2.11) on the highest
level j
and we attempt to solve (3.3) by preconditioned cg-iterations. The preconditioner
will be constructed by means of the natural splitting of N j into the
standard nonconforming part CR j := CR h j and the "bubble" part B
and a further multilevel preconditioning of BPX-type for the nonconforming
part. For that purpose we introduce the bilinear form a CR j
a CR j (u CR
aj K (u CR
is the standard bilinear form associated with
the primal variational formulation (2.1)
Z
\Omega
(aru
In the sequel we will refer to A : H 1
0(\Omega\Gamma as the operator associated
with the bilinear form a.
Further, we define the bilinear form a
(w B
Z
a -
for all w B
. Denoting by AD j , D g, the operators
associated with aD j
, we will prove the spectral equivalence of AN j and ACR
. To this end we need the following technical lemmas:
Lemma 3.1 For all u CR
there holds
Proof. For the reference triangle -
K with vertices (0; 0), (1; 0) and (0; 1) it is
easy to establish
(3.7) can be deduced by the affine equivalence of the Crouzeix-Raviart elements

Lemma 3.2 For all w B
there holds
Proof. Since
are the barycentric
coordinates of K, we have
Denoting by - i
the local basis of -
h and by ( -
the
matrix representation of -aj K in case
stands for the vertex opposite to e i , by Green's
3: (3.11)
If we consider the reference triangle -
K where the vertices are given by (0; 0),
refers to the usual partial order on the set of symmetric, positive
definite matrices. Moreover, taking advantage of the affine equivalence of the
Raviart-Thomas elements it is easy to show
Using (3.1), (3.11) and (3.12) in (3.10) and observing (3.9) it follows that
We assume a and b to be locally constant, i.e., a ij
and we denote by ff 0;K and ff 1;K the lower and
upper bounds arising in (2.2) when restricting a to K. We further suppose
that a and b are such that
min
- 0: (3.13)
Note, that only for simplicity we have chosen the strong inequality (3.13).
All results can be extended to the more general case that a constant c ? 0,
independent of K exists such that for all K 2 ch 2
holds.
Under the assumption (3.13) there holds:
Theorem 3.3 Under the assumption (3.13) there exist constants
depending on the local bounds ff l;K , l 2 f0; 1g, K 2 T j , such that for all
with
a CR j
(w B
a CR j (u CR
Proof. For the proof of the preceding result we use the following lemma which
can easily established.
Lemma 3.4 For all there hold
ff 0;K
(aru CR
(3.15 a)
Proof. Using the Cauchy-Schwarz inequality we obtain
j as well as the orthogonality (ru CR
we obtain
ff 0;K
(aru CR
The following inequality deduces (3.15 b)
On the other hand, in view of kP h j u CR
0;k we have
(bu CR
Combining (3.15 a) and (3.16 a) gives the upper bound in (3.14) with c
ff 1;K
ff 0;K
. Further, by Young's inequality,
(bu CR
(bu CR
Consequently, using (3.15 b), (3.16 b), (3.13) and
ff 1;K
(aru CR
(bu CR
which yields the lower bound in (3.14) with c
ff 0;K
ff 1;K
We note that the bilinear form a B j gives rise to a diagonal matrix which thus
can be easily used in the preconditioning process. On the other hand, the bi-linear
form a CR j
corresponds to the standard nonconforming approximation of
(2.1) by the lowest order Crouzeix-Raviart elements. Multilevel preconditioner
for such nonconforming finite element discretizations have been developed by
Oswald [39, 40], Zhang [53] and the authors [30, 49]. Here we will use a BPX-
type preconditioner based on the use of a pseudo-interpolant which allows to
identify with a closed subspace of the standard conforming ansatz space
with respect to the next higher level. More precisely, we denote by T j+1 , the
triangulation obtained from T h by regular refinement of all elements
and we refer to S k ae H 1
as the standard conforming ansatz
space generated by continuous, piecewise linear finite elements with respect to
the triangulation T k . Denoting by N 0
j+1 the set of interior vertices of T j+1 and
recalling that the midpoints m e of interior edges
j correspond to vertices
j+1 , we define a mapping P CR
are the midpoints of those interior edges
having
j+1 as a common vertex. We note that this pseudo-interpolant has been
originally proposed by Cowsar [15] in the framework of related domain decomposition
techniques. The following result will lay the basis for the construction
of the multilevel preconditioner:
Lemma
j be the pseudo-interpolant given by (3.17). Then there
exist constants depending only on the local geometry of T 0 such
that for all
Proof. The assertion follows by arguing literally in the same way as in [15]
(Theorem 2) and taking advantage of the local quasiuniformity of the triangulations

It follows from (3.18) that ~
represents a closed subspace of
being isomorphic to CR j . Based on this observation we may now use the
well known BPX-preconditioner for conforming discretizations with respect to
the hierarchy (S k ) j+1
k=0 of finite element spaces (cf. e.g. [10], [11], [16], [41], [50],
[52], and [53]). We remark that for a nonvanishing Helmholtz term in (2.1) the
initial triangulation T 0 should be chosen in such a way that the magnitude of
the coefficients of the principal part of the elliptic operator is not dominated
by the magnitude of the Helmholtz coefficient times the square of the maximal
diameter of the elements in T 0 (cf. e.g. [37], [51]).
Denoting by \Gamma k := fOE (k)
, the set of nodal basis functions
of the BPX-preconditioner is based on the following
structuring of the nodal bases of varying index k:
We introduce the Hilbert space
Y
equipped with the inner product
where
, and we consider the bilinear form
denoting by ~
the operator associated with ~ b. We further define a
mapping
and refer to R
V as its adjoint in the sense that (R V u; v)
. Then the BPX-preconditioner is given by
satisfying
with constants depending only on the local geometry of T 0 and
on the bounds for the data a, b in (2.2).
The condition number estimates (3.23) have been established by various authors
(cf. [10], [16], [38]). They can be derived using the powerful Dryja-
theory [18] of additive Schwarz iterations. Another approach due to
Oswald [41] is based on Nepomnyaschikh's fictitious domain lemma:
Lemma 3.6 Let S and V be two Hilbert spaces with inner products (\Delta; \Delta) S and
generated by symmetric, positive definite operators A S : S 7\Gamma! S and ~
. Assume that there exist a linear operator
necessarily linear) operator
R
a S (Rv; Rv) - c 1
Then there holds
c 0 a S (u; u) - a S (R ~
is the adjoint of R in the sense that (Rv; u)
Proof. See e.g. [35].
In the framework of BPX-preconditioner with a S being the bilinear
form in (3.5) while V , ~ b and R are given by (3.19), (3.20) and
(3.21), respectively. The estimate (3.24 b) is usually established by means of a
strengthened Cauchy-Schwarz inequality. Further, is an appropriately
chosen decomposition operator such that the P.L. Lions type estimate (3.24 c)
holds true (cf. e.g. [41] Chapter 4).
Now, returning to the nonconforming approximation we define I S
CR j by
I S
. Note that in view of (3.17)
the operators I S
j corresponds to the identity on CR j . Then, with C as
in (3.22) the operator
is an appropriate BPX-preconditioner for the nonconforming discretization of
(2.1). In particular, we have:
Theorem 3.7 Let CNC be given by (3.25). Then there exist positive constants
depending only on the local geometry of T 0 and on the bounds for the
coefficients a, b in (2.2) such that for all u 2 CR j
Proof. In view of the fictitious domain lemma we choose
as in (3.4) and V and ~ b according to (3.19), (3.20). Furthermore, we specify
with T S as the decomposition operator in the conforming setting.
Obviously
Moreover, using the obvious inequality
a CR j
and (3.24 b), for all v 2 V we have
a CR j
(Rv;
Finally, using again (3.18) and (3.24 c) for
get
~
1 a CR j (P CR
(3.27 c)
In terms of (3.27 a-c) we have verified the hypotheses of the fictitious domain
lemma which gives the assertion.
4 A posteriori error estimation.
Efficient and reliable error estimators for the total error providing indicators
for local refinement of the triangulations are an indispensable tool for efficient
adaptive algorithms. Concerning the finite element solution of elliptic
boundary value problems we mention the pioneering work done by Babuska
and Rheinboldt [2, 3] which has been extended among others by Bank and
Weiser [7] and Deuflhard, Leinen Yserentant [17] to derive element-oriented
and edge-oriented local error estimators for standard conforming approxima-
tions. We remark that these concepts have been adapted to nonconforming
discretizations by the authors in [29, 30] and [49]. The basic idea is to discretize
the defect problem for the available approximation with respect to a
finite element space of higher accuracy. For a detailed representation of the
different concepts and further references we refer to the monographs of Johnson
[31], Szabo and Babuska [45] and Zienkiewicz and Taylor [54] (cf. also the
recent survey articles by Bornemann et al. [9] and Verf?rth [47, 48]).
In this section we will derive an error estimator for the L 2 -norm of the total
error in the primal variable u based on the superconvergence result (2.9). As
we shall see this estimator does not require the solution of an additional defect
problem and hence is much more cheaper than the estimators mentioned
above. We note, however, that an error estimator for the total error in the
flux based on the solution of localized defect problems has been developed by
the first author in [28].
We suppose that ~
is an approximation of the solution / h 2 N h
of (2.11) obtained, for example, by the multilevel iterative solution process
described in the preceding section. Then, in view of (2.7) and (2.10) we
get an approximation ( ~
h \Theta W h \Theta M h of the unique solution
(j
h \Theta W h \Theta M h of (2.7) by means of
~
and
The last equality is obtained due to:
R
R
R
R
Further, we denote by - ~ u h 2 CR h the nonconforming extension of ~ - h .
In lights of the superconvergence result (2.9) we assume the existence of a
constant
In other words, (4.3) states that the nonconforming extension - u h of - h does
provide a better approximation of the primal variable u than the piecewise
constant approximation u h .
It is easy to see that (4.3) yields
(4.
Observing (2.10) and (4.1), we have
K2T h3
K2T h3
ii
q3
Using (4.5), (4.6) in (4.4) we get
q3
We note that k/ h \Gamma ~
represents the L 2 - norm of the iteration error whose
actual size can be controlled by the iterative solution process. Therefore, the
provides an efficient and reliable error estimator for the L 2 -
norm of the total error whose local contributions k~u
be used as indicators for local refinement of T h . Moreover, the estimator can
be cheaply computed, since it only requires the evaluation of the available
approximations ~ u h 2 W h and ~
For a better understanding of the estimator the rest of this section will be
devoted to show that it is equivalent to a weighted sum of the squares of the
jumps of ~
across the edges e 2 E h . For that purpose we introduce the jump
and the average of piecewise continuous functions v along edges In
particular, for e
h we denote by K in and K out the two adjacent triangles
and by n e the unit normal outward from K in . On the other hand, for
we refer to n e as the usual outward normal. Then, we define the average [v] A
of v on e 2 E h and the jump [v] J of v on e 2 E h according to
(vj K in
It is easy to see that for piecewise continuous functions u, v there holds
Z
Z
e
(uj K in \Delta vj K in
Z
e
(uj K in \Delta vj K in \Gamma uj Kout \Delta vj Kout
Z
e
Further, we observe that for vector fields q the quantity [n \Delta q] J is independent
of the choice of K in and K out .
In terms of the averages [n e \Delta q h ] A and the jumps [n e \Delta q h ] J we may decompose
the nonconforming Raviart-Thomas space -
into the sum
where the subspaces -
h and -
H are given by
Obviously, we have -
g. As the main result of this section we will prove:
Theorem 4.1 Let ( ~ j h
h \Theta W h \Theta M h be an approximation of the
unique solution of (2.6) obtained according to (4.1), (4.2) and let - ~
be the nonconforming extension of ~
- h . Then there exist constants
depending only on the shape regularity of T h and the ellipticity constants in
(2.2) such that
The proof of the preceding result will be provided in several steps. Firstly, due
to the shape regularity of T h we have:
Lemma 4.2 Under the assumptions of Theorem 4.1 there holds3
'i
'i
Proof. By straightforward computation
~
jK in j
~
~
which easily gives (4.12) by taking advantage of (3.1).
As a direct consequence of Lemma 4.1 we obtain the lower bound in (4.11)
with oe
However, the proof of the upper bound is more elaborate. In
view of (4.12) it is sufficient to show that
ii
holds true with an appropriate positive constant c. As a first step in this direction
we will establish the following relationship between ~ - h and the averages
and jumps of ~
Lemma 4.3 Under the assumptions of Theorem 4.1 for all q h
V h there holds
ii
where P c denotes the projection onto V h with respect to the weighted L 2 -inner
product (\Delta; \Delta) 0;c .
Proof. We denote by ~
OE h the unique element in B h satisfying
Z
~
Z
~
In view of (4.2) we thus have
K2T h@
Z
Z
By Green's formula, observing ~
Z
Z
and hence
Z
\Omega
c
ar ~
which shows that
Consequently, for q h
Z
Z
cP c (ar ~
Z
r ~
Z
~
Z
~
Z
~
It follows from (4.2) that
Z
~
which by (4.8 b) is clearly equivalent to the assertion.
For a particular choice of q h 2 -
V h in Lemma 4.2 we obtain an explicit representation
of ~
- h on e
h . We choose q
(- K in
and - K in
e are the standard basis vector fields in -
h with support in K in
resp. K out given by
Corollary 4.4 Let the assumptions of Lemma 4.2 be satisfied and let - e 2 -
h , be given by (4.15). Then there holds
~
Proof. Observing [n e 0
the assertion is a direct consequence of (4.14).
Moreover, with regard to (4.13) we get:
Corollary 4.5 Under the assumption of Lemma 4.2 there holds@ X
ii
Proof. Since for each - h 2 M h (E h ) there exists a unique q h 2 -
A
satisfying
by means of (4.14) we get
ii
A
A
which gives (4.17) by the Schwarz inequality.
The preceding result tells us that for the proof of (4.13) we have to verify
A
Since (4.18) obviously holds true for q h 2 -
, it is sufficient to show:
Lemma 4.6 Let the assumptions of Lemma 4.2 satisfied. Then there holds
A
Proof. We refer to A, -
A and P c as the matrix representations of the operators
A and P c . With respect to the standard bases of V h and -
we may identify
vectors q
respectively. We remark that q h 2 -
A
iff q K in
h , and q K
h . It follows that for q h
A
(q K in
Obviously
stands for the spectral radius of P c \Delta P T
c . Denoting by S the
natural embedding of V h into -
h and by S its matrix representation, it is easy
to see that
A
whence
ASA
ASA
We further refer to AK and -
AK as the local stiffness matrices. Using (2.2)
and (3.12), we get
AKg
with
. Consequently, introducing the local vectors
it follows that
(jK in
e
(jK in
Using (4.24), (4.25) in (4.23) we find
which gives (4.19) in view of (4.20), (4.21) and (4.22).
Summarizing the preceding results it follows that the upper estimate in (4.11)
holds true with oe
s
5 Numerical results.
In this section, we will present the numerical results obtained by the application
of the adaptive multilevel algorithm to some selected second order elliptic
boundary problems. In particular, we will illustrate the refinement process
as well as the performance of both the multilevel preconditioner and the a
posteriori error estimator. The following model problems from [5] and [20]
have been chosen as test examples:
Problem 1. Equation (2.1) with on the unit square
with the right-hand side f and the Dirichlet boundary conditions
according to the solution u(x;
which has a boundary layer along the lines Fig. 5.1).
Problem 2. Equation (2.1) with the right-hand side f j 0 and a hexagon
\Omega with corners (\Sigma1; 0),
. The coefficients are chosen according
to b j 0 and a(x; y) being piecewise constant with the values 1 and 100
on alternate triangles of the initial triangulation (cf. Fig. 5.2). The solution
given by u(x; continuous with a jump discontinuity of
the first derivatives at the interfaces.
Starting from the initial coarse triangulations depicted in Figures 5.1 and 5.2,
on each refinement level l the discretized problems are solved by preconditioned
cg-iterations with a BPX-type preconditioner as described in Section 3. The
iteration on level l stopped when the estimated iteration error " l+1 is
l
l
, with the safety factor l denotes the
estimated error on level l, the number of nodes on level l and l +1 are given by
N l and N l+1 , respectively. Denoting by ( ~ j l ; ~
resulting approximation
and by -
~ u l the nonconforming extension of ~ - l , for the local refinement of T l the
error contributions ffl 2
~
and the weighted
mean value
K2T l
K are computed. Then, an element K 2 T l is
marked for refinement if jKj
oe is a safety factor which is
chosen as 0:95. The interpolated values of the level l approximation are
used as startiterates on the next refinement level. For the global refinement
process we use -ffl 2
0;\Omega as stopping criteria, where ff is a safety
factor which is chosen as and tol is the required accuracy,
Level 0,

Figure

5.1: Initial triangulation T 0 and final triangulation T 6 (Problem 1)
Level 0,

Figure

5.2: Initial triangulation T 0 and final triangulation T 5 (Problem 2)

Figures

5.1 and 5.2 represent the initial triangulations T 0 and the final triangulations
T 6 and T 5 for Problems 1 and 2, respectively. For Problem 1 we observe
a pronounced refinement in the boundary layer (cf. Fig. 5.1). For Problem 2
there is a significant refinement in the areas where the diffusion coefficient is
large with a sharp resolution of the interfaces between the areas of large and
small diffusion coefficient (cf. Fig. 5.2).1.11.31.51.710 100 1000 10000
estimated
error/true
error
Number of nodes
Boundary layer
Discontinuous coefficients

Figure

5.3: Error Estimation for Problem 1 and 25152535450 10000 20000 30000 40000 50000
Number
of
cg-iterations
Number of nodes
Boundary layer
Discontinuous coefficients

Figure

5.4: Preconditioner for Problem 1 and 2
The behaviour of the a posteriori L 2 -error estimator is illustrated in Figure
5.3 where the ratio of the estimated error and the true error is shown as a
function of the total number of nodes. The straight and the dashed lines refer
to Problem 1 (boundary layer) and Problem 2 (discontinuous coefficients),
respectively. In both cases we observe a slight overestimation at the very beginning
of the refinement process, but the estimated error rapidly approaches
the true error with increasing refinement level.
Finally, the performance of the preconditioner is depicted in Figure 5.4 displaying
the number of preconditioned cg-iterations as a function of the total
number of nodal points. Note that for an adequate representation of the performance
we use zero as initial iterates on each refinement level and iterate
until the relative iteration error is less than In both cases, we observe
an increase of the number of iterations at the beginning of the refinement
process until we get into the asymptotic regime where the numerical results
confirm the theoretically predicted O(1) behaviour.



--R

Mixed and nonconforming finite element methods: im- plementation
estimates for adaptive finite element compu- tations
A posteriori error estimates for the finite element method.


Refinement algorithm and data structures for regular local mesh refinement.
Some a posteriori error estimators for elliptic partial differential equations.
A class of iterative methods for solving saddle point problems.

A basic norm equivalence for the theory of multilevel methods.
Parallel multilevel preconditioners.
Mixed and Hybrid Finite Element Methods.
Two dimensional exponential fitting and application to drift-diffusion models
Domain decomposition and mixed finite elements for the neutron diffusion equation.
Domain decomposition methods for nonconforming finite element spaces of Lagrange-type

Concepts of an adaptive hierarchical finite element code.
Towards a unified theory of domain decomposition alogrithms for elliptic problems.
Version 2.0.

Local refinement via domain decomposition techniques for mixed finite element methods with rectangular Raviart-Thomas elements
Convergence analysis of an approximation of miscible displacement in porous media by mixed finite elements and a modified method of characteristics.
The Schwarz algorithm and multilevel decomposition iterative techniques for mixed finite element methods.
Analysis of multilevel decomposition iterative methods for mixed finite element methods.
Analysis of the Schwarz algorithm for mixed finite element methods.

Mixed finite element discretization of continuity equations arising in semiconductor device simulation.
Adaptive mixed finite elements methods using flux- based a posteriori error estimators

Adaptive multilevel iterative techniques for nonconforming finite elements discretizations.
Numerical Solutions of Partial Differential Equations by the Finite Element Method.
Schwarz alternating and iterative refinement methods for mixed formulations of elliptic problems.
Schwarz alternating and iterative refinement methods for mixed formulations of elliptic problems.

Decomposition and fictitious domain methods for elliptic boundary value problems.
Some aspects of mixed finite elements methods for semiconductor simulation.
Two remarks on multilevel preconditioners.
On a hierarchical basis multilevel method with nonconforming P1 elements.
On discrete norm estimates related to multilevel preconditioner in the finite elements methods.
On a BPX-preconditioner for P1 elements
Multilevel finite element approximation: Theory and Application.
Multigrid applied to mixed finite elements schemes for current continuity equations.
A preconditioned iterative method for saddle point problems.
Mixed finite elements methods for flow through unstructed porous media.
Finite Element Analysis.


Multilevel approaches to nonconforming finite elements discretizations of linear second order elliptic boundary value problems.
Iterative methods by space decomposition and subspace correction.
Hierarchical bases in the numerical solution of parabolic problems.
Old and new convergence proofs for multigrid methods.

The Finite Element Method
--TR
mixed finite elements;multilevel preconditioned CG iterations;a posteriori error estimator