Institut f¨ ur Numerische Mathematik
Kompaktkurs
Lineare Gleichungssysteme Hierarchische Matrizen M. Bebendorf, O. Steinbach
O. Steinbach
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 1 / 16
Institut f¨ ur Numerische Mathematik
Numerische Simulation I station¨ are und instation¨are partielle Differentialgleichungen I I I
I
Kopplung verschiedener physikalischer Felder I I I
I
Fluid–Struktur–Akustik elektromechanische Felder Mehrk¨ orpersysteme
Diskretisierung I I I I
I
Potentialgleichung (Temperatur, elektromagnetische Felder,. . . ) Festk¨ orpermechanik (lineare Elastostatik, Elastoplastizit¨ at,. . . ) Str¨ omungsmechanik (Stokes, Navier–Stokes)
Finite Element Methoden Finite Differenzen Methoden Randelementmethoden Zeitdiskretisierung
Familie von linearen (nichtlinearen) Gleichungssystemen Ax = f ,
O. Steinbach
A ∈ Rn×n , n → ∞
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 2 / 16
Institut f¨ ur Numerische Mathematik
Lineares Gleichungssystem Ax = f ,
I
A ∈ Rn×n ,
n→∞
Beschreibung und Anwendung der Systemmatrix A I
I
I
I
O. Steinbach
A heißt vollbesetzt, falls Anzahl der Nichtnullelemente von A von der Gr¨ oßenordnung O(n2 ) ist. A heißt schwachbesetzt, falls Anzahl der Nichtnullelemente von A von der Gr¨ oßenordnung O(n logα n) ist. A heißt data sparse, falls A n 2 Nichtnullelemente besitzt, diese aber durch O(n logα n) Eintr¨ age beschrieben werden kann. 0 1 1 1 0 0 a1 a2 an 1 ... n 1 B an a1 a2 C B B . C B . C .. C A = B C , A = @ .. A (1 . . . n) = @ .. .. . A @ A . n n . . . n2 a2 an a1 zirkulante Matrix (n Eintr¨ age) Struktur der Matrix A
Rang 1 Matrix (2n Eintr¨ age)
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 3 / 16
Institut f¨ ur Numerische Mathematik
Lineares Gleichungssystem Ax = f ,
A ∈ Rn×n ,
n→∞
L¨osungsverfahren I direkte Verfahren I I I I
I
Frage: Wird tats¨achlich die exakte L¨osung x = A−1 f ben¨otigt? I I I
I
Gauß Elimination mit O(n 3 ) Multiplikationen LU Zerlegung (Cholesky) mit O(n 3 ) Multiplikationen schnelle direkte Methoden f¨ ur schwachbesetzte Matrizen strukurierte Matrizen (zirkulante Matrizen, FFT) Modellierungsfehler Diskretisierungsfehler Verfahrensfehler (numerische Integration, Rundungsfehler)
Hinreichend genaue L¨osung ist ausreichend! Iterationsverfahren I I I
klassische Iterationsverfahren (Jacobi, Gauß–Seidel, SOR) Gradientenverfahren (steilster Abstieg, minimaler Defekt) Verfahren orthogonaler Richtungen (CG, GMRES, BiCGStab)
Ziel: Effiziente L¨osung von Au = f f¨ ur n → ∞. O. Steinbach
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 4 / 16
Institut f¨ ur Numerische Mathematik
Beispiel: Bestimmung einer st¨ uckweise linearen Approximation n P @ uh (x) = uk ϕk (x) @ k=0 @ ϕk (x) @ @ @ @ @ Minimierungsproblem xk #2 Z1 " n X ∂ uk ϕk (x) dx → min , u(x) − F (u) = F (u) = 0 u0 ,...,un ∂u` k=0
0
Variationsproblem n X k=0
uk
Z1
ϕk (x)ϕ` (x)dx =
Z1
u(x)ϕ` (x)dx
f¨ ur ` = 0, . . . , n
0
0
Lineares Gleichungssystem Mh u = f ,
Mh [`, k] =
Z1
ϕk (x)ϕ` (x)dx
0
O. Steinbach
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 5 / 16
Institut f¨ ur Numerische Mathematik
Systemmatrix
2 1 h Mh = 6 I I I I I
1 4 1
1 .. .
..
.
..
..
.
.
1
∈ R(n+1)×(n+1) 1 4 1 1 2
Mh ist symmetrisch Mh ist streng diagonal dominant (nur in einer Raumdimension) Mh ist tridiagonal (nur in einer Raumdimension) Mh ist positiv definit (gilt f¨ ur beliebige Ansatzfunktionen) Mh ist schwach besetzt (3n + 1 Nichtnulleintr¨age)
O. Steinbach
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 6 / 16
Institut f¨ ur Numerische Mathematik
Beispiel f¨ ur n = 7 2 B 1 B B 1 B B Mh = B 42 B B B @ 0
0
Mh−1
B B B 14 B B = B 2911 B B B @
5042 −1351 362 −97 26 −7 2 −1
−1351 2702 −724 194 −52 14 −4 2
1 4 1
362 −724 2534 −679 182 −49 14 −7
1 4 1
1
1 4 1
−97 194 −679 2522 −676 182 −52 26
1 4 1
1 4 1 26 −52 182 −676 2522 −679 194 −97
1 4 1
C C C C C C C C C 1 A 2 −7 14 −49 182 −679 2534 −724 362
2 −4 14 −52 194 −724 2702 −1351
−1 2 −7 26 −97 362 −1351 5042
1 C C C C C C C C C A
Matrix Mh schwach besetzt, aber inverse Matrix Mh−1 vollbesetzt Frage: Gibt es eine Struktur in der Darstellung von Mh−1 ? O. Steinbach
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 7 / 16
Institut f¨ ur Numerische Mathematik
0
5042 −1351 B B B B„ « B ` −362 B B 97 B 14 B −1 B M = B h 2911 B B 0 B 26 B B B −7 B B B @ 2 @ −1 0
„
−1351 2702 2
−1
1
C` C A
1
−1 2
´
−2
1
1
2 B B 1 B B B B B B 1 B B Mh = B 42 B B B B B B B B @
4 1
1 4
1
1
4 1
«
`
−362
97
1 4
1
1
4
1
1
4
1
1
2
´
0
C C C C C C C C C C C C C C C C C C C A
1 1 ` −2 C C 26 A 7 −26
1
C ´ C C −7 2 −1 C C 2534 −679 C C −679 2522 C C C C „ « C ` ´C 2522 −679 97 2 −1 C −679 2534 −362 C C ´ C 7 −26 C „ « ` ´ 2 2702 −1351 A 97 −362 −1 −1351 5042 B B @
Mh−1 ist exakt als Hierarchische Matrix darstellbar! O. Steinbach
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 8 / 16
Institut f¨ ur Numerische Mathematik
Beispiel: Finite Element Methoden Dirichlet Randwertproblem f¨ ur Poisson Gleichung −∆u(x) = −
d X ∂2 u(x) = f (x) ∂xi2
f¨ ur x ∈ Ω ⊂ Rd ,
u(x) = 0
f¨ ur x ∈ ∂Ω
i=1
Variationsproblem Z Z ∇u(x)∇v (x)dx = f (x)v (x)dx, Ω
u(x), v (x) = 0 f¨ ur x ∈ ∂Ω
Ω
Diskretisierung mit st¨ uckweise linearen finiten Elementen Z Kh [`, k] = ∇ϕk (x)∇ϕ` (x)dx Kh u = f , Ω
mit FEM Steifigkeitsmatrix Kh I Kh ist symmetrisch und positiv definit I Kh ist schwach besetzt O. Steinbach
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 9 / 16
Institut f¨ ur Numerische Mathematik
Beispiel: d = 1 0
2 B −1 B B 1B B Kh = B hB B B @
−1 2 −1
1
−1 .. . .. .
..
C C C C C C ∈ R(n−1)×(n−1) C C C −1 A 2
.
..
. −1
−1 2 −1
n=9 0
−1
Kh
B B B 1B B = B 9B B B @
„
« 8 7 „ 7 « 14 ` ´ 6 1 2 5 1 0 4 B 3 C` @ 2 A 1 1
„
2
«
1 „2
18 15
3
`
6
4
´
15 20
5 «
´
1 1 ´ B 2 C` @ 3 A 4 3 2 1 4 « « „ „ ` ´ 5 20 15 2 1 6 15 18 „ « „ « ` ´ 2 14 7 5 6 1 7 8 0
1 C C C C C C C C C A
Kh−1 ist exakt als Hierarchische Matrix darstellbar!
O. Steinbach
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 10 / 16
Institut f¨ ur Numerische Mathematik
Neumann Randwertproblem f¨ ur Laplace Gleichung −∆u(x) = 0
∂ u(x) = g (x) ∂n
f¨ ur x ∈ Ω,
f¨ ur x ∈ ∂Ω
u = 1 ist L¨osung des homogenen Neumann Randwertproblems −∆u(x) = 0
f¨ ur x ∈ Ω,
∂ u(x) = 0 ∂n
f¨ ur x ∈ ∂Ω
Randwertproblem ist nicht eindeutig l¨osbar! Variationsproblem Z Z ∇u(x)∇v (x)dx = g (x)v (x)dsx Ω
∂Ω
L¨osbarkeitsbedingung v =1:
Z
g (x)dsx = 0
∂Ω
O. Steinbach
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 11 / 16
Institut f¨ ur Numerische Mathematik
FEM Diskretisierung (d = 1) Kh u = f ,
Z
Kh [`, k] =
∇ϕk (x)∇ϕ` (x)dx
Ω
mit FEM Steifgkeitsmatrix 0
B B B 1B B Kh = B hB B B @
−1
2 −1
−1 2 −1
−1 .. . .. .
1 ..
.
..
. −1
I
Kh ist symmetrisch
I
Kh ist schwachbesetzt Kh ist singul¨ar (Eigenvektor u = 1)
I
O. Steinbach
−1 2 −1
Lineare Gleichungssysteme
−1 2
−1
C C C C C C ∈ R(n+1)×(n+1) C C C A
SIMNET Kurs 24.–27.4.2006 12 / 16
Institut f¨ ur Numerische Mathematik
Erweitertes Variationsproblem Z Z Z Z ∇u(x)∇v (x)dx + u(x)dx v (x)dx = g (x)v (x)dsx Ω
Ω
Ω
v = 1 liefert Skalierungsbedingung Z
∂Ω
u(x)dx = 0
Ω
modifzierte Steifigkeitsmatrix e h = Kh + a a > , K I I
ak =
eh ist symmetrisch und positiv definit K eh ist data sparse K
O. Steinbach
Z
ϕk (x)dx
Ω
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 13 / 16
Institut f¨ ur Numerische Mathematik
Beispiel: Stokes System −∆u(x)+∇p(x) = f (x),
div u(x) = 0
f¨ ur x ∈ Ω,
u(x) = 0
f¨ ur x ∈ ∂Ω
mit Geschwindigkeitsfeld u und Druck p. Bemerkung: Druck ist nur eindeutig bis auf additive Konstante! Erweitertes Variationsproblem Z Z ∇u(x)∇v (x)dx + ∇p(x)v (x)dx Z Ω
O. Steinbach
Ω
div u(x)q(x)dx +
Z Ω
Ω
p(x)dx
Z
q(x)dx
=
Z
f (x)v (x)dx
Ω
=
0
Ω
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 14 / 16
Institut f¨ ur Numerische Mathematik
Erweitertes Variationsproblem Z Z ∇u(x)∇v (x)dx + ∇p(x)v (x)dx Z
Ω
div u(x)q(x)dx +
Ω
Z
Ω
p(x)dx
Ω
Z
=
Z
f (x)v (x)dx
Ω
q(x)dx
=
0
!
f 0
!
Ω
FEM Diskretisierung Ah Bh> I I I
−Bh a a>
!
u p
=
Voraussetzung: diskrete Stabilit¨atsbedingung Systemmatrix positiv definit, aber Block schief symmetrisch Erweiterung auf nichtlineares Navier–Stokes System
O. Steinbach
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 15 / 16
Institut f¨ ur Numerische Mathematik
Literatur 1. R. Barrett et. al.: Templates for the Solution of Linear Systems. Building Blocks for Iterative Methods. SIAM, Philadelphia, 1993. 2. A. Meister: Numerik linearer Gleichungssysteme. Eine Einf¨ uhrung in moderne Verfahren. Vieweg, Braunschweig, 1999. 3. Y. Saad: Iterative Methods for Sparse Linear Systems. PWS Publishing, 1995. 4. O. Steinbach: L¨osungsverfahren f¨ ur lineare Gleichungssysteme. Algorithmen und Anwendungen. B. G. Teubner, Stuttgart, Leipzig, Wiesbaden, 2005.
O. Steinbach
Lineare Gleichungssysteme
SIMNET Kurs 24.–27.4.2006 16 / 16