<<

. 15
( 16)



>>

’= for h ’ 0 .
uM = = ’M
1’» ’»
M+1 » » 3

So the numerical solution does not converge to the exact solution at the
node aM .
Again this is no contradiction to possible convergence results for the ¬nite
di¬erence method in the discrete maximum norm, since now the di¬usion
coe¬cient is not ¬xed, but rather the discretization is to be viewed as
belonging to the limit case k = 0, with an arti¬cial di¬usion part in the
discretization (see (9.8) below).


Finite Di¬erence Methods with Symmetric and One-Sided Dif-
ference Quotients
The oscillations in Example 9.1 show that in this case no comparison prin-
ciple as in Corollary 1.13 is valid. Such a comparison principle, or more
strongly a maximum principle, will lead to nonnegative solutions in the
case of nonnegative right-hand side and Dirichlet data. This avoids for a
homogeneous right-hand side an undershooting, as observed in Example 9.1,
i.e., negative solution values in this case, and also an overshooting, i.e., so-
lution values larger than the maximum of the Dirichlet data, provided that
condition (1.32) (6)* holds.
In the following we will examine how the convective part in¬‚uences the
matrix properties (1.32) and thus the validity of a maximum or comparison
principle and also conclude a ¬rst simple remedy.
We consider the model problem (9.2), for simplicity on a rectangle
„¦ = (0, a) — (0, b), with constant, scalar K = kI and equipped with an
equidistant grid „¦h . To maintain the order of consistency 2 of a spatial
discretization of ’∇ · (K∇u) = ’k∆u by the ¬ve-point stencil, the use of
the symmetric di¬erence quotient for the discretization of

(c · ∇u)(x) = c1 (x)‚1 u(x) + c2 (x)‚2 u(x)

suggests itself, i.e., for a grid point x ∈ „¦h ,

1
c1 (x)‚1 u(x) ∼ c1 (x) (ui+1,j ’ ui’1,j ), (9.3)
2h

and similarly for c2 (x)‚2 u(x) (cf. (1.7) for the notation). This leads to
˜
the following entries of the system matrix Ah , for example in a rowwise
372 9. Discretization of Convection-Dominated Problems

numbering (compare (1.13)):
c1 (x) k
’ ’ 2;
left secondary diagonal:
2h h
c1 (x) k
’ 2;
right secondary diagonal: +
2h h
c2 (x) k
’ ’ 2;
l + 1 positions to the left:
2h h
c2 (x) k
’ 2;
l + 1 positions to the right: +
2h h
4k
diagonal: .
h2
Condition (1.32) (1) and (1.32) (6)* obviously hold.
We check the conditions su¬cient for a comparison principle (Corol-
lary 1.13). To satisfy condition (1.32) (2) we require
|c1 (x)|
k
’ + < 0,
h2 2h
|c2 (x)|
k
’ 2+ < 0.
h 2h
Denoting the grid P´clet number by
e
c ∞h
Peh := , (9.4)
2k
the above conditions are satis¬ed if
Peh < 1 (9.5)
is satis¬ed. Under this assumption also the conditions (1.32) (5) and (7)
are satis¬ed, and thus also (3), i.e., (9.5), is su¬cient for the validity of a
comparison principle. In Example 9.1 this is just the condition h < 2k.
The grid P´clet number is obviously related to the global P´clet number
e e
from (9.1) by
h
Peh = Pe .
2 diam(„¦)
The requirement (9.4) can always be met by choosing h su¬ciently small,
but for large Pe this may be a severe requirement, necessary for the sake of
stability of the method, whereas for the accuracy desired a larger step size
may be su¬cient. A simple remedy to ensure condition (1.32) (2) is to use
a one-sided (upwind) discretization of c1 ‚1 u and c2 ‚2 u, which is selected
against the stream direction de¬ned by c1 and c2 , respectively:
1
For c1 (x) ≥ 0 : c1 (x)‚1 u(x) ∼ c1 (x) (ui,j ’ ui’1,j ) , (9.6)
h
1
c1 (x)‚1 u(x) ∼ c1 (x) (ui+1,j ’ ui,j ) ,
for c1 (x) < 0 :
h
9.1. Standard Methods 373

and analogously for c2 ‚2 u.
Due to this choice there are only additional nonnegative addends to the
diagonal position and nonpositive ones to the o¬-diagonal positions com-
pared to the ¬ve-point stencil or another discretization of a di¬usive part.
Thus all properties (1.32) (1)“(7), (4)*, and (6)* remains una¬ected; i.e.,
the upwind discretization satis¬es all qualitative properties of Section 1.4
from the inverse monotonicity to the strong maximum principle, without
any restrictions to the local P´clet number.
e
The drawback lies in the reduced accuracy, since the one-sided di¬erence
quotient has only order of consistency 1. In Section 9.3 we will develop
more re¬ned upwind discretizations.
Due to
u(x, y) ’ u(x ’ h, y)
c1 = (9.7)
h
c1 h ’u(x ’ h, y) + 2u(x, y) ’ u(x + h, y) u(x + h, y) ’ u(x ’ h, y)
+ c1 ,
h2
2 2h
and analogously for the forward di¬erence quotient, the upwind dis-
cretization can be perceived as a discretization with symmetric di¬erence
quotients if a step-size-dependent di¬usive part, also discretized with ‚ ’ ‚ + ,
is added with the di¬usion coe¬cient
|c1 (x)|
h 0
Kh (x) := . (9.8)
|c2 (x)|
0
2
Therefore, one also speaks of adding arti¬cial di¬usion (or viscosity). The
disadvantage of this full upwind method is that it recognizes the ¬‚ow di-
rection only if the ¬‚ow is aligned to one of the coordinate axes. This will
be improved in Section 9.2.

Error Estimates for the Standard Finite Element Method
In order to demonstrate the theoretical de¬ciencies, we will again reproduce
the way for obtaining standard error estimates for a model problem. So let
K(x) ≡ µI with a constant coe¬cient µ > 0, c ∈ C 1 („¦, Rd ), r ∈ C(„¦),
f ∈ L2 („¦). Furthermore, assume that the following inequality is valid in
„¦, where r0 > 0 is a constant: r ’ 1 ∇ · c ≥ r0 .
2
Then the bilinear form a : V — V ’ R, V := H0 („¦), corresponding to
1

the boundary value problem (9.2), reads as (cf. (3.23))

[µ∇u · ∇v + c · ∇u v + r uv ] dx , u, v ∈ V .
a(u, v) := (9.9)
„¦

To get an ellipticity estimate of a, we set u = v ∈ V in (9.9) and take the
relation 2v(c · ∇v) = c · ∇v 2 into account. Then, by partial integration of
the middle term, we obtain

µ|v|2 + c · ∇v, v
a(v, v) = + rv, v
1 0 0
374 9. Discretization of Convection-Dominated Problems

1 1
µ|v|2 ’ ∇ · c, v 2 = µ|v|2 + r ’ ∇ · c, v 2
= + rv, v .
1 1
0
2 2
0 0

Introducing the so-called µ-weighted H 1 -norm by
2 1/2
:= µ|v|2 + v
v , (9.10)
µ 1 0

we immediately arrive at the estimate
a(v, v) ≥ µ|v|2 + r0 v ≥ ± v 2,
2
˜ (9.11)
1 0 µ

where ± := min{1, r0 } does not depend on µ.
˜
Due to c · ∇u = ∇ · (cu) ’ (∇ · c)u, partial integration yields for arbitrary
u, v ∈ V the identity
c · ∇u, v = ’ u, c · ∇v ’ (∇ · c)u, v .
0 0 0

So we get the continuity estimate
|a(u, v)| ¤ µ|u|1 |v|1 + c 0,∞ u 0|v|1 + (|c|1,∞ + r 0,∞ ) u 0 v 0
√ √
¤ ( µ |u|1 + u 0 ) {( µ + c 0,∞ ) |v|1 + (|c|1,∞ + r 0,∞ ) v 0 }
¤ M u µ v 1,
˜
(9.12)
˜ := 2 max{√µ + c 0,∞ , |c|1,∞ + r 0,∞ }.
where M
Since we are interested in the case of small di¬usion µ > 0 and present
˜
convection (i.e., c 0,∞ > 0), the continuity constant M can be bounded
independent of µ. It is not very surprising that the obtained continuity
estimate is nonsymmetric, since also the di¬erential expression L behaves
like that. Passing over to a symmetric estimate results in the following
relation:
˜
M
√ u µ v µ.
|a(u, v)| ¤
µ
Now, if Vh ‚ V denotes a ¬nite element space, we can argue as in the proof
of C´a™s lemma (Theorem 2.17) and get an error estimate for the corre-
e
sponding ¬nite element solution uh ∈ Vh . To do this, the nonsymmetric
continuity estimate (9.12) is su¬cient. Indeed, for arbitrary vh ∈ Vh , we
have
± u ’ uh ¤ a(u ’ uh, u ’ uh) = a(u ’ uh , u ’ vh ) ¤ M u ’ uh u ’ vh
˜
2
˜ .
µ 1
µ

Thus
˜
M
u ’ uh ¤ inf u ’ vh 1.
µ
± vh ∈Vh
˜
˜˜
Here the constant M /± does not depend on µ, h, and u. This estimate
is weaker than the standard estimate, because the µ-weighted H 1 -norm is
weaker than the H 1 -norm. Moreover, the error of the best approximation is
9.2. The Streamline-Di¬usion Method 375

not independent of µ, in general. For example, if we apply continuous, piece-
wise linear elements, then, under the additional assumption u ∈ H 2 („¦),
Theorem 3.29 yields the estimate
u ’ vh ¤ u ’ Ih (u) ¤ Ch|u|2 ,
inf 1 1
vh ∈Vh

where the constant C > 0 does not depend on µ, h, and u. So, we ¬nally
arrive at the relation
u ’ uh ¤ Ch|u|2 . (9.13)
µ

However, the H 2 -seminorm of the solution u depends on µ in a disadvan-
tageous manner; for example, it may be (cf. also [27, Lemma III.1.18])
that
|u|2 = O(µ’3/2 ) (µ ’ 0) .
This result is sharp, since for examples of boundary value problems for
ordinary linear di¬erential equations the error of the best approximation
already exhibits this asymptotic behaviour.
So the practical as well as the theoretical problems mentioned above
indicate the necessity to use special numerical methods for solving
convection-dominated equations. In the next sections, a small collection
of these methods will be depicted.


9.2 The Streamline-Di¬usion Method
The streamline-di¬usion method is the prevalent method in the numerical
treatment of stationary convection-dominated problems. The basic idea
is due to Brooks and Hughes [49], who called the method the streamline
upwind Petrov“Galerkin method (SUPG method).
We describe the idea of the method for a special case of boundary value
problem (9.2) under consideration. Let the domain „¦ ‚ Rd be a bounded
polyhedron. We consider the same model as in the preceding section, that
is, K(x) ≡ µI with a constant coe¬cient µ > 0, c ∈ C 1 („¦, Rd ), r ∈ C(„¦),
f ∈ L2 („¦). We also assume that the inequality r ’ 1 ∇ · c ≥ r0 is valid in „¦,
2
where r0 > 0 is a constant. Then the variational formulation of (9.2) reads
as follows:
Find u ∈ V such that
for all v ∈ V ,
a(u, v) = f, v (9.14)
0

where a is the bilinear form (9.9).
Given a regular family of triangulations {Th }, let Vh ‚ V denote the set
of continuous functions that are piecewise polynomial of degree k ∈ N and
satisfy the boundary conditions, i.e.,
Vh := vh ∈ V vh |K ∈ Pk (K) for all K ∈ Th . (9.15)
376 9. Discretization of Convection-Dominated Problems

If in addition the solution u ∈ V of (9.14) belongs to the space H k+1 („¦),
we have, by (3.87), the following error estimate for the interpolant Ih (u):
u ’ Ih (u) ¤ cint hk+1’l |u|k+1,K (9.16)
l,K K

for 0 ¤ l ¤ k+1 and all K ∈ Th . Since the spaces Vh are of ¬nite dimension,
a so-called inverse inequality can be proven (cf. Theorem 3.43, (2) and
Exercise 9.3):
cinv
∆vh 0,K ¤ |vh |1,K (9.17)
hK
for all vh ∈ Vh and all K ∈ Th . Here it is important that the constants
cint , cinv > 0 from (9.16) and (9.17), respectively, do not depend on u or vh
and on the particular elements K ∈ Th .
The basic idea of the streamline-di¬usion method consists in the addition
of suitably weighted residuals to the variational formulation (9.14). Because
of the assumption u ∈ H k+1 („¦), k ∈ N, the di¬erential equation can be
interpreted as an equation in L2 („¦). In particular, it is valid on any element
K ∈ Th in the sense of L2 (K), i.e.,
’µ∆u + c · ∇u + ru = f almost everywhere in K and for all K ∈ Th .
Next we take an elementwise de¬ned mapping „ : Vh ’ L2 („¦) and multiply
the local di¬erential equation in L2 (K) by the restriction of „ (vh ) to K.
Scaling by a parameter δK ∈ R and summing the results over all elements
K ∈ Th , we obtain
δK ’µ∆u + c · ∇u + ru, „ (vh ) = δK f, „ (vh ) .
0,K 0,K
K∈Th K∈Th

If we add this relation to equation (9.14) restricted to Vh , we see that the
weak solution u ∈ V © H k+1 („¦) satis¬es the following variational equation:
for all vh ∈ Vh ,
ah (u, vh ) = f, vh h

where
δK ’µ∆u + c · ∇u + ru, „ (vh )
ah (u, vh ) := a(u, vh ) + ,
0,K
K∈Th

f, vh := f, v + δK f, „ (vh ) .
h 0 0,K
K∈Th

Then the corresponding discretization reads as follows:
Find uh ∈ Vh such that
for all vh ∈ Vh .
ah (uh , vh ) = f, vh (9.18)
h

Corollary 9.2 Suppose the problems (9.14) and (9.18) have a solution u ∈
V © H k+1 („¦) and uh ∈ Vh , respectively. Then the following error equation
is valid:
ah (u ’ uh , vh ) = 0 for all vh ∈ Vh . (9.19)
9.2. Streamline-Di¬usion Method 377

In the streamline-di¬usion method (sdFEM), the mapping „ used in (9.18)
is chosen as „ (vh ) := c · ∇vh .
Without going into details, we mention that a further option is to set
„ (vh ) := ’µ∆vh + c · ∇vh + rvh . This results in the so-called Galerkin/least
squares“FEM (GLSFEM) [54].
Especially with regard to the extension of the method to other ¬nite ele-
ment spaces, the discussion of how to choose „ and δK is not yet complete.

Interpretation of the Additional Term in the Case of Linear Ele-
ments
If the ¬nite element spaces Vh are formed by piecewise linear functions (i.e.,
in the above de¬nition (9.15) of Vh we have k = 1), we get ∆vh |K = 0 for
all K ∈ Th . If in addition there is no reactive term (i.e., r = 0), the discrete
bilinear form is


µ∇uh ·∇vh dx+ c · ∇uh , vh 0 + δK c · ∇uh , c · ∇vh
ah (uh , vh ) = .
0,K
„¦ K∈Th


Since the scalar product appearing in the sum can be rewritten as
c · ∇uh , c · ∇vh 0,K = K (ccT ∇uh ) · ∇vh dx , we obtain the following
equivalent representation:


(µI + δK ccT )∇uh · ∇vh dx + c · ∇uh , vh
ah (uh , vh ) = .
0
K
K∈Th


This shows that the additional term introduces an element-dependent extra
di¬usion in the direction of the convective ¬eld c (cf. also Exercise 0.3),
which motivates the name of the method. In this respect, the streamline-
di¬usion method can be understood as an improved version of the full
upwind method, as seen, for example, in (9.6).

Analysis of the Streamline-Di¬usion Method
To start the analysis of stability and convergence properties of the
streamline-di¬usion method, we consider the term ah (vh , vh ) for arbitrary
vh ∈ Vh .
As in Section 3.2.1, the structure of the discrete bilinear form ah allows
us to derive the estimate

ah (vh , vh ) ≥ µ|vh |2 +r0 vh ’µ∆vh + c · ∇vh + rvh , c · ∇vh
2
0+ δK .
1 0,K
K∈Th


Furthermore, neglecting for a moment the second term in the sum and
using the elementary inequality ab ¤ a2 + b2 /4 for arbitrary a, b ∈ R, we
378 9. Discretization of Convection-Dominated Problems

get

δK ’µ∆vh + rvh , c · ∇vh 0,K
K∈Th

¤ ’µ |δK | ∆vh , |δK | c · ∇vh
0,K
K∈Th

|δK | rvh , |δK | c · ∇vh
+
0,K

¤ µ2 |δK | ∆vh + |δK | r
2 2 2
vh
0,K 0,∞,K 0,K
K∈Th
|δK |
c · ∇vh 2
+ .
0,K
2
By means of the inverse inequality (9.17) it follows that

c2
δK ’µ∆vh + rvh , c · ∇vh ¤ µ2 |δK | 2 |vh |1,K
inv 2
0,K
hK
K∈Th K∈Th
|δK |
+|δK | r c · ∇vh 2
2 2
vh + .
0,K
0,∞,K 0,K
2
Putting things together, we obtain
c2
ah (vh , vh ) ≥ µ ’ µ |δK | inv |vh |2
2 2
1,K vh 0,K
h2
K
K∈Th
|δK |
+ r0 ’ |δK | r + δK ’ c · ∇vh
2 2
.
0,∞,K 0,K
2
The choice
h2
1 r0
¤ min K
0 < δK 2 , r2 (9.20)
2 µcinv 0,∞,K

leads to
µ r0 1
ah (vh , vh ) ≥ |vh |2 + δK c · ∇vh
2 2
vh + .
1 0 0,K
2 2 2
K∈Th

Therefore, if the so-called streamline-di¬usion norm is de¬ned by
1/2

δK c · ∇v v∈V ,
µ|v|2 + r0 v 2 2
v := + ,
sd 1 0 0,K
K∈Th

then the choice (9.20) implies the estimate
1
¤ ah (vh , vh ) for all vh ∈ Vh .
2
vh (9.21)
sd
2
9.2. Streamline-Di¬usion Method 379

Obviously, the streamline-di¬usion norm · sd is stronger than the
µ-weighted H 1 -norm (9.10); i.e.,

min{1, r0 } v µ ¤ v sd for all v ∈ V .

Now an error estimate can be proven. Since estimate (9.21) holds only on
the ¬nite element spaces Vh , we consider ¬rst the norm of Ih (u) ’ uh ∈ Vh
and make use of the error equation (9.19):
1
Ih (u) ’ uh ¤ ah (Ih (u) ’ uh , Ih (u) ’ uh ) = ah (Ih (u) ’ u, Ih (u) ’ uh ) .
2
sd
2
In particular, under the assumption u ∈ V © H k+1 („¦) the following three
estimates are valid:

∇(Ih (u) ’ u) · ∇(Ih (u) ’ uh ) dx ¤ µ |Ih (u) ’ u|1 Ih (u) ’ uh
µ sd
„¦

¤ cint µ hk |u|k+1 Ih (u) ’ uh ,
sd



[c · ∇(Ih (u) ’ u) + r(Ih (u) ’ u)](Ih (u) ’ uh ) dx
„¦

(r ’ ∇ · c)(Ih (u) ’ u)(Ih (u) ’ uh ) dx
=
„¦

’ (Ih (u) ’ u) c · ∇(Ih (u) ’ uh ) dx
„¦

¤ r’∇·c Ih (u) ’ u Ih (u) ’ uh
0,∞ 0 0

+ Ih (u) ’ u c · ∇(Ih (u) ’ uh )
0 0
±
 1/2

¤ Ih (u) ’ u 2
C
 0,K
K∈Th

1/2 
’1
δK Ih (u) ’ u Ih (u) ’ uh
2
+ sd

0,K
K∈Th
1/2
’1
¤ 1 + δK h2 |u|2 Ih (u) ’ uh
Chk ,
sd
K k+1,K
K∈Th

and

δK ’µ∆(Ih (u) ’ u) + c · ∇(Ih (u) ’ u)
K∈Th

+ r(Ih (u) ’ u), c · ∇(Ih (u) ’ uh ) (9.22)
0,K
380 9. Discretization of Convection-Dominated Problems

¤ δK µhk’1 + c k+1
k
cint 0,∞,K hK +r 0,∞,K hK
K
K∈Th
— |u|k+1,K δK c · ∇(Ih (u) ’ uh ) 0,K


1/2
2
¤C |u|2 Ih (u) ’ uh
δK µhk’1 + hk + hk+1 sd .
K k+1,K
K K
K∈Th

Condition (9.20), which was already required for estimate (9.21), implies
that
h2
¤ 2K ,
µδK
cinv
and so the application to the ¬rst term of the last bound leads to

δK ’ µ∆(Ih (u) ’ u) + c · ∇(Ih (u) ’ u)
K∈Th
+ r(Ih (u) ’ u), c · ∇(Ih (u) ’ uh ) 0,K

1/2

¤ Ch δK ] |u|2 Ih (u) ’ uh
k
[µ + .
sd
k+1,K
K∈Th

Collecting the estimates and dividing by Ih (u) ’ uh sd , we obtain the
relation
1/2
h2
Ih (u) ’ uh ¤ Chk µ + K + h2 + δK |u|2 .
sd K k+1,K
δK
K∈Th

Finally, the terms in the square brackets will be equilibrated with the help
of condition (9.20). We rewrite the µ-dependent term in this condition as
h2 2
K
=2 PeK hK
2
µcinv cinv c ∞,K

with
c ∞,K hK
PeK := . (9.23)

This local P´clet number is a re¬nement of the de¬nition (9.4).
e
The following distinctions concerning PeK are convenient:
PeK ¤ 1 and PeK > 1 .
In the ¬rst case, we choose
h2 2
= δ1 K ,
δK = δ0 PeK hK δ0 = δ1 ,
µ c ∞,K
9.2. Streamline-Di¬usion Method 381

with appropriate constants δ0 > 0 and δ1 > 0, respectively, which are
independent of K and µ. Then we have
h2 1 2PeK
hK ¤ C(µ + hK ) ,
+ h2 + δK = µ + h2 + δ1
K
µ+ 1+
K K
δK δ1 c 0,∞,K
where C > 0 is independent of K and µ. In the second case, it is su¬cient to
choose δK = δ2 hK with an appropriate constant δ2 > 0 that is independent
of K and µ. Then
δ2 c 0,∞,K h2
δ2 K
δK = PeK hK =
PeK 2PeK µ
and
h2 1
+ δ2 hK + h2 ¤ C(µ + hK ) ,
µ + K + h2 + δK = µ +
K K
δK δ2
with C > 0 independent of K and µ. Note that in both cases the constants
can be chosen su¬ciently small, independent of PeK , that the condi-
tion (9.20) is satis¬ed. Now we are prepared to prove the following error
estimate.
Theorem 9.3 Let the parameters δK be given by
±
 h2
δ1 K , PeK ¤ 1 ,
δK =
 δ hµ , Pe > 1 ,
2K K

where δ1 , δ2 > 0 do not depend on K and µ and are chosen such that
condition (9.20) is satis¬ed. If the weak solution u of (9.14) belongs to
H k+1 („¦), then


u ’ uh sd ¤ C µ + h hk |u|k+1 ,

where the constant C > 0 is independent of µ, h, and u.

Proof: By the triangle inequality, we get
u ’ uh ¤ u ’ Ih (u) + Ih (u) ’ uh sd .
sd sd

An estimate of the second addend is already known. To deal with the ¬rst
term, the estimates of the interpolation error (9.16) are used directly:


u ’ Ih (u) 2
sd

µ|u ’ Ih (u)|2 + r0 u ’ Ih (u) δK c · ∇(u ’ Ih (u))
2 2
= +
1 0 0,K
K∈Th
2(k+1)
¤ |u|2
c2 µh2k + r0 hK 2 2k
+ δK c 0,∞,K hK
int K k+1,K
K∈Th
382 9. Discretization of Convection-Dominated Problems

¤ µ + h2 + δK |u|2
k+1,K ¤ C(µ + h)hK |u|k+1 .
Ch2k 2k 2
K K
K∈Th

2

Remark 9.4 (i) In the case of large local P´clet numbers, we have
e
µ¤ 1 c ∞,K hK and thus
2

1/2

u ’ uh hK c · ∇(u ’ uh ) ¤ Chk+1/2 |u|k+1 .
2
+ δ2
0 0,K
K∈Th


So the L2 -error of the solution is not optimal in comparison with the
estimate of the interpolation error

u ’ Ih (u) ¤ Chk+1 |u|k+1 ,
0


whereas the L2 -error of the directional derivative of u in the direction
of c is optimal.

(ii) In general, the seminorm |u|k+1 depends on negative powers of µ.
Therefore, if h ’ 0, the convergence in Theorem 9.3 is not uniform
with respect to µ.

Comparing the estimate from Theorem 9.3 for the special case of continuous
linear elements with the estimate (9.13) for the corresponding standard
method given at the end of the introduction, i.e.,

u ’ uh ¤ Ch|u|2 ,
µ


we see that the error of the streamline-di¬usion method is measured in a
stronger norm than the · µ -norm and additionally, that the error bound
is asymptotically better in the interesting case µ < h . A further advan-
tage of the streamline-di¬usion method is to be seen in the fact that its
implementation is not much more di¬cult than that of the standard ¬nite
element method.
However, there are also some disadvantages: Since the error bound in-
volves the H k+1 -seminorm of the solution u, it may depend on negative
powers of µ. Furthermore, there is no general rule to determine the param-
eters δ1 , δ2 . Usually, they are chosen more or less empirically. This may be
a problem when the streamline-di¬usion method is embedded into more
complex programs (for example, for solving nonlinear problems). Finally,
in contrast to the ¬nite volume methods described in the next section, the
property of inverse monotonicity (cf. Theorem 6.19) cannot be proven in
general.
9.3. Finite Volume Methods 383

Exercises
9.1 (a) Given a constant di¬usion coe¬cient µ > 0, rewrite the ordinary
boundary value problem
(’µu + u) = 0 in „¦ := (0, 1) ,
u(0) = u(1) ’ 1 = 0 ,
into an equivalent form but with nonnegative right-hand side and
homogeneous Dirichlet boundary conditions.
(b) Compute the H 2 (0, 1)-seminorm of the solution of the transformed
problem and investigate its dependence on µ.

9.2 Prove the error equation of the streamline-di¬usion method (Corollary
9.2).

9.3 Given an arbitrary, but ¬xed, triangle K with diameter hK , prove the
inequality
cinv
∆p 0,K ¤ |p|1,K
hK
for arbitrary polynomials p ∈ Pk (K), k ∈ N, where the constant cinv > 0
is independent of K and p.

·
9.4 Verify that the streamline-di¬usion norm is indeed a norm.
sd




9.3 Finite Volume Methods
In the convection-dominated situation, the ¬nite volume method intro-
duced in Chapter 6 proves to be a very stable, but not so accurate, method.
One reason for this stability lies in an appropriate asymptotic behaviour of
the weighting function R for large absolute values of its argument.
Namely, if we consider the examples of nonconstant weighting functions
given in Section 6.2.2, we see that
(P4) lim R(z) = 0 , lim R(z) = 1 .
z’’∞ z’∞

In the general case of the model problem (6.5) with k = µ > 0, (P4)
γij dij
’1 the term rij ui + (1 ’ rij )uj in the bilinear
implies that for
µ
γij dij
form bh e¬ectively equals uj , whereas in the case 1 the quantity
µ
ui remains.
In other words, in the case of dominating convection, the approximation
bh evaluates the “information” (uj or ui ) upwind, i.e., just at that node (aj
or ai ) from which “the ¬‚ow is coming”.
384 9. Discretization of Convection-Dominated Problems

This essentially contributes to the stabilization of the method and makes
it possible to prove properties such as global conservativity or inverse mono-
tonicity (cf. Section 6.2.4) without any restrictions on the size of the local
γij dij
P´clet number
e and thus without any restrictions on the ratio of h
µ
and µ. This local P´clet number (note the missing factor 2 in comparison
e
to (9.23)) also takes the direction of the ¬‚ow compared to the edge ai aj
into account.

The Choice of Weighting Parameters
In order to motivate the choice of the weighting parameters in the case of
the Voronoi diagram, we recall the essential step in the derivation of the
¬nite volume method, namely the approximation of the integral

[µij (νij · ∇u) ’ γij u] dσ .
Iij :=
“ij

It ¬rst suggests itself to apply a simple quadrature rule, for example
Iij ≈ qij mij ,
where qij denotes the value of the expression to be integrated at the
point aij of the intersection of the boundary segment “ij with the edge
bounded by the vertices ai and aj (i.e., 2aij = ai + aj ). Next, if this edge
is parametrised according to
11
„∈ ’ ,
x = x(„ ) = aij + „ dij νij , ,
22
and if we introduce the composite function w(„ ) := u(x(„ )) , then we can
write
µij dw
µij (νij · ∇u) ’ γij u = q(0) with q(„ ) := („ ) ’ γij w(„ ) .
dij d„
The relation de¬ning the function q can be interpreted as a linear ordinary
di¬erential equation for the unknown function w : ’ 2 , 1 ’ R . Provided
1
2
that q is continuous on the interval ’ 2 , 1 , the equation can be solved
1
2
exactly:

dij γij dij 1 1
q(s) exp ’ ds + w ’
w(„ ) = s+
µij µij 2 2
’1/2

γij dij 1
— exp „+ .
µij 2
Approximating q by a constant qij , we get in the case γij = 0,
qij γij dij 1 1
≈ 1 ’ exp ’ +w ’
w(„ ) „+
γij µij 2 2
γij dij 1
— exp „+ .
µij 2
9.3. Finite Volume Methods 385

In particular,
1 qij γij dij 1 γij dij
≈ 1 ’ exp ’ +w ’
w exp ; (9.24)
2 γij µij 2 µij
that is, the approximation qij of q(0) can be expressed by means of the
values w(± 1 ):
2

γij dij
w( 1 ) ’ w(’ 1 ) exp
2 2 µij
qij ≈ γij . (9.25)
γij dij
’1
exp µij

In the case γij = 0, it immediately follows from the exact solution and the
approximation q ≈ qij that

w( 1 ) ’ w(’ 1 )
qij ≈ µij 2 2
.
dij
Since this is equal to the limit of (9.25) for γij ’ 0, we can exclusively
work with the representation (9.25).
If we de¬ne the weighting function R : R ’ [0, 1] by
1 z
R(z) := 1 ’ 1’ , (9.26)
ez ’ 1
z
γij dij
then with the choice rij := R , (9.25) can be written as
µij

uj ’ ui
qij ≈ µij ’ [rij ui + (1 ’ rij ) uj ] γij .
dij
A simple algebraic manipulation shows that this is exactly the approxima-
tion scheme given in Section 6.2.
The use of the weighting function (9.26) yields a discretization method
that can be interpreted as a generalization of the so-called Il™in“ Allen“
Southwell scheme. However, in order to avoid the comparatively expensive
computation of the function values rij of (9.26), often simpler functions
R : R ’ [0, 1] are used (see Section 6.2.2), which are to some extent
approximations of (9.26) keeping the properties (P1) to (P4).
At the end of this paragraph we will illustrate the importance of the prop-
erties (P1) to (P3), especially for convection-dominated problems. Property
(P2) has been used in the proof of the basic stability estimate (6.20). On
the other hand, we have seen at several places (e.g., in Section 1.4 or in
Chapter 5) that the matrix Ah of the corresponding system of linear al-
gebraic equations should have positive diagonal entries. For example, if
in the di¬erential equation from (9.2) the reaction term disappears, then
properties (P1) and (P3) guarantee that the diagonal entries are at least
nonnegative. This can be seen as follows:
386 9. Discretization of Convection-Dominated Problems

From (6.9) we conclude the following formula:
µij µij γij dij
i ∈ Λ.
(Ah )ii = + γij rij mij = 1+ rij mij ,
dij dij µij
If we replace in property (P3) the number z by ’z, then we get, by property
(P1),
0 ¤ 1 + [1 ’ R(’z)] z = 1 + zR(z) .
Therefore, if the weighting function R satis¬es (P1) and (P3), then we have
that (Ah )ii ≥ 0 for all i ∈ Λ .
The simple choice rij ≡ 1 does not satisfy property (P3). In this case,
2
the condition (Ah )ii ≥ 0 leads to the requirement
γij dij
’ ¤ 1,
2µij
which in the case γij ¤ 0, i.e., for a local ¬‚ow from aj to ai , is a restriction
to the ratio of h and µ, and this is analogous to the condition (9.5) on the
grid P´clet number, where only the sizes of K, c, and h enter.
e
Similarly, it can be shown that property (P3) implies the nonpositivity
of the o¬-diagonal entries of Ah .

An Error Estimate
At the end of this section an error estimate will be cited, which can be
derived similarly to the corresponding estimate of the standard method.
The only special aspect is that the dependence of the occurring quantities
on µ is carefully tracked (see [40]).
Theorem 9.5 Let {Th }h be a regular family of conforming triangula-
tions, all triangles of which are nonobtuse. Furthermore, in addition to the
assumptions on the coe¬cients of the bilinear form (9.9), let f ∈ C 1 („¦).
If the exact solution u of the model problem belongs to H 2 („¦) and if
uh ∈ Vh denotes the approximative solution of the ¬nite volume method
(6.11), where the approximations γij , respectively ri , are chosen according
¯
to (6.7), respectively (6.8), then for su¬ciently small h > 0 the estimate
h
¤C√ [ u h ∈ (0, ¯ ,
u ’ uh + |f |1,∞ ] , h]
µ 2
µ
¯
holds, where both the constant C > 0 and h > 0 do not depend on µ.
In special, but practically not so relevant, cases (for example, if the trian-
gulations are of Friedrichs“Keller type), it is possible to remove the factor
1
√ in the bound above.
µ
Comparing the ¬nite volume method with the streamline-di¬usion
method, we see that the ¬nite volume method is less accurate. However, it
is globally conservative and inverse monotone.
9.4. The Lagrange“Galerkin Method 387

Exercise
9.5 Using an equidistant grid, formulate both the streamline-di¬usion
method and the ¬nite volume method for a one-dimensional model problem
(d = 1, „¦ = (0, 1), r = 0) with constant coe¬cients and compare the re-
sulting discretizations. Based on that comparison, what can be said about
the choice of the parameters in the streamline-di¬usion method?



9.4 The Lagrange“Galerkin Method
In the previous sections, discretization methods for stationary di¬usion-
convection equations were presented. In conjunction with the method of
lines, these methods can also be applied to parabolic problems. However,
since the method of lines decouples spatial and temporal variables, it can-
not be expected that the peculiarities of nonstationary di¬usion-convection
equations are re¬‚ected adequately.
The so-called Lagrange“Galerkin method attempts to bypass this prob-
lem by means of an intermediate change from the Eulerian coordinates
(considered up to now) to the so-called Lagrangian coordinates. The latter
are chosen in such a way that the origin of the coordinate system (i.e., the
position of the observer) is moved with the convective ¬eld, and in the new
coordinates no convection occurs.
To illustrate the basic idea, the following initial-boundary value problem
will be considered, where „¦ ‚ Rd is a bounded domain with Lipschitz
continuous boundary and T > 0:
For given functions f : QT ’ R and u0 : „¦ ’ R , ¬nd a function
u : QT ’ R such that
‚u
+ Lu = f in QT ,
‚t
(9.27)
u=0 on ST ,
„¦ — {0} ,
u = u0 on
where

(Lu) (x, t) := ’∇·(K(x) ∇u(x, t))+c(x, t)·∇u(x, t)+r(x, t)u(x, t), (9.28)

with su¬ciently smooth coe¬cients

K : „¦ ’ Rd,d , c : Q T ’ Rd , r : QT ’ R .

As usual, the di¬erential operators ∇ and ∇· act only with respect to the
spatial variables.
The new coordinate system is obtained by solving the following
parameter-dependent auxiliary problem:
388 9. Discretization of Convection-Dominated Problems

Given (x, s) ∈ QT , ¬nd a vector ¬eld X : „¦ — [0, T ]2 ’ Rd such that
d
t ∈ (0, T ) ,
X(x, s, t) = c(X(x, s, t), t) ,
(9.29)
dt
X(x, s, s) = x .
The trajectories X(x, s, ·) are called characteristics (through (x, s)). If c
is continuous on QT and, for ¬xed t ∈ [0, T ], Lipschitz continuous with
respect to the ¬rst argument on „¦, then there exists a unique solution
X = X(x, s, t). Denoting by u the su¬ciently smooth solution of (9.27)
and setting
u(x, t) := u(X(x, s, t), t) for ¬xed s ∈ [0, T ] ,
ˆ
then the chain rule implies that
‚u
ˆ ‚u
+ c · ∇u (X(x, s, t), t) .
(x, t) =
‚t ‚t
The particular value
‚u
ˆ ‚u
(x, s) + c(x, s) · ∇u(x, s)
(x, s) =
‚t ‚t
is called the material derivative of u at (x, s). Thus the di¬erential equation
reads as
‚u ˆ
’ ∇ · (K ∇u) + ru = f ;
‚t
i.e., it is formally free of any convective terms.
Now the equation will be semidiscretized by means of the horizontal
method of lines. A typical way is to approximate the time derivative by
backward di¬erence quotients. So let an equidistant partition of the time
interval (0, T ) with step size „ := T /N, N ∈ N (provided that T < ∞), be
given.
Tracking the characteristics backwards in time, in the strip „¦ —
[tn , tn+1 ), n ∈ {0, 1, . . . , N ’ 1}, with x = X(x, tn+1 , tn+1 ) the following
approximation results:
‚u ˆ 1 1
≈ [ˆ(x, tn+1 ) ’ u(x, tn )] = [u(x, tn+1 ) ’ u(X(x, tn+1 , tn ), tn )] .
u ˆ
‚t „ „
Further, if Vh denotes a ¬nite-dimensional subspace of V in which we want
to ¬nd the approximations to u(·, tn ), the method reads as follows:
Given u0h ∈ Vh , ¬nd an element U n+1 ∈ Vh , n ∈ {0, . . . , N ’ 1}, such
that
1
U n+1 ’ U n (X(·, tn+1 , tn )), vh 0

+ K∇U n+1 · ∇vh , 1 0 + r(·, tn+1 )U n+1 , vh 0 = f (·, tn+1 ), vh 0
for all vh ∈ Vh ,
U 0 = u0h .
(9.30)
9.4. Lagrange“Galerkin Method 389

A possible extension of the method is to use time-dependent subspaces;
that is, given a sequence of subspaces Vh ‚ V, n ∈ {0, . . . , N }, the
n

approximations U n to u(·, tn ) are chosen from Vh .n

So the basic idea of the Lagrange“Galerkin method, namely, the elimi-
nation of convective terms by means of an appropriate transformation of
coordinates, allows the application of standard discretization methods and
makes the method attractive for situations where convection is dominating.
In fact, there exists a whole variety of papers dealing with error esti-
mates for the method in the convection-dominated case, but often under
the condition that the system (9.29) is integrated exactly.
In practice, the exact integration is impossible, and the system (9.29) has
to be solved numerically (cf. [61]). This may lead to stability problems, so
there is still a considerable need in the theoretical foundation of Lagrange“
Galerkin methods.
Only recently it has been possible for a model situation to prove order
of convergence estimates uniformly in the P´clet number for (9.30) (see
e
[43]). The key is the consequent use of Lagrangian coordinates, reveal-
ing that (9.30) is just the application of the implicit Euler method to an
equation arising from a transformation by characteristics de¬ned piecewise
backward in time. This equation is a pure di¬usion problem, but with a
coe¬cient re¬‚ecting the transformation. In conjunction with the backward
Euler method this is not visible in the elliptic part to be discretized. Thus
tracking the characteristics backward in time turns out to be important.
A
Appendices




A.1 Notation
C set of complex numbers
N set of natural numbers
N0 := N ∪ {0}
Q set of rational numbers
R set of real numbers
R+ set of positive real numbers
Z set of integers
z real part of the complex number z
z imaginary part of the complex number z
transpose of the vector x ∈ Rd , d ∈ N
T
x
1/p
d
|x|p j=1 |xj | , x = (x1 , . . . , xd )T ∈ Rd , d ∈ N, p ∈ [1, ∞)
p
:=
|x|∞ := maxj=1,...,d |xj | maximum norm of the vector x ∈ Rd , d ∈ N
|x| := |x|2 Euclidean norm of the vector x ∈ Rd , d ∈ N
d
x·y := xT y = j=1 xj yj scalar product of the vectors x, y ∈ Rd
:= y T Ax = y · Ax energy product of the vectors x, y ∈ Rd w.r.t.
x, y A
a symmetric, positive de¬nite matrix A
|±| := |±|1 order (or length) of the multi-index ± ∈ Nd , d ∈ N 0
I identity matrix or identity operator
jth unit vector in Rm , j = 1, . . . , m
ej
diagonal matrix in Rm,m with diagonal
diag(»i ) = diag(»1 , . . . , »m )
entries »1 , . . . , »m ∈ C
A.1. Notation 391

AT transpose of the matrix A
A’T transpose of the inverse matrix A’1
det A determinant of the square matrix A
»min (A) minimum eigenvalue of a matrix A with real eigenvalues
»max (A) maximum eigenvalue of a matrix A with real eigenvalues
σ(A) set of eigenvalues (spectrum) of the square matrix A
(A) spectral radius of the square matrix A
m(A) bandwidth of the symmetric matrix A
Env (A) hull of the square matrix A
p(A) pro¬le of the square matrix A
:= {x : x ’ x0 < } open ball in a normed space
B (x0 )
:= {x : x ’ x0 ¤ } closed ball in a normed space
B (x0 )
diameter of the set G ‚ Rd
diam (G)
|G|n n-dimensional (Lebesgue) measure of the G ‚ Rn , n ∈ {1, . . . , d}
|G| := |G|d d-dimensional (Lebesgue) measure of the set G ‚ Rd
vol (G) length (d = 1), area (d = 2), volume (d = 3) of “geometric bodies”
G ‚ Rd
int G interior of the set G
‚G boundary of the set G
G closure of the set G
span G linear hull of the set G
conv G convex hull of the set G
|G| cardinal number of the discrete set G
outer unit normal w.r.t. the set G ‚ Rd
ν
domain of Rd , d ∈ N
„¦
:= ‚„¦ boundary of the domain „¦ ‚ Rd

supp • support of the function •
’1
f inverse of the mapping f
f [G] image of the set G under the mapping f
’1
f [G] preimage of the set G under the mapping f
f |K restriction of f : G ’ R to a subset K ‚ G
vX norm of the element v of the normed space X
dim X dimension of the ¬nite-dimensional linear space X
L[X, Y ] set of linear, continuous operators acting from the normed space
X in the normed space Y
:= L[X, R] dual space of the real normed space X
X
O(·), o(·) Landau symbols of asymptotic analysis
(i, j ∈ N0 ) Kronecker symbol, i.e., δii = 1 and δij = 0 if i = j
δij

Di¬erential expressions
(l ∈ N) symbol for the partial derivative w.r.t. the lth variable
‚l
(t ∈ R) symbol for the partial derivative w.r.t. the variable t
‚t
(± ∈ Nd multi-index) ±th partial derivative
‚± 0
∇ := (‚1 , . . . , ‚d )T Nabla operator (symbolic vector)
392 A. Appendices

∆ Laplace operator
:= µ · ∇ directional derivative w.r.t. the vector µ
‚µ
:= ‚¦ := (‚j ¦i )m
D¦ Jacobi matrix or functional matrix
i,j=1
‚x
of a di¬erentiable mapping ¦ : Rm ’ Rm
Coe¬cients in di¬erential expressions
K di¬usion coe¬cient (a square matrix function)
c convection coe¬cient (a vector function)
r reaction coe¬cient

Discretization methods
Vh ansatz space
Xh extended ansatz space without any homogeneous Dirichlet
boundary conditions
ah approximated bilinear form
bh approximated linear form

Function spaces (see also Appendix A.5)
Pk (G) set of polynomials of maximum degree k on G ‚ Rd
C(G) = C 0 (G) set of continuous functions on G
(l ∈ N) set of l-times continuously di¬erentiable functions on G
C l (G)

C (G) set of in¬nitely often continuously di¬erentiable functions on G
0
C(G) = C (G) set of bounded and uniformly continuous functions on G
(l ∈ N) set of functions with bounded and uniformly continuous
C l (G)
derivatives up to the order l on G

C (G) set of functions, all partial derivatives of which are bounded and
uniformly continuous on G
0
C0 (G) = C0 (G) set of continuous functions on G with compact support
(l ∈ N) set of l-times continuously di¬erentiable functions on G
l
C0 (G)
with compact support

C0 (G) set of in¬nitely often continuously di¬erentiable functions on G
with compact support
(p ∈ [1, ∞)) set of Lebesgue-measurable functions whose pth power
p
L (G)
of their absolute value is Lebesgue-integrable on G

L (G) set of measurable, essentially bounded functions
scalar product in L2 (G) †
·, · 0,G
norm in L2 (G) †
· 0,G
· 0,p,G (p ∈ [1, ∞]) norm in Lp (G) †
· ∞,G norm in L∞ (G) †
(l ∈ N, p ∈ [1, ∞]) set of l-times weakly di¬erentiable functions
l
Wp (G)
from Lp (G), with derivatives in Lp (G)
· l,p,G (l ∈ N, p ∈ [1, ∞]) norm in Wp (G) †
l

(l ∈ N, p ∈ [1, ∞]) seminorm in Wp (G) †
| · |l,p,G l
A.2. Basic Concepts of Analysis 393

:= W2 (G) (l ∈ N)
H l (G) l

(l ∈ N) scalar product in H l (G) †
·, · l,G
(l ∈ N) norm in H l (G) †
· l,G
(l ∈ N) seminorm in H l (G) †
| · |l,G
·, · 0,h discrete L2 („¦)-scalar product
· 0,h discrete L2 („¦)-norm
L2 (‚G) set of square Lebesgue-integrable functions on the boundary ‚G
1
set of functions from H 1 (G) with vanishing trace on ‚G
H0 (G)
C([0, T ], X) = C 0 ([0, T ], X) set of continuous functions on [0, T ] with values
in the normed space X
C ([0, T ], X) (l ∈ N) set of l-times continuously di¬erentiable functions on
l

[0, T ] with values in the normed space X
Lp ((0, T ), X) (p ∈ [1, ∞]) Lebesgue-space of functions on [0, T ] with values
in the normed space X

Convention: In the case G = „¦, this speci¬cation is omitted.



A.2 Basic Concepts of Analysis
A subset G ‚ Rd is called a set of measure zero if, for any number µ > 0, a
countable family of balls Bj with d-dimensional volume µj > 0 exists such
that
∞ ∞
G‚
µj < µ and Bj .
j=1 j=1

Two functions f, g : G ’ R are called equal almost everywhere (in short:
equal a.e., notation: f ≡ g) if the set {x ∈ G : f (x) = g(x)} is of measure
zero.
In particular, a function f : G ’ R is called vanishing almost everywhere
if it is equal to the constant function zero almost everywhere.
A function f : G ’ R is called measurable if there exists a sequence
(fi )i of step functions fi : G ’ R such that fi ’ f for i ’ ∞ almost
everywhere.
In what follows, G denotes a subset of Rd , d ∈ N.
(i) A point x = (x1 , x2 , . . . , xd )T ∈ Rd is called a boundary point of G
if every open neighbourhood (perhaps an open ball) of x contains a
point of G as well as a point of the complementary set R \ G.
(ii) The collection of all boundary points of G is called the boundary of
G and is denoted by ‚G.
(iii) The set G := G ∪ ‚G is called the closure of G.
(iv) The set G is called closed if G = G.
394 A. Appendices

(v) The set G is called open if G © ‚G = ….
(vi) The set G \ ‚G is called the interior of G and is denoted by int G.
A subset G ‚ Rd is called connected if for arbitrary distinct points
x1 , x2 ∈ G there exists a continuous curve in G connecting them.
The set G is called convex if any two points from G can be connected by
a straight-line segment in G.
A nonempty, open, and connected set G ‚ Rd is called a domain in Rd .
By ± = (±1 , . . . , ±d )T ∈ Nd a so-called multi-index is denoted. Multi-
0
indices are a popular tool to abbreviate some elaborate notation. For
example,
d d d
|±| :=
±
±
‚i i
‚ := , ±! := ±i ! , ±i .
i=1 i=1 i=1

The number |±| is called the order (or length) of the multi-index ±.
For a continuous function • : G ’ R, the set supp • := {x ∈ G : •(x) = 0}
denotes the support of •.


A.3 Basic Concepts of Linear Algebra
A square matrix A ∈ Rn,n with entries aij is called symmetric if aij = aji
holds for all i, j ∈ {1, . . . , n}.
A matrix A ∈ Rn,n is called positive de¬nite if x · Ax > 0 for all x ∈
Rn \ {0}.
Given a polynomial p ∈ Pk , k ∈ N0 , of the form
k
with aj ∈ C, j ∈ {0, . . . , k}
aj z j
p(z) =
j=0

and a matrix A ∈ Cn,n , then the following matrix polynomial of A can be
established:
k
aj Aj .
p(A) :=
j=0


Eigenvalues and Eigenvectors
Let A ∈ Cn,n . A number » ∈ C is called an eigenvalue of A if
det(A ’ »I) = 0 .
If » is an eigenvalue of A, then any vector x ∈ Cn \ {0} such that
(” (A ’ »I)x = 0)
Ax = »x
is called an eigenvector of A associated with the eigenvalue ».
A.3. Linear Algebra 395

The polynomial pA (») := det(A ’ »I) is called the characteristic
polynomial of A.
The set of all eigenvalues of a matrix A is called the spectrum of A,
denoted by σ(A).
If all eigenvalues of a matrix A are real, then the numbers »max (A) and
»min (A) denote the largest, respectively smallest, of these eigenvalues.
The number (A) = max»∈σ(A) |»| is called the spectral radius of A.

Norms of Vectors and Matrices
The norm of a vector x ∈ Rn , n ∈ N, is a real-valued function x ’ |x|
satisfying the following three properties:
(i) |x| ≥ 0 for all x ∈ Rn , |x| = 0 ” x = 0 ,
(ii) |±x| = |±| |x| for all ± ∈ R , x ∈ Rn ,
(iii) |x + y| ¤ |x| + |y| for all x, y ∈ Rn .
For example, the most frequently used vector norms are
(a) the maximum norm:
|x|∞ := max |xj | . (A3.1)
j=1...n

p ∈ [1, ∞):
(b) the p -norm,

1/p
n
|x|p := |xj |p . (A3.2)
j=1

The important case p = 2 yields the so-called Euclidean norm:
1/2
n
|x|2 := x2 . (A3.3)
j
j=1

The three most important norms (that is, p = 1, 2, ∞) in Rn are equivalent
in the following sense: The inequalities

1
√ |x|2 ¤ |x|∞ ¤ |x|2 ¤ n |x|∞ ,
n
1
|x|1 ¤ |x|∞ ¤ |x|1 ¤ n |x|∞ ,
n

1
√ |x|1 ¤ |x|2 ¤ |x|1 ¤ n |x|2
n
are valid for all x ∈ Rn .
The norm of the matrix A ∈ Rn,n is a real-valued function A ’ A
satisfying the following four properties:
A ≥ 0 for all A ∈ Rn,n , A = 0 ” A = 0,
(i)
±A = |±| A for all ± ∈ R , A ∈ Rn,n ,
(ii)
A+B ¤ A + B for all A, B ∈ Rn,n ,
(iii)
AB ¤ A B for all A, B ∈ Rn,n .
(iv)
396 A. Appendices

In comparison with the de¬nition of a vector norm, we include here an
additional property (iv), which is called the submultiplicative property. It
restricts the general set of matrix norms to the practically important class
of submultiplicative norms.
The most common matrix norms are
(a) the total norm:
:= n max |aik | ,
A (A3.4)
G
1¤i,k¤n

(b) the Frobenius norm:
1/2
n
a2
A := , (A3.5)
F ik
i,k=1

(c) the maximum row sum:
n
|aik | ,
A := max (A3.6)

1¤i¤n
k=1

(d) the maximum column sum:
n
|aik | .
A := max (A3.7)
1
1¤k¤n
i=1

All these matrix norms are equivalent. For example, we have
1
¤A ¤A ¤n A p ∈ {1, ∞} ,
A ,
G p G p
n
or
1
¤A ¤A ¤n A
A .
G F G F
n

<<

. 15
( 16)



>>