<< . .

. 22
( : 59)

. . >>

plying (3.112) assumes the evaluability of v in ˆi , which is in the following
ensured by the continuity of v . This implies the same assumption for the
coe¬cients, since the shape functions Ni and their derivatives are continu-
ous. In order to ensure the numerical stability of a quadrature formula, it
is usually required that
ωi > 0 for all i = 1, . . . , R ,
ˆ (3.113)
which we will also do. Since all the considered ¬nite elements are such
that their faces with the enclosed degrees of freedom represent again a ¬-
nite element (in Rd’1 ) (see (3.42)), the boundary integrals are included
in a general discussion. In principle, di¬erent quadrature formulas can be
applied for each of the above integrals, but here we will disregard this pos-
sibility (with the exception of distinguishing between volume and boundary
integrals because of their di¬erent dimensions).
A quadrature formula on K generates a quadrature formula on a general
element K, recalling

v (ˆ) dˆ | det(B)|
v(x) dx = ˆx x

ωi,K v(bi,K ) ,
152 3. Finite Element Methods for Linear Elliptic Problems

:= F (ˆi ) are dependent on
where ωi = ωi,K = ωi | det(B)| and bi = bi,K
ˆ b
K. The positivity of the weights is preserved. Here, again F (ˆ) = B x + d
x ˆ
ˆ to K. The errors of the
denotes the a¬ne-linear transformation from K
quadrature formulas
ωi v (ˆi ) ,
v (ˆ) dˆ ’
E(ˆ) := ˆx x ˆ ˆb
K i=1
v(x) dx ’
EK (v) := ωi v(bi )
K i=1

are related to each other by
EK (v) = | det(B)|E(ˆ) .
ˆv (3.115)
The accuracy of a quadrature formula will be de¬ned by the requirement
that for l as large as possible,
E(ˆ) = 0 for p ∈ Pl (K)
ˆp ˆ
is satis¬ed, which transfers directly to the integration over K. A quadrature
formula should further provide the desired accuracy by using quadrature
nodes as less as possible, since the evaluation of the coe¬cient functions is
often expensive. In contrast, for the shape functions and their derivatives
a single evaluation is su¬cient. In the following we discuss some exam-
ples of quadrature formulas for the elements that have been introduced in
Section 3.3.
The most obvious approach consists in using nodal quadrature formu-
ˆ ˆˆ
las, which have the nodes a1 , . . . , aL of the reference element (K, P , Σ) as
ˆ ˆ
quadrature nodes. The requirement of exactness in P is then equivalent to

ωi =
ˆ Ni (ˆ) dˆ ,
xx (3.116)

so that the question of the validity of (3.113) remains.
We start with the unit simplex K de¬ned in (3.47). Here, the weights
of the quadrature formulas can be given directly on a general simplex K: If
the shape functions are expressed by their barycentric coordinates »i , the
integrals can be computed by
±1 !±2 ! · · · ±d+1 ! vol (K)
»±1 »±2 · · · »d+1 (x) dx =
(±1 + ±2 + · · · + ±d+1 + d)! vol (K)
1 2 ˆ

(see Exercise 3.28).
If P = P1 (K) and thus the quadrature nodes are the vertices, it follows
ωi = »i (x) dx = vol (K) for all i = 1, . . . , d + 1 . (3.118)
3.5. The Implementation of the Finite Element Method: Part 2 153

For P = P2 (K) and d = 2 we get, by the shape functions »i (2»i ’ 1), the
weights 0 for the nodes ai and, by the shape functions 4»i »j , the weights
ωi = vol (K) for bi = aij , i, j = 1, . . . , 3 , i > j ,
so that we have obtained here a quadrature formula that is superior to
(3.118) (for d = 2). However, for d ≥ 3 this ansatz leads to negative weights
and is thus useless. We can also get the exactness in P1 (K) by a single
quadrature node, by the barycentre (see (3.52)):
ω1 = vol (K) and b1 = aS = ai ,
d+1 i=1

which is obvious due to (3.117).
As a formula that is exact for P2 (K) and d = 3 (see [53]) we present
R = 4, ωi = 1 vol (K), and the bi are obtained by cyclic exchange of the
barycentric coordinates:
√ √ √ √
5’ 5 5’ 5 5’ 5 5+3 5
, , , .
20 20 20 20

On the unit cuboid K we obtain nodal quadrature formulas, which are
exact for Qk (K), from the Newton“Cˆtes formulas in the one-dimensional
situation by
i1 id
for ˆi1 ...id =
ωi1 · · · ωid
ωi1 ...id
ˆ = ˆ ˆ b ,..., (3.119)
k k
for ij ∈ {0, . . . , k} and j = 1, . . . , d .
Here the ωij are the weights of the Newton“Cˆtes formula for 0 f (x)dx
ˆ o
(see [30, p. 128]). As in (3.118), for k = 1 we have here a generalization
of the trapezoidal rule (cf. (2.38), (8.31)) with the weights 2’d in the 2d
vertices. From k = 8 on, negative weights arise. This can be avoided and
the accuracy for a given number of points increased if the Newton“Cˆtes o
integration is replaced by the Gauss“(Legendre) integration: In (3.119), ij /k
has to be replaced by the jth node of the kth Gauss“Legendre formula
(see [30, p. 156] there on [’1, 1]) and analogously ωij . In this way, by
(k + 1)d quadrature nodes the exactness in Q2k+1 (K), not only in Qk (K),
is obtained.
Now the question as to which quadrature formula should be chosen arises.
For this, di¬erent criteria can be considered (see also (8.29)). Here, we re-
quire that the convergence rate result that was proved in Theorem 3.29
should not be deteriorated. In order to investigate this question we have
to clarify which problem is solved by the approximation uh ∈ Vh based on
quadrature. To simplify the notation, from now on we do not consider
boundary integrals, that is, only Dirichlet and homogeneous Neumann
154 3. Finite Element Methods for Linear Elliptic Problems

boundary conditions are allowed. However, the generalization should be
clear. Replacing the integrals in (3.111) and (3.111) by quadrature formu-
las R ωi v (ˆi ) leads to some approximation Ah of the sti¬ness matrix
i=1 ˆ ˆ b
and q h of the load vector in the form
¯ ¯
Ah = (ah (•j , •i ))i,j , q h = (bh (•i ))i ,
for i, j = 1, . . . , M . Here the •i are the basis functions of Xh (see (3.101))
without taking into account the Dirichlet boundary condition and
ωl,K (K∇v · ∇w)(bl,K )
ah (v, w) :=
K∈Th l=1
ωl,K (c · ∇vw)(bl,K ) +
+ ωl,K (rvw)(bl,K )
K∈Th l=1 K∈Th l=1
for v, w ∈ Xh , (3.120)

for v ∈ Xh .
bh (v) := ωl,K (f v)(bl,K )
K∈Th l=1

The above-given mappings ah and bh are well-de¬ned on Xh — Xh and Xh ,
respectively, if the coe¬cient functions can be evaluated in the quadrature
nodes. Here we take into account that for some element K, ∇v for v ∈
Xh can have jump discontinuities on ‚K. Thus, for the quadrature nodes
bl,K ∈ ‚K in ∇v(bl,K ) we have to choose the value “belonging to bl,K ” that
corresponds to the limit of sequences in the interior of K. We recall that
in general ah and bh are not de¬ned for functions of V . Obviously, ah is
bilinear and bh is linear. If we take into account the analysis of incorporating
the Dirichlet boundary conditions in (3.99)“(3.106), we get a system of
equations for the degrees of freedom ξ = (ξ1 , . . . , ξM1 )T , which is equivalent
M1 ¯
to the variational equation on Vh for uh = i=1 ξi •i ∈ Vh :
ah (¯h , v) = bh (v) ’ ah (wh , v) for all v ∈ Vh
u (3.121)
with wh according to (3.105). As has been shown in (3.109), (3.121) is
equivalent, in the sense of the total approximation uh + wh of u + w, to the
¯h ∈ Vh ,
variational equation for u
ah (uh , v) = ¯h (v) := bh (v) ’ ah (Ih (w), v) for all v ∈ Vh ,
¯ b (3.122)
if this system of equations is uniquely solvable.

3.28 Prove equation (3.117) by ¬rst proving the equation for K = K
and then deducing from this the assertion for the general simplex by
Exercise 3.18.
3.6. Convergence Rate Results in Case of Quadrature and Interpolation 155

3.29 Let K be a triangle with vertices a1 , a2 , a3 . Further, let a12 , a13 , a23
denote the corresponding edge midpoints, a123 the barycenter and |K| the
area of K. Check that the quadrature formula
Qh (u) := 3 u(ai ) + 8 u(aij ) + 27u(a123 )
60 i=1 i<j

computes the integral Q(u) := u dx exactly for polynomials of third

3.6 Convergence Rate Results in the Case of
Quadrature and Interpolation
The purpose of this section is to analyze the approximation quality of a
solution uh + Ih (w) according to (3.122) and thus of uh + wh according to
(3.121) of the boundary value problem (3.12), (3.18)“(3.20).
Hence, we have left the ¬eld of Galerkin methods, and we have to
investigate the in¬‚uence of the errors
a ’ ah , b ’ a(w, ·) ’ bh + ah (Ih (w), ·).

To this end, we consider in general the variational equation in a normed
space (V, · )
u ∈ V satis¬es for all v ∈ V ,
a(u, v) = l(v) (3.123)
and the approximation in subspaces Vh ‚ V for h > 0,
uh ∈ Vh satis¬es for all v ∈ Vh .
ah (uh , v) = lh (v) (3.124)
Here a and ah are bilinear forms on V — V and Vh — Vh , respectively, and
l, lh are linear forms on V and Vh , respectively. Then we have the following
Theorem 3.38 (First Lemma of Strang)
Suppose there exists some ± > 0 such that for all h > 0 and v ∈ Vh ,
¤ ah (v, v) ,
±v (3.125)
and let a be continuous in V — V .
Then, there exists some constant C independent of Vh such that
|a(v, w) ’ ah (v, w)|
u ’ uh ¤ u ’ v + sup
C inf
v∈Vh w∈Vh
|l(w) ’ lh (w)|
+ sup .
156 3. Finite Element Methods for Linear Elliptic Problems

Proof: Let v ∈ Vh be arbitrary. Then it follows from (3.123)“(3.125) that
± uh ’ v ¤ ah (uh ’ v, uh ’ v)

= a(u ’ v, uh ’ v) + a(v, uh ’ v) ’ ah (v, uh ’ v)
+ lh (uh ’ v) ’ l(uh ’ v)
and moreover, by the continuity of a (cf. (3.2)),
|a(v, w) ’ ah (v, w)|
± uh ’ v ¤ M u ’ v + sup
|lh (w) ’ l(w)|
for v ∈ Vh .
+ sup

By means of u ’ uh ¤ u ’ v + uh ’ v and taking the in¬mum over
all v ∈ Vh , the assertion follows. 2

For ah = a and lh = l the assertion reduces to C´a™s lemma (Theo-
rem 2.17), which was the initial point for the analysis of the convergence
rate in Section 3.4. Here we can proceed analogously. For that purpose, the
following conditions must be ful¬lled additionally:
• The uniform Vh -ellipticity of ah according to (3.125) must be ensured.
• For the consistency errors
|a(v, w) ’ ah (v, w)|
Ah (v) := sup (3.127)

for an arbitrarily chosen comparison function v ∈ Vh and for
|l(w) ’ lh (w)|

the behavior in h must be analyzed.
The ¬rst requirement is not crucial if only a itself is V -elliptic and Ah
tends suitably to 0 for h ’ 0 :
Lemma 3.39 Suppose the bilinear form a is V -elliptic and there exists
some function C(h) with C(h) ’ 0 for h ’ 0 such that
Ah (v) ¤ C(h) v for v ∈ Vh .
Then there exists some h > 0 such that ah is uniformly Vh -elliptic for
h ¤ h.

Proof: By assumption, there exists some ± > 0 such that for v ∈ Vh ,
¤ ah (v, v) + a(v, v) ’ ah (v, v)
|a(v, v) ’ ah (v, v)| ¤ Ah (v) v ¤ C(h) v 2
3.6. Convergence Rate Results in Case of Quadrature and Interpolation 157

¯ ¯
Therefore, for instance, choose h such that C(h) ¤ ±/2 for h ¤ h. 2

We concretely address the analysis of the in¬‚uence of numerical quadra-
ture, that is, ah is de¬ned as in (3.120) and lh corresponds to ¯h in (3.122)
with the approximate linear form bh according to (3.120). Since this is an
extension of the convergence results (in · 1 ) given in Section 3.4, the as-
sumptions about the ¬nite element discretization are as summarized there
at the beginning. In particular, the triangulations Th consist of elements
that are a¬ne equivalent to each other. Furthermore, for a simpli¬cation of
the notation, let again d ¤ 3 and only Lagrange elements are considered. In
particular, let the general assumptions about the boundary value problems
which are speci¬ed at the end of Section 3.2.1 be satis¬ed.
According to Theorem 3.38, the uniform Vh -ellipticity of ah must be
ensured and the consistency errors (for an appropriate comparison element
v ∈ Vh ) must have the correct convergence behavior. If the step size h is
small enough, the ¬rst proposition is implied by the second proposition
by virtue of Lemma 3.39. Now, simple criteria that are independent of this
restriction will be presented. The quadrature formulas satisfy the properties
(3.112), (3.113) introduced in Section 3.5; in particular, the weights are
Lemma 3.40 Suppose the coe¬cient function K satis¬es (3.16) and let
c = 0 in „¦, let |“3 |d’1 > 0, and let r ≥ 0 in „¦. If P ‚ Pk (K) for the
ansatz space and if the quadrature formula is exact for P2k’2 (K), then ah
is uniformly Vh -elliptic.

Proof: Let ± > 0 be the constant of the uniform positive de¬niteness of
K(x). Then we have for v ∈ Vh :
ah (v, v) ≥ ± ωl,K |∇v|2 (bl,K ) = ± |∇v|2 (x) dx = ±|v|2 ,
K∈Th l=1

since |∇v|2 ∈ P2k’2 (K). The assertion follows from Corollary 3.14. 2

Further results of this type can be found in [9, pp. 194]. To investigate
the consistency error we can proceed similarly to the estimation of the
interpolation error in Section 3.4: The error is split into the sum of the errors
over the elements K ∈ Th and there transformed by means of (3.115) into
the error over the reference element K. The derivatives (in x) arising in the
error estimation over K are backtransformed by using Theorem 3.26 and
Theorem 3.27, which leads to the desired hK -factors. But note that powers
of B ’1 or similar terms do not arise. If the powers of det(B) arising in
both transformation steps cancel each other (which will happen), in this
way no condition about the geometric quality of the family of triangulations
arises. Of course, these results must be combined with estimates for the
158 3. Finite Element Methods for Linear Elliptic Problems

approximation error of Vh , for which, in particular, both approaches of
Section 3.4 (either regularity or maximum angle condition) are admissible.
For the sake of simplicity, we restrict our attention in the following to the
case of the polynomial ansatz space P = Pk (K). More general results of
similar type, in particular for triangulations with the cuboid element and
ˆ ˆ
P = Qk (K) as reference element, are summarized in [9, p. 207].
We recall the notation and the relations introduced in (3.114), (3.115)
for the local errors. In the following theorems we make use of the Sobolev
spaces W∞ on „¦ and on K with the norms · l,∞ and · l,∞,K , respectively,

and the seminorms | · |l,∞ and | · |l,∞,K , respectively. The essential local
assertion is the following:

<< . .

. 22
( : 59)

. . >>