Some multilevel methods on graded meshes.
We consider Yserentant's hierarchical basis method and multilevel diagonal scaling method on a class of refined meshes used in the numerical approximation of boundary value problems on polygonal domains in the presence of singularities. We show, as in the uniform case, that the stiffness matrix of the first method has a condition number bounded by (ln(1/h))2, where h is the meshsize of the triangulation. For the second method, we show that the condition number of the iteration operator is bounded by ln(1/h), which is worse than in the uniform case but better than the hierarchical basis method. As usual, we deduce that the condition number of the BPX iteration operator is bounded by ln(1/h). Finally, graded meshes fulfilling the general conditions are presented and numerical tests are given which confirm the theoretical bounds.
Introduction
The solution of boundary value problems (b.v.p.) in non-smooth domains presents singularities
in the neighbourhood of singular points of the boundary, e.g. in the neighbourhood
of re-entrant corners. Consequently, the use of uniform finite element meshes yields a poor
rate of convergence. Many authors proposed to build graded meshes in the neighbourhood
of these singular points in order to restore the optimal convergence order (see, e.g.
[13, 16]). Roughly speaking, such meshes consist in moving the nodal points by some
coordinate transformation in order to compensate the singular behaviour of the solution,
i.e. that the nodes accumulate near the singular point.
As usual the finite element discretization leads to the resolution of large-scale systems
of linear algebraic equations, where the system matrices in the nodal basis have a large
condition number. This implies that the resolution by iterative methods requires a large
number of iterations. Using preconditioners based on multilevel techniques one can reduce
this number of iterations drastically. The first obstacle is that the graded meshes proposed
in [13, 16] are actually not nested. Consequently, we propose here to build a sequence of
nested graded meshes T two-dimensional domains which are also appropriate
for the approximation of singularities. A similar algorithm was proposed in [12].
For uniform meshes standard multilevel methods, e.g. the hierarchical basis method
[20] and BPX-like preconditioners [3, 4, 5, 8, 10, 14, 15, 19, 21] allow the reduction of the
condition number to the order O((ln h respectively, for two-dimensional
problems.
Similar results were obtained in the case of nonuniformly refined meshes (see, e.g.,
[4, 5, 8, 15, 19, 20]). But these meshes are different from the above graded meshes.
Therefore, our goal is to extend this kind of results to our new meshes. The main idea is
to prove that our graded meshes satisfy the conditions
with positive constants - 1 , - 2 , fi, and fl; hK k and hK l are the exterior diameter of the
triangles l. Using this property, we can prove
that the condition number of the stiffness matrix in the hierarchical basis is of the order
and that the condition number of a (j 1)-level additive Schwarz operator
with multilevel diagonal scaling (MDS method) is of the order O(ln h
The outline of the paper is the following one: In Section 2, we present our model
problem and describe its finite element discretization. In Section 3, we analyse the condition
number of the stiffness matrix in the hierarchical basis by showing the equivalence
between the H 1 -norm and the standard discrete one, and in Section 4, we derive estimates
of the condition number of the MDS method by adapting Zhang's arguments [21]. Section
5 is devoted to the building of the nested graded meshes. We also check that these
meshes are regular and fulfil the conditions (1). Finally, numerical tests are presented in
Section 6 which confirm our theoretical estimates.
2 The model problem
be a bounded domain of the plane with a polygonal boundary \Gamma (i.e. the
union of a finite number of linear segments). On \Omega\Gamma we shall consider usual Sobolev spaces
of norm and semi-norm denoted by k \Delta k
s;\Omega , respectively (we
refer to [11] for more details). As usual, a
s(\Omega\Gamma is the closure in H
of C 1
0(\Omega\Gamma7 the
space of C 1 functions with compact support in \Omega\Gamma
Consider the boundary value problem
whose variational formulation is: Find u 2
a
1(\Omega\Gamma such that
a
where we have set
Z
\Omega
r T urv dx and
Z
\Omega
(\Omega\Gamma3 It is well known that
if\Omega is convex then
2(\Omega\Gamma and consequently
the use of uniform meshes in standard finite element methods yields an optimal order of
convergence h. On the contrary,
if\Omega is not convex then u 62 H
in general and uniform
meshes yield a poor rate of convergence. Many authors [13, 16, 18] have shown that local
mesh grading allows to restore the optimal order. But such meshes are not uniform in the
sense used in standard multilevel techniques. Hereabove and later on, by uniform meshes
we mean either regular refinements (partition of triangles of level k into four congruent
subtriangles of level nonuniformly refinements (PLTMG package of [2]), see for
instance Section 4 of [15] and the references cited there. For this reason, as in [20, 21],
we relax the conditions of the meshes in the following way (graded meshes that fulfil
these conditions are built in Section 5). We suppose that we have a sequence of nested
triangulations fT k g k2IN such that any triangle of T k is divided into four triangles of T k+1 .
We assume that the triangulations are regular in Ciarlet's sense [6], i.e., the ratios hK =ae K
between the exterior diameters hK and the interior diameters ae K of elements
are uniformly bounded from above and the maximal mesh size h tends to
zero as k goes to infinity. We further assume (see Section 3 of [20] and Section 2 of [21])
that there exist positive constants fi; positive constants such that for all
l with K k ae K l , we have
For regular refinements we have We shall see later on
that our graded meshes satisfy (4) with is the
grading parameter.
In each triangulation T k , we use the approximation space
a
is the set of polynomials of order - 1 on K. We consider the Galerkin
approximation solution of
Let us remark that with the mesh T k built in Section 5 and an appropriate parameter -,
we have the error estimate
1;\Omega . 2 \Gammak kfk
where here and in the sequel a . b means that there exists a positive constant C independent
of k and of the above constants fi; fl such that a - Cb. In Section 5, the constant
will also be independent of the grading parameter -.
3 Yserentant's hierarchical basis method
The goal of this section is to show that the stiffness matrix of the Galerkin method in the
hierarchical basis on meshes T k of the previous section has a condition number bounded
by (ln( 1
as in the uniform case. The same result was already underlined by Yserentant
in Section 3 of [20] for nonuniformly refined meshes (in the above sense) by introducing
the condition (4) and by showing that the results for uniformly refined meshes proved in
Section 2 of [20] could be adapted to this kind of meshes satisfying (4). We then follow
the arguments of Section 2 of [20], underline the differences with the standard refinement
rule and also give the dependence with respect to the parameters fi; fl.
Let N k be the set of vertices of the triangles of T k and S k be the space of continuous
functions on -
\Omega and linear on the triangles of T k . For a continuous function u in -
be the function in S k interpolating u at the nodes of T k , i.e.,
I k u
For further use, let us also denote by V k the subspace of S k of functions vanishing at the
nodes of level k \Gamma 1, in other words, V k is the range of I k \Gamma I k\Gamma1 .
On the finite element space S j , define the semi-norm j \Delta j as follows:
The proof of the equivalence of norms we have in mind is based on the two following
preliminary lemmas. The first one concerns equivalence of semi-norms (cf. Lemma 2.4 of
[20]).
Lemma 3.1 For all u
1;\Omega . juj
Proof: In view of Lemma 2.4 of [20], we simply need to show that the following estimates
1;K .
for all To prove this estimate, we remark that K 2 T k\Gamma1 is
divided into four triangles K l 4, such that v is linear in each K l and
are the vertices of K (see

Figure

1). Due to the fact that the triangulation T k is regular, by an affine coordinate
transformation (reducing to the reference element -
we prove that
1;K l .
where I(K l ) is the set of vertices of K l which are not vertex of K. Summing these
equivalences on l = 4, we obtain (9).
The second ingredient is a Cauchy-Schwarz type inequality already proved in Lemma
2.7 of [20] in the case of regularly refined meshes and that we easily extend to the case of
our mesh as suggested in Section 3 of [20].

Figure

1: Triangle divided in four subtriangles K l
Lemma 3.2 For all l , we have
a(u; v) . fl jk\Gammalj
Proof: Similar to the proof of Lemma 2.7 of [20] with the following slight modification:
if K is a fixed triangle of T l and S the boundary strip of K consisting of all triangles of
are subsets of K and meet the boundary of K then due to (4), we
have
meas (K) .
In view to the proof of Lemma 2.7 of [20], this yields the assertion.
Now we can formulate the equivalence between the H 1 norm and the discrete one (see
Theorem 2.2 of [20]).
Theorem 3.3 For all u
1;\Omega .
Proof: For the lower bound, we remark that the assumption (4) and Lemmas 2.2 and
2.3 of [20] imply that
for every K 2 Summing these inequalities on all K 2 T k , we get
1;\Omega . (1
0;\Omega . (1
Therefore by Lemma 3.1 and the triangular inequality, we get
1;\Omega . kI 0 uk 2
By the estimates (12) and (13), we then obtain the lower bound in (11).
Let us now pass to the upper bound. First, Lemma 3.2 and the arguments of Lemma
2.8 of [20] yield
1;\Omega .
On the other hand, the assumption (4), the fact that our triangulation is regular and the
arguments of Lemma 2.9 of [20] lead to
0;\Omega . kI 0 uk 2
The sum of the two above estimates gives the upper bound in (11).
Using a hierachical basis of V j and the former results, we directly get the
Corollary 3.4 The Galerkin stiffness matrix A j of the approximated problem (5) in the
hierarchical basis has a spectral condition number which grows at most quadratically
with the number of levels j, more precisely
4 Multilevel diagonal scaling method
In this section, we analyse the multilevel diagonal scaling method and the BPX algorithm
in the spirit of [21]. Here the main difficulty relies on the fact that our meshes are not
quasi-uniform (quasi-uniform meshes means that hK - h k , for all triangles K 2 T k , for
leading to the fact that the assumption 2.1.c of [21] is violated.
Let us recall that the multilevel diagonal scaling method consists in the following
algorithm: First we represent as a sum
i is the nodal basis function of V k associated with the
interior vertex
being the number of interior vertices of T k . Define
the operator A from V j to V j by
where (\Delta; \Delta) means the L 2
inner product. Let us further define the preconditioner B
MDS
and the j 1-level multilevel diagonal scaling operator PMDS by
The multilevel diagonal scaling algorithm consists in finding u of the Galerkin
problem (5) by solving iteratively (using for instance the conjugate gradient method) the
equation
As usual to solve iteratively (16), the crucial point is to estimate the condition number
of the iteration operator PMDS . For quasi-uniform meshes, it was shown by X. Zhang in
Theorem 3.1 and Section 4 of [21] that this condition number is uniformly bounded (with
respect to the level j). The same result was extended to the case of nonuniformly refined
meshes [8, x5], [15, x4.2.2]. Our goal is to extend this type of results to meshes satisfying
only (4) (actually only the upper bound is sufficient) which can be non quasi-uniform.
Analysing carefully the proof of Theorem 3.1 of [21] we remark that the upper bound is
valid under the assumption (4) (only the upper bound) and is fully independent of the
quasi-uniformity of the meshes. On the contrary the proof of the lower bound uses this
last property. The key point in our proof of this lower bound is the use of Scott-Zhang's
interpolation operator that we recall now for convenience [17]. For a fixed k 2 f0;
with each g, we associate the macro-element
which is actually the support of OE k
i . For any triangle K 2 T k , let us further denote by
S(K) the union of all macro-elements containing K, i.e.,
The following well known facts result from the regularity of the family fT k g k2IN : There
exists a positive integer M (independent of k) such that
card
hK . hK 0 ; for any K; K
A direct consequence of these two properties is that the diameter of S(K) is equivalent
to hK , indeed from the triangular inequality we have
diam
Using the properties (18) and (17), we get
With any nodal point p k
i , we associate one edge oe k
i of one triangle K 2 T k such that
. We now fix a dual basis f/ k
i g of the nodal one fOE k
in the sense that
Z
Then for all v 2
a
Scott-Zhang's interpolation operator - k v on T k is defined by
Z
Note that the operator - k is actually linear continuous from a
a projection on V k (i.e. - k that it enjoys the following local
interpolation property (see Section 4 of [17]): for all triangles K 2 T k and or 1, we
a
Let us notice that Cl'ement's interpolation operator [7, 9] also satisfies (20) but unfortunately
is not a projection on V k .
Now we are able to prove the estimate of the condition number -(PMDS ) of the iteration
operator PMDS .
Theorem 4.1 The multilevel diagonal scaling operator PMDS
Consequently we have
which means that -(PMDS ) grows at most linearly with the number of levels j + 1.
Proof: As already mentioned, the upper bound was proved by X. Zhang in Lemmas 3.2
to 3.5 in [21]. To prove the lower bound instead of using the H 1 -projection on V k , for
which has a global approximation property which is not convenient for
non quasi-uniform meshes, we take advantage of the local interpolation property (20) of
Scott-Zhang's interpolation operator. Indeed for any u 2
a
1(\Omega\Gamma8 we set
with the convention - \Gamma1 Consequently any u 2 V j may be written
Then for all triangles K 2 T k and or 1, we have:
where M(K) is the unique triangle in T k\Gamma1 containing K if
Owing to (20) and (18), we deduce that
Now we decompose u k in the nodal basis, in other words we write
where
. Consequently we get
. ju k (p k
This last estimate being obtained using the equivalence of norms in finite dimensional
spaces on the reference element -
K and an affine coordinate transformation. Using now
the estimate (24) we arrive at
1;\Omega .
Summing this last estimate on and using the property (17), we obtain
1;\Omega .
1;K . juj 2
The sum on
1;\Omega . (j
With the help of Lemma 3.1 of [21] (see also Remark 3.1 of [21]) and the definition of the
bilinear form a, we conclude
The lower bound directly follows.
Let us finish this section by looking at the BPX algorithm. As the BPX preconditioner
is defined by
the BPX operator
equivalent to 1 (uniformly with respect to k), the condition numbers of
PBPX and PMDS are equivalent. This means that the following holds.
Corollary 4.2 The BPX operator enjoys the property
5 Graded nested meshes
The triangulations T k
of\Omega are graded according to Raugel's procedure [11, 16]. But here
since we need a nested sequence of triangulations this procedure is slightly modified. As
a consequence we need to check the regularity of the meshes. In a second step we shall
show that this family satisfies the condition (4).
Let us first describe the construction of the meshes:
Divide\Omega into a coarse triangular mesh T 0 such that each triangle has either one or
no singular point
(of\Omega\Gamma as vertex. If a triangle has a singular point as vertex (i.e. the
interior angle at this point is ? -), it is called a singular triangle and we suppose that all
its angles are acute and the edges hitting the singular point have the same length (this is
always possible by eventual subdivisions).
ii) Any non singular triangle T of T 0 is divided using the regular refinement procedure,
i.e., divide any triangle of T k included in T into four congruent subtriangles of T k+1 , see

Figure

2.

Figure

2: Triangle K 2 T k divided into four congruent subtriangles
iii) Any singular triangle T of T 0 is refined iteratively as follows: Fix a grading parameter
(that for simplicity we take identical for all singular triangles; if there exists more
than one singular point, then we simply need to take the same parameter for triangles
containing the same singular point). In order to make understandable our procedure we
describe explain how to pass from T " T k to T " T k+1 . For
convenience we first recall Raugel's grading procedure.
Introduce barycentric coordinates - in T such that the singular point of T has
the coordinate - 1. For all n 2 IN   , define vertices
in T whose
coordinates are
Raugel's grading procedure consists in defining as the set of triangles described
by their three vertices as follows:
First simply defined by Raugel's procedure, i.e., it is the set of four triangles
described by (26) with Figure 3).
Secondly, the triangulation T " T 2 is built as follows (see Figure 4): The part below
the line - 1 is identical with Raugel's one, namely it is described by the four
triangles of vertices:
On the contrary the part above the line - 1 is modified in order to guarantee
the nested property. More precisely, the set of triangles in this zone is described by
~
~
where for the points ~
i;j are identical with p (4)
i;j except in the case (i;
and (i;
1;2 ) as the intersection between the line
(2)
(2)
(2)
(2)
(2)
(2)

Figure

3:
and the line joining the points p (2)

Figure

4.
Notice that these points ~
i;j are actually on one edge of a triangle of T " T 1 . We now
remark that in this procedure the three triangles K l above the line
are divided into four triangles in the following way: determine the two
points which are intersection between the line and the edges of K l ;
determine the mid point of the third edge (uniform subdivision in two parts). Using these
three points on the edges of K l and the vertices of K l , we divide K l into four triangles in
a standard way (see Figure 1). This will be the general rule.
~
~

Figure

4: our procedure
Now we can describe the passage from . The triangle of
containing the singular corner is divided into four triangles in Raugel's way: these triangles
are described by their three vertices
Any triangle K 2 T "T k above the line - 1
divided into four triangles in the
following way: First there exists i - 1 such that K is between the lines
Two vertices are on one line that we denote by and the
third one denoted by p 1 is on the other line. Secondly determine the two points p 0
0which are intersection between the line - 1
and the edges of K; determine
the mid point
1 of the third edge. Now the four triangles K l
are described by their three vertices (see Figure 5):
Remark that the triangle of containing the singular corner is also refined with
the same rule.
Let us finally notice that the above procedure guarantees the conformity of the meshes.
Now we want to show that this family of meshes is regular.
Lemma 5.1 The above family is regular in the sense that
Proof: To prove the assertion it suffices to look at the triangles of T " T k for any singular
triangle T of T 0 . Now we remark that our procedure preserves the acute property of the
angles. Therefore if we show that for all K 2
are the lengths of the edges of K in increasing order, then we deduce that
the smallest angle ff K of K satisfiesq
. sin(ff K
By Zl'amal's condition [22], we then deduce
ae K
which yields (27).
It then remains to prove (28). We now remark that if we apply a similarity of center at
the singular point and of ratio 2 \Gamma1=- to the triangulation we obtain the part of the
triangulation of below the line - 1 . This means that we are reduced
to prove (28) for the triangles above that line . Therefore we say that
only if K is between the lines
with
For any triangle K 2 ~
us denote by p K the length of the edge parallel to the
line
when K is between the lines
We first prove that
hK . p K . e 3( 1
Indeed we shall establish inductively that
r l
~ hK . p K .
l
1. It is clear that (30) holds for
Consequently to prove (30) for all k, it suffices to show that if (30) holds for k, it also
holds 1. Fix any K 2 ~
then as already explained it is divided into four
triangles K l 4, of ~
Two geometrical cases can be distinguished: either
1 is on the line - 1
1 is on the line - 1
us first show
that (30) holds for the triangles K l 4, in the first case. With the notation from

Figure

5, we deduce from the construction of the mesh that p 0
h is the similarity of center p 1 and ratio
This implies that
(a)
Figure

5: Definition of the nodes p i and p 0
Since by assumption K satisfies (30), K 1 and K 3 directly satisfies
r l
l
leading to (30) for K 1 (with 1. For the triangle K 3 ,
the above estimate yields
r K
r l
~ hK 3 . p K 3 . r K
l
where r
. This leads to (30) for K 3 because
due to the fact that i - 2 k\Gamma1 .
For K 2 and K 4 , we have p K Therefore by the inductive assumption
and the fact that ~ hK 2
r l
~ hK l . p K l
l
Again this leads to (30) for K 2 and K 4 because we easily check that (note that r - 1=2)
The second case is treated similarly, for K 3 we have the same estimate than before
K instead of r K that is the reason of the factor r \Gamma1
k+2 on the right-hand side. For
simply remark that the ratio ~
r of the second similarity is use
the fact that r k+2 - 2~r.
The proof of (30) is then complete.
Now (29) follows from (30) because using the fact that
log a
with
l=3
l
l=3
l - e 3( 1
Let us now come back to (28). For any K 2 ~
by construction of the mesh, we
clearly have
with the above notation for the vertices of K. On the other hand, since all the angles of
K are acute, if t denotes the orthogonal projection of p 1 on the edge we have by
3:
. e 6( 1
is the angle between the lines
Using the estimates (29), (31) and (32), we conclude that
This yields (28).
Remark 5.2 It was shown by Raugel in [16, p.96] that Raugel's graded meshes satisfy
Let us now show that our meshes satisfy the condition (4).
Lemma 5.3 The above family satisfies the condition (4) with
positive constant C independent of -.
Proof: As before it suffices to prove the assertion for the triangles in a fixed singular
triangle T of T 0 (since the remainder of the triangulation is quasi-uniform). By the
estimates (29), (31) and (32), we can claim that
Consequently we are reduced to estimate the quotient
~ hK l
when k - l, for any triangle K k l with K k ae K l . This quotient is
now easily estimated from above and from below by using the mean value theorem and
by distinguishing the case when K l contains the singular corner or not.
Remark 5.4 With our meshes, we have by Corollaries 3.4 and 4.2 and the two above
Lemmas that
(j
where C(-) is a positive constant which depends on e 6( 1
and then can blow
up as - tends to 0. This fact is confirmed by the numerical tests given in the next section.
6 Numerical results
In this Section, we present some numerical results which confirm our theoretical results
derived in Sections 3 and 4.
Let us consider boundary value problem (2) in a
domain\Omega with a re-entrant corner
(see

Figure

6).
It is well-known that the weak solution u of such a problem admits in the neighbourhood
of the singular point, i.e. in the neighbourhood of the re-entrant corner, the singular
representation
2(\Omega\Gamma3 the singular function
3 - in our example), and the stress intensity factor c (see, e.g., [13, 16]). Here, (r; ')
are polar coordinates with x
Using
graded meshes with a grading parameter -
one gets the optimal convergence order
of the finite element solution of problem (2). Figure 6 shows the mesh T 0 and the mesh

Figure

Domain\Omega with mesh T 0 and T 3
3 resulting from the mesh generation procedure described in Section 5 with the grading
parameter
Next we want to show by our experiments the dependence of the condition number
of the Galerkin stiffness matrix A j in the hierarchical basis on the number
levels used (Figure 7). In the experiments we use different values of the grading parameter
-. On can observe that is nearly a constant, and consequently, the
experiments confirm the theoretical estimate given in Corollary 3.4.0.51.52.53.54.51 2 3 4 5 6 7
(j

Figure

7: as a function of j

Figure

8 shows the behaviour of -(PMDS ) in dependence on the number
used. The numerical experiments confirm the statement given in Theorem 4.1.

Figure

8: -(PMDS ) as a function of j

Acknowledgement

We thank Dr. T. Apel for many discussions on this topics.



--R

spaces.
A Software Package for Solving Elliptic Partial Differential Equations - Users' Guide 7.0
A basic norm equivalence for the theory of multilevel methods.
New estimates for multilevel algorithms including the V-cycle
Parallel multilevel preconditioners.
The Finite Element Method for Elliptic Problems.
Approximation by finite element functions using local regularization.
Multilevel preconditioning.
Finite Element Methods for Navier-Stokes Equations: Theory and Algorithms
Multilevelmethoden als Iterationsverfahren - uber Erzeugendensystemen
Elliptic Problems in Nonsmooth Domains.
On adaptive grids in multilevel methods.

On discrete norm estimates related to multilevel preconditioners in the finite element method.
Multilevel Finite Element Approximation: Theory and Applications.
R'esolution num'erique de probl'emes elliptiques dans des domaines avec coins.
Finite element interpolation of nonsmooth functions satisfying boundary conditions.
An Analysis of the Finite Element Method.
Iterative methods by space decomposition and subspace correction.
On the multi-level splitting of finite element spaces
Multilevel Schwarz methods.
On the finite element method.
--TR
On the multi-level splitting of finite element spaces
Iterative methods by space decomposition and subspace correction
Finite Element Method for Elliptic Problems
graded meshes;finite element discretizations;multilevel methods;mesh refinement