skript-mehrgitter-sander.tex 119 KB
Newer Older
1
% The following line selects the page layout for printing
2
\documentclass{scrbook}
3
4
5

% The following lines select a page layout without vertical margins,
% and a very large right margin.  This is much better for online teaching.
6
7
8
9
% \documentclass[paper=310mm:297mm]{scrbook}
% \KOMAoptions{twoside=false}
%
% \usepackage[width=12cm, left=3cm, top=0.2cm, bottom=0.5cm]{geometry}
10
11
12
13

\usepackage[utf8]{inputenc}

\usepackage{longtable,tabularx,tabulary,booktabs}
14
\usepackage{aligned-overset}
15
\usepackage{array}
16
\usepackage[linesnumbered,inoutnumbered]{algorithm2e}
17
18
\usepackage{caption}
\usepackage{url}
Sander, Oliver's avatar
Sander, Oliver committed
19
\usepackage{overpic}
20
21
\usepackage{graphicx}
\usepackage{tikz}
22
\usepackage{pgfplots}
23
\pgfplotsset{compat=1.17}
24
25
26
27

\usepackage{amsfonts}
\usepackage{lmodern}
\usepackage[ngerman]{babel}
28
\usepackage{csquotes}
29
30
31
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{amsthm}
32
\usepackage{bm}
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
\usepackage{mathtools}
\usepackage{colonequals}
\usepackage{enumitem}

\usepackage[defernumbers=true,
            maxbibnames=99,  % Never write 'et al.' in a bibliography
            giveninits=true, % Abbreviate first names
            sortcites=true,  % Sort citations when there are more than one
            natbib=true,     % Allow commands like \citet and \citep
            backend=biber]{biblatex}
\addbibresource{skript-mehrgitter-sander.bib}

\usepackage{microtype}
\usepackage{todonotes}
\usepackage{esdiff}
\usepackage{hyperref}

\allowdisplaybreaks[0]
\setlength{\parindent}{0pt}  % No indent of first line of paragraph

\newcommand{\R}{\mathbb R}
\newcommand{\N}{\mathbb N}
\newcommand{\curl}{\operatorname{curl}}
\renewcommand{\div}{\operatorname{div}}
57
\newcommand{\grad}{\operatorname{grad}}
58
59
\newcommand{\tr}{\operatorname{tr}}
\DeclareMathOperator*{\argmin}{arg\,min}
Sander, Oliver's avatar
Sander, Oliver committed
60
\newcommand{\Gammatight}[1]{{\Gamma\hspace{-0.8mm}_{#1}}}  % A \Gamma with a subscript, and extra kerning
61
62
63

% Bold letters
\newcommand{\bn}{\mathbf{n}}
64
\newcommand*\bu{\bm{u}}
65
66
\newcommand{\bv}{\mathbf{v}}

67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
\newcommand*\bA{\bm{A}}
\newcommand*\bH{\bm{H}}
\newcommand*\bI{\bm{I}}
\newcommand*\bL{\bm{L}}
\newcommand*\bP{\bm{P}}
\newcommand*\bQ{\bm{Q}}
\newcommand*\bV{\bm{V}}
\newcommand*\bHd{\bH(\div)}
\newcommand*\til{\widetilde}
\newcommand{\stackrela}[2]{\overset{#1}&{#2}}

\newcommand{\X}{\mathcal{X}}				 % Triangulierung
\newcommand{\V}{\mathcal{V}}				 % Triangulierung
\newcommand{\E}{\mathcal{E}}				 % Triangulierung
\newcommand{\F}{\mathcal{F}}				 % Triangulierung
\newcommand{\T}{\mathcal{T}}				 % Triangulierung

\DeclarePairedDelimiter{\set}{\lbrace}{\rbrace}
\DeclarePairedDelimiter{\norma}{\Vert}{\Vert_{A}}
86
87
\DeclarePairedDelimiter{\abs}{\lvert}{\rvert}
\DeclarePairedDelimiter{\norm}{\lVert}{\rVert}
88
89
90

\theoremstyle{plain} %Text ist Kursiv
\newtheorem{theorem}{Satz}[chapter]
91
\newtheorem{hilfssatz}[theorem]{Hilfssatz}
92
93
94
95
96
97
98
99
100
\newtheorem{lemma}[theorem]{Lemma}
\newtheorem{proposition}[theorem]{Proposition}
\newtheorem{corollary}[theorem]{Korollar}
\newtheorem{definition}[theorem]{Definition}
\newtheorem{problem}[theorem]{Problem}

\theoremstyle{definition} %Text ist \"upright"
\newtheorem{remark}[theorem]{Bemerkung}
\newtheorem{example}[theorem]{Beispiel}
101
\newtheorem{exercise}[theorem]{Übung}
102
103
104
105
106
107
108
109
110
111
112
113
114
115

\swapnumbers

% Suche auch in diesem Verzeichnis nach Bildern
\graphicspath{{gfx/}}

\begin{document}

\begin{titlepage}

\begin{center}

\vspace*{0.3\textheight}

116
{\Huge Mehrgitterverfahren}
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
\vspace*{0.05\textheight}

{\Large Oliver Sander}
\vspace*{0.1\textheight}

\today
\vspace{0.4\textheight}

\end{center}

\end{titlepage}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%  Copyright notice
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\thispagestyle{empty}

\vspace*{\fill}
\noindent
\includegraphics[width=0.2\textwidth]{license-logo-cc-by-sa}
Oliver Sander, 2021

\bigskip

\noindent
Copyright 2021 by Oliver Sander.
This work is licensed under the Creative Commons Attribution-ShareAlike 4.0 International License.
To view a copy of this license, visit \url{http://creativecommons.org/licenses/by-sa/4.0/}.

\pagebreak

\tableofcontents
\thispagestyle{empty}
\setcounter{page}{0}

Sander, Oliver's avatar
Sander, Oliver committed
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
\chapter{Finite Elemente Methoden}

In diesem Kapitel geben wir eine kurze Einführung in die Methode der Finiten Elemente.
Damit legen wir die Grundlagen, und motivieren die weiteren Auführungen.

\section{Schwache Formulierung}

\begin{figure}
 \begin{center}
  \begin{overpic}[width=0.3\textwidth]{gfx/gebiet}
   \put(50,82){$\Gamma_D$}
   \put(25,2){$\Gamma_N$}
   \put(30,40){$\Omega$}
  \end{overpic}
 \end{center}
\caption{Einfaches Gebiet $\Omega$ mit Dirichletrand $\Gamma_D$ und Neumannrand
$\Gamma_N$.}
\label{fig:gebiet}
\end{figure}

Sei $\Omega \subset \R^d$ ein Gebiet und $f : \Omega \to \R$.
Der Rand $\partial \Omega$ von $\Omega$ sei zerlegt in zwei disjunkte Teile $\Gammatight{D}$ und $\Gammatight{N}$,
und wir bezeichnen die äußere Einheitsnormale mit $\nu$.
175
Wir betrachten als Beispiel die \emph{Poisson-Gleichung}
Sander, Oliver's avatar
Sander, Oliver committed
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
\begin{equation}
\label{eq:starke_formulierung}
 - \Delta u \colonequals - \sum_{i=1}^d \frac{\partial^2 u}{\partial x_i^2} = f
\qquad \text{auf $\Omega$},
\end{equation}
mit den Randbedingungen
\begin{align}
\nonumber
 u & = 0 \qquad \text{auf $\Gammatight{D}$}, \\
 \label{eq:neumann_randbedingungen}
\frac{\partial u}{\partial \nu} & = g \qquad \text{auf $\Gammatight{N}$}.
\end{align}
Wir bezeichnen mit $H^1(\Omega)$ den Raum aller skalarer $L^2$-Funktionen auf
$\Omega$ mit schwacher Ableitung in $L^2(\Omega)$.
Als Raum von Testfunktionen führen wir
\begin{equation*}
 H^1_{\Gammatight{D}} (\Omega)
=
\{ v \in H^1(\Omega) \; | \; v|_{\Gammatight{D}} = 0 \}
\end{equation*}
ein, wobei $v|_{\Gammatight{D}}$ im Sinne des Spursatzes zu verstehen ist.
Zur Abkürzung setzen wir $H \colonequals H^1_{\Gammatight{D}} (\Omega)$.

Multiplikation von \eqref{eq:starke_formulierung} mit einer Testfunktion
$v \in H^1_{\Gammatight{D}} (\Omega)$ und Integration über $\Omega$ ergibt
\begin{equation*}
 - \int_\Omega \Delta u \cdot v \, dx = \int_\Omega f v \, dx.
\end{equation*}
Wir wenden darauf die Greensche Formel an, und erhalten
\begin{equation*}
 \int_\Omega \nabla u \nabla v \, dx
=
\int_{\partial \Omega} \frac{\partial u}{\partial \nu} v \, ds
+
\int_\Omega f v \, dx.
\end{equation*}
Nach Voraussetzung besteht $\partial \Omega$ aus den zwei Teilen $\Gammatight{D}$ und $\Gammatight{N}$.
Die Testfunktion $v$ wurde so konstruiert, dass $v|_{\Gammatight{D}} = 0$.  Mithin bleibt
nur ein Integral über $\Gammatight{N}$ übrig.  Dort können wir die Neumannrandbedingungen~\eqref{eq:neumann_randbedingungen}
einsetzen und erhalten
\begin{equation*}
 \int_\Omega \nabla u \nabla v \, dx
=
\int_{\Gamma_N} g v \, ds
+
\int_\Omega f v \, dx.
\end{equation*}
Wir schreiben dies als
\begin{equation}
\label{eq:schwache_formulierung}
226
 a(u,v) = \ell(v)
Sander, Oliver's avatar
Sander, Oliver committed
227
228
229
\end{equation}
mit
\begin{equation*}
230
231
232
 a(v,w) \colonequals \int_\Omega \nabla v \nabla w \, dx
\qquad \text{und} \qquad
 \ell(v) \colonequals \int_{\Gamma_N} g v \, ds
Sander, Oliver's avatar
Sander, Oliver committed
233
234
235
+
\int_\Omega f v \, dx.
\end{equation*}
236
237
Man nennt \eqref{eq:schwache_formulierung} die \emph{schwache Formulierung}
(oder \emph{Variationsformulierung})
Sander, Oliver's avatar
Sander, Oliver committed
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
von \eqref{eq:starke_formulierung}.

$a(\cdot, \cdot)$ ist eine symmetrische Bilinearform.  Hat $\Gammatight{D}$
positives $d-1$-dimensionales Maß, so ist $a(\cdot,\cdot)$
$H^1_{\Gamma_D}$-elliptisch (siehe z.B.\ \cite{braess:2013}), d.h.\ es existieren $\gamma, \Gamma > 0$
so daß
\begin{equation*}
 \gamma \norm{v}^2_H \le a(v,v)
\qquad
\abs{a(v,w)} \le \Gamma \norm{v}_H \norm{w}_H
\qquad
\text{für alle $v,w \in H^1_{\Gamma_D}(\Omega)$}.
\end{equation*}


\begin{theorem}[Lax-Milgram-Lemma \cite{braess:2013}]
Sei $a(\cdot,\cdot)$ $H$-elliptisch (nicht notwendig symmetrisch!).
Dann hat die Variationsgleichung
\begin{equation*}
 u \in H
\quad : \quad
259
a(u,v) = \ell(v)
Sander, Oliver's avatar
Sander, Oliver committed
260
261
\qquad v \in H
\end{equation*}
262
for jedes $\ell \in H'$ eine eindeutig bestimmte Lösung.
Sander, Oliver's avatar
Sander, Oliver committed
263
264
265
266
267
268
269
270
271
\end{theorem}




\section{Diskretisierung (Galerkin-Verfahren)}

Sei $V_h$ ein endlichdimensionaler Teilraum von $H^1_{\Gammatight{D}}$ (mit Dimension $n$).
Der Index $h$
272
bezeichne einen Ordnungsparameter.  Die \emph{diskrete Variationsformulierung}
Sander, Oliver's avatar
Sander, Oliver committed
273
274
275
lautet
\begin{equation}
\label{eq:diskrete_formulierung}
276
 u_h \in V_h \quad : \quad a(u_h,v_h) = \ell(v_h)
Sander, Oliver's avatar
Sander, Oliver committed
277
278
279
280
281
282
283
284
285
\qquad \text{für alle $v_h \in V_h$}.
\end{equation}
Anwendung des Lax--Milgram-Lemmas ergibt wieder Existenz und Eindeutigkeit
von Lösungen.

Sei $\{\varphi_i\}, 0 \le i < n$, eine Basis von $V_h$. Dann ist
\eqref{eq:diskrete_formulierung} äquivalent zu
\begin{equation*}
%\label{eq:diskrete_formulierung}
286
 u_h \in V_h \quad : \quad a(u_h, \varphi_i) = \ell(\varphi_i)
Sander, Oliver's avatar
Sander, Oliver committed
287
288
289
\qquad \text{für alle $0 \le i < n$}.
\end{equation*}
Dies wiederum entspricht dem linearen Gleichungssystem
290
291
\begin{equation}
\label{eq:fe_lineares_gleichungssystem}
292
 Ax = b,
293
\end{equation}
Sander, Oliver's avatar
Sander, Oliver committed
294
295
296
wobei
\begin{alignat*}{2}
A &\in \R^{n \times n},
297
\qquad &A_{ij} & \colonequals a(\varphi_i, \varphi_j)
Sander, Oliver's avatar
Sander, Oliver committed
298
299
300
301
302
=
\int_\Omega \nabla \varphi_i \nabla \varphi_j \, dx \\
%
b &\in \R^n,
\qquad &
303
b_i & \colonequals \ell(\varphi_i)
Sander, Oliver's avatar
Sander, Oliver committed
304
305
306
307
308
=
\int_{\Gamma_N} g \varphi_i \, ds
+
\int_\Omega f \varphi_i \, dx
\end{alignat*}
309
und $x \in \R^n$ die Koeffizienten von $u_h$ bzgl.\ der Basis $\{\varphi_i\}$ sind.
Sander, Oliver's avatar
Sander, Oliver committed
310

311
312
313
\bigskip

Die Matrix $A$ nennt man \emph{Steifigkeitsmatrix}.  Dieser Ausdruck kommt aus den
Sander, Oliver's avatar
Sander, Oliver committed
314
315
316
Anfangstagen der Finiten Elemente in der Ingenieurswelt, wo insbesondere mechanische
Probleme gelöst wurden.  Die Matrix $A$ beschreibt dabei die Steifigkeit des simulierten
Objekts, und der Vektor $b$ beschreibt die von außen wirkenden Kräfte.  Den Vektor $b$
317
nennt man deshalb manchmal auch \emph{Lastvektor}.
Sander, Oliver's avatar
Sander, Oliver committed
318
319
320
321
322
323
324
325


\section{Finite Elemente}
\label{sec:finite_elemente}

Finite Elemente sind eine bestimmte Möglichkeit, den Teilraum
$V_h$ zu wählen.

326
327
328
329
\medskip

Sei $\Omega$ von jetzt ab \emph{polygonal berandet}.  Wir füllen
$\Omega$ mit einer Triangulierung / einem \emph{Gitter}.  Ein Beispielgitter für
Sander, Oliver's avatar
Sander, Oliver committed
330
331
332
333
334
335
336
337
338
339
340
341
342
343
ein zweidimensionales Gebiet sieht man in Abbildung~\ref{fig:2d-Beispielgitter}.

\begin{figure}
\begin{center}
\includegraphics[width=0.5\textwidth]{dreiecksgitter}
\end{center}
\caption{Ein einfaches zweidimensionales Gitter}
\label{fig:2d-Beispielgitter}
\end{figure}

\begin{definition}
\label{def:triangulierung}
 Es sei $\Omega \subset \R^2$ polygonal berandet, und $\Gammatight{D}$ auf dem
Rand aufgelöst.  Eine Menge $\mathcal{T}$
344
von abgeschlossenen Dreiecken $T$ heißt \emph{Triangulierung} von $\Omega$, falls gilt
Sander, Oliver's avatar
Sander, Oliver committed
345
346
347
\begin{enumerate}
 \item Die Menge der Dreiecke überdeckt den Abschluss des Gebiets
\begin{equation*}
348
 \overline{\Omega} = \bigcup_{T \in \mathcal{T}} T.
Sander, Oliver's avatar
Sander, Oliver committed
349
350
351
352
353
354
355
356
\end{equation*}
%
\item Der Schnitt zweier Dreiecke aus $\mathcal T$ ist entweder eine
gemeinsame Kante, ein gemeinsamer Eckpunkt, oder leer.
\end{enumerate}
\end{definition}
Der Diskretisierungsparameter
\begin{equation*}
357
 h \colonequals \max_{T \in \mathcal{T}} \operatorname{diam} T,
Sander, Oliver's avatar
Sander, Oliver committed
358
359
\end{equation*}
der als Index in $V_h$ auftaucht, beschreibt traditionell die Feinheit des Gitters.
360
Oft nennt man $h$ auch \emph{Gitterweite}.
Sander, Oliver's avatar
Sander, Oliver committed
361
362
363
364

\begin{definition}
Der Raum
\begin{equation*}
365
 V_h^{p}
Sander, Oliver's avatar
Sander, Oliver committed
366
367
368
\colonequals
\{ v \in C(\overline{\Omega})
\; | \;
369
\text{$v$ ist Polynom von Grad $\le p$ auf $T$ für alle $T \in \mathcal{T}$, und $v|_{\Gamma_D} = 0$} \}
Sander, Oliver's avatar
Sander, Oliver committed
370
\end{equation*}
371
heißt Raum der Lagrange-Elemente $p$-Ordnung
Sander, Oliver's avatar
Sander, Oliver committed
372
373
374
375
376
377
378
379
bezüglich der Triangulierung $\mathcal{T}$.
\end{definition}

\begin{figure}
\begin{center}
\includegraphics[width=0.4\textwidth]{generalP1function}
\includegraphics[width=0.4\textwidth]{nodalP1basis}
\end{center}
380
\caption{Links: eine Finite-Elemente Funktion.  Rechts: ein Element der Knotenbasis erster Ordnung}
Sander, Oliver's avatar
Sander, Oliver committed
381
382
383
\label{fig:beispiel-fe-funktion}
\end{figure}

384
385
386
Als Basis von $V_h^{p}$ wählt man gewöhnlich die \emph{Knotenbasis}.  Sei $\mathcal L$
die Menge der Lagrange-Knoten des Gitters.  Dann ist die Knotenbasis die Menge
der Funktionen $\varphi_i \in V_h^{p}$ mit
Sander, Oliver's avatar
Sander, Oliver committed
387
388
389
\begin{equation*}
 \varphi_i (v_j) = \delta_{ij}
\qquad
390
\text{für alle $v_j \in \mathcal{L}$}.
Sander, Oliver's avatar
Sander, Oliver committed
391
\end{equation*}
392
Es ist also insbesondere $\dim V_h^{p} = \abs{\mathcal L \setminus \mathcal{L}_{\Gamma_D}}$
Sander, Oliver's avatar
Sander, Oliver committed
393
394
(Abbildung~\ref{fig:beispiel-fe-funktion}).

395
396
397
398
\bigskip

Der Diskretisierungsfehler hängt hauptsächlich von der Gitterauflösung $h$,
der Ordnung $p$ der Finite-Elemente-Funktionen, und der Glattheit von $u$ ab.
Sander, Oliver's avatar
Sander, Oliver committed
399
400

\begin{theorem}
401
402
403
 Es sei $h$ klein genug, $u \in H^1_0(\Omega) \cap H^{p+1}(\Omega)$ die Lösung
 der schwachen Formulierung, und $u_h$ die entsprechende Finite-Elemente-Lösung in $V_h^p$.
 Dann gilt die a priori Abschätzung
Sander, Oliver's avatar
Sander, Oliver committed
404
\begin{equation*}
405
 \norm{u - u_h}_1 \le C h^p \abs{u}_{p+1}
Sander, Oliver's avatar
Sander, Oliver committed
406
407
408
409
410
411
412
413
\end{equation*}
mit einer Konstanten $C$.
\end{theorem}

\begin{proof}
 \citet{braess:2013}
\end{proof}

414
\chapter{Mehrgitterverfahren}
Sander, Oliver's avatar
Sander, Oliver committed
415

416
\section{Iterative Lösungsverfahren für dünnbesetzte lineare Gleichungssysteme}
417

418
419
Wir beschäftigen uns jetzt mit der Frage, wie die linearen
Gleichungssysteme~\eqref{eq:fe_lineares_gleichungssystem} gelöst werden können.
420

421
\subsection{Eigenschaften von Finite-Elemente-Matrizen}
422
423
424

Finite-Elemente-Matrizen haben bestimmte Eigenschaften, die das Lösen schwierig machen.

425
\subsubsection{Größe}
426
Finite-Elemente-Matrizen können sehr groß werden:
427
\begin{itemize}
428
429
430
431
432
433
434
 \item Für jeden Lagrangeknoten enthält~\eqref{eq:fe_lineares_gleichungssystem} eine Gleichung.

 \item Für $d$-dimensionale Gebiete und Lagrange-Elemente erster Ordnung
       hat man etwa $n \approx \operatorname{Vol}( \Omega ) \cdot h^{-d}$ Knoten.

 \item Je feiner das Gitter, desto präziser die Lösung, desto größer aber auch die Matrix.

435
 \item Eine Zahl in doppelter Genauigkeit braucht $8$~Byte.  In ein Gigabyte Hauptspeicher
436
437
438
439
440
441
442
443
444
   passen also ca.\ $1{,}25\cdot 10^{8}$ Zahlen.
\end{itemize}

\begin{exercise}
 Wie hängt die Problemgröße von der Ordnung der Ansatzfunktionen ab?
\end{exercise}

\subsubsection{Kondition}

445
446
Die Matrizen zu elliptischen Problemen sind positiv definit, aber sie
sind schlecht konditioniert.
447

448
\medskip
449

450
451
452
453
454
455
456
457
458
459
460
461
462
\emph{Beispiel:} Die Poisson-Gleichung auf $\Omega = [0,1]^2$, Finite-Elemente
mit Lagrange-Elementen auf einem strukturierten Dreiecksgitter mit Kantenlänge $h = \frac{1}{n}$.
\begin{itemize}
 \item Die Eigenwerte der Steifigkeitsmatrix $A$ sind~\cite[Kapitel~12.3.3]{dahmen_reusken:2008}
  \begin{align*}
   \lambda_{\nu \mu}
   & =
   \frac{4}{h^2} \Big[\sin^2 \Big(\frac{1}{2} \nu \pi h \Big)+\sin^2 \Big(\frac{1}{2} \mu \pi h \Big) \Big],
   \qquad
   1 \le \nu,\mu < n.
  \end{align*}

 \item Die Kondition ist deshalb ($A$ ist symmetrisch)
463
  \begin{equation*}
464
   \kappa(A)
465
466
467
   =
   \frac{\lambda_\text{max}}{\lambda_\text{min}}
   =
468
469
470
471
472
473
474
   \frac{\lambda_{n-1,n-1}}{\lambda_{1,1}}
   =
   \frac{\sin^2 \Big(\frac{1}{2} (n-1) \pi h \Big)}{\sin^2 \Big(\frac{1}{2} \pi h \Big)}
   =
   \frac{\cos^2 \big(\frac{1}{2} \pi h \big)}{\sin^2 \big(\frac{1}{2} \pi h \big)}
   =
   \frac{4}{(\pi h)^2} \big(1+\mathcal{O} (h^2) \big).
475
476
477
478
479
480
481
482
483
484
485
486
487
  \end{equation*}

 \item Die Kondition hängt außerdem ab von:
  \begin{itemize}
   \item der Qualität der Elemente~\cite{shewchuk:2002}

   \item Bei einer Gleichung mit ortsabhängigen Koeffizienten wie
    \begin{equation*}
     - \div (D(x) \grad u) = f
    \end{equation*}
    hängt sie auch von den Eigenwerten von $D$ ab, und wie unterschiedlich
    diese auf $\Omega$ werden.
  \end{itemize}
488
489
\end{itemize}

490
491
\emph{Konsequenz:} Wir brauchen Lösungsverfahren, deren Geschwindigkeit
möglichst unabhängig von $\kappa(A)$ ist.
492

493
494
495
496
497

\begin{exercise}
 Wie hängt die Kondition bei einem dreidimensionalen Gitter von $h$ ab?
\end{exercise}

498
499

\subsubsection{Dünnbesetztheit}
500

501
Mehr zu dünnen Matrizen findet man z.B.\ bei \citet{pissanetzky:1984}.
502

503
\begin{itemize}
504
\item Ein Matrixeintrag $A_{ij} \colonequals a(\varphi_i, \varphi_j)$ ist genau dann
505
506
507
508
509
510
511
512
513
514
515
  nicht Null, wenn sich die Träger von $\varphi_i$ und $\varphi_j$ überlappen.

\item D.h.\ der überwiegende Teil der Matrixeinträge ist Null.

\item \glqq Eine Matrix heißt dünnbesetzt wenn sie so viele Nullen enthält dass es sich lohnt,
   dies auszunutzen.\grqq

\item Die Anzahl der Nicht-Null-Einträge pro Zeile ist durch eine (kleine) Konstante beschränkt,
  die nicht von der Gitterauflösung abhängt.

\item Die Matrix enthält also nur $O(n)$ Einträge
516
517
518
519
520
521
522
523
524
525
526
527
528
529
\end{itemize}

\bigskip

\emph{Konsequenz}: Wir brauchen Algorithmen und Datenstrukturen, die die Dünnbesetztheit ausnutzen.

\medskip Beispiel: Matrix--Vektor-Multiplikation $v = Aw$.

  \begin{algorithm}[H]
   \SetAlgoLined
   % \caption{...}
   \For{alle Zeilen $i$}
   {
     $v_i = 0$\\
530
     \For{alle Spalten $j$ in denen $A_{ij} \neq 0$}
531
532
533
534
535
     {
       $v_i = v_i + A_{ij} w_j$
     }
   }
  \end{algorithm}
536
Das braucht nur $O(\# \text{Nichtnulleinträge})$ Operationen (statt $O(n^2)$).
537

538
539
\bigskip

540
Direkte Verfahren wie Gauß-Elimination oder Cholesky-Zerlegung sind i.A. zu teuer.
541

542
543
544
545
546
547
548
\begin{itemize}
 \item \emph{Erinnerung:} Gauß-Elimination braucht $O(n^3)$ Rechenoperationen!

 \item Wendet man das Gauß-Verfahren auf solch eine Matrix an, so entstehen bei den Zwischenschritten in der Matrix eine beträchtliche Anzahl von zusätzlichen Einträgen ("`fill-in"').

 \item Das Gauß-Verfahren ist deshalb nicht nur zu langsam, es braucht auch zu viel Speicher.
\end{itemize}
549

550
\subsection{Direkte Verfahren für dünnbesetzte Gleichungen}
551

552
553
\begin{itemize}
 \item Es gibt direkte Verfahren zum Lösen von dünnbesetzten Gleichungssystemen.
554

555
 \item Sehr kompliziert / sehr ausgefeilt.
556

557
 \item Die Geschwindigkeit solcher Löser hängt nicht von der Kondition der Matrix ab.
558

559
560
561
562
563
564
565
566
567
568
569
570
 \item Manchmal gibt es Probleme mit Rundungsfehler, wenn Matrizen sehr schlecht
   konditioniert sind.

 \item Für kleine Probleme (In 2021: ca.\ 50\,000 Unbekannte oder weniger)
   sind sie konkurrenzlos schnell.

 \item Danach steigt die Rechenzeit schnell an.

 \item Limitierender Faktor ist aber häufig der Speicherverbrauch.

 \item Dennoch in der Praxis meist die Methode der Wahl.
\end{itemize}
571
572
573
574


\subsection{Lineare iterative Verfahren}

575
576
577
Iterative Verfahren konstruieren eine Folge $\{ x^k \}_{k \in \N}$,
die (hoffentlich) gegen die Lösung $x^*$ von $Ax=b$ konvergiert.

578
\bigskip
579

580
581
Für \emph{lineare} iterative Verfahren~\cite{dahmen_reusken:2008} schreibt man
$Ax=b$ als Fixpunktgleichung
582
\begin{equation*}
583
 x=x+C(b-Ax)
584
\end{equation*}
585
mit einer nichtsingulären Matrix $C \in \R^{n \times n}$.
586

587
\medskip
588

589
Die Größe $r(x) \colonequals b - Ax$
590
nennt man \emph{Residuum}.  Die Matrix $C$ heißt \emph{Vorkonditionierer}.
591

592
593
\bigskip

594
595
Dafür machen wir jetzt eine Fixpunktiteration:
\begin{equation*}
596
  x^{k+1}\colonequals x^k+C (b-Ax^k)=(I-CA)x^k+Cb,
597
598
599
600
  \qquad
  k = 0, 1, 2, \dots
\end{equation*}

601
Der Fehler im $k$-ten Iterationsschritt ist $e^k \colonequals x^k-x^*$. Es gilt
602
603
604
\begin{equation*}
 e^{k+1}=x^{k+1}-x^*
 =
605
\underbrace{(I-CA)}_{\text{Iterationsmatrix}}e^k.
606
\end{equation*}
607
Das erklärt den Namen: Die Abbildung $e^k \mapsto e^{k+1}$ ist linear.
608
609
610
611
612

\begin{definition}
 Die Matrix $I - CA$ heißt \emph{Iterationsmatrix} der Methode.
\end{definition}

613
\subsubsection{Konvergenz:}
614
615

\begin{theorem}
616
 Sei $\norm{\cdot}$ eine submultiplikative Matrixnorm.
617
618
 Das Verfahren konvergiert für jeden Startwert $x^0 \in \R^n$ gegen die Lösung
 von $Ax = b$, wenn $\norm{I - CA} < 1$.
619
620
\end{theorem}

621
622
623
Sei $\rho$ der Spektralradius einer Matrix.  Für jede submultiplikative Matrixnorm
gilt $\rho (B) \leq \norm{B}$, denn für jeden Eigenvektor $v$ von $B$
zum Eigenwert $\lambda$ gilt
624
625
626
627
628
629
630
631
632
\begin{equation*}
 \abs{\lambda}\norm{v}
 =
 \norm{\lambda v}
 =
 \norm{B v}
 \le
 \norm{B}\norm{v}.
\end{equation*}
633
634

Damit erhält man:
635

636
637
638
\begin{theorem}
Sei $\rho (I-CA)$ der Spektralradius von $I-CA$. Das Verfahren konvergiert für jeden Startwert $x^0 \in \R$ gegen die Lösung von $Ax=b$ genau dann, wenn $\rho(I-CA)<1$.
\end{theorem}
639

640
641
642
643
644
645
646
647
648
Je nach Wahl des Vorkonditionierers $C$ erhält man unterschiedliche Verfahren.
Es ergibt sich folgendes Dilemma:
\begin{enumerate}
    \item $C$ soll $A^{-1}$ möglichst gut approximieren,
    \item die Operation $y \mapsto Cy$ soll möglichst billig sein.
\end{enumerate}

Seien $D$, $L$, $R$ die Diagonale bzw.\ der linke und rechte Dreiecksteil von $A$.
Man erhält die folgenden Standardverfahren:
649
\begin{itemize}
650
 \item Das Jacobi-Verfahren
651
  \begin{equation*}
652
653
   C = D^{-1}
  \end{equation*}
654

655
656
657
658
 \item Das Gauß--Seidel-Verfahren
  \begin{equation*}
   C = (D+L)^{-1}
  \end{equation*}
659

660
 \item Das symmetrische Gauß--Seidel-Verfahren
661
  \begin{equation*}
662
   C = (D+L)^{-1}+(D+R)^{-1}-(D+R)^{-1}A(D+L)^{-1}
663
664
665
666
667
668
669
  \end{equation*}
\end{itemize}

\emph{Beachte:} Die Matrix $C$ wird nie explizit ausgerechnet!


\subsection{Das Jacobi-Verfahren}
670
671
672
673

Wir betrachten dieses Verfahren etwas genauer.  Man kann es auch
\glqq direkt\grqq motivieren.

674
Wir nehmen im Folgenden an, dass $A_{ii} \neq 0$ für alle $i=1,\ldots,n$.
675
676
677
678
679

\medskip

Betrachte das lineare Gleichungssystem:
\begin{align*}
680
681
682
    A_{11} x_1+A_{12} x_2+\ldots +A_{1m} x_m & = b_1 \\
    A_{12} x_1+A_{22} x_2+\ldots +A_{2m} x_m & = b_2 \\
    & \;\; \vdots
683
\end{align*}
684
\emph{Idee:} Löse $i$-te Zeile nach $x_i$ auf für $i=1,\ldots,n$:
685
\begin{align*}
686
687
688
    x_1 & =\frac{1}{A_{11}} (b_1-A_{12} x_2-\ldots -A_{1n} x_n ) \\
    x_2 & =\frac{1}{A_{22}} (b_2-A_{21} x_1-\ldots -A_{2n} x_n ) \\
    & \;\; \vdots
689
690
691
\end{align*}
Mache daraus ein iteratives Verfahren:
\begin{align*}
692
693
694
    x_1^{k+1} & =\frac{1}{A_{11}} (b_1-A_{12}x_2^k-\ldots -A_{1n}x_n^k) \\
    x_2^{k+1} & =\frac{1}{A_{22}} (b_2-A_{21}x_1^k-\ldots -A_{2n}x_n^k) \\
    & \;\; \vdots
695
696
697
\end{align*}
Für $i=1,\ldots,n$
\begin{equation*}
698
699
 x_i^{k+1}
 =
700
 \frac{1}{A_{ii}} \bigg(b_i-\sum_{\substack{j=1\\j\neq i}}^n A_{ij}x_j^k \bigg)
701
 =
702
 x^k_i + \frac{1}{A_{ii}} \bigg(b_i-\sum_{j=1}^n A_{ij}x_j^k \bigg).
703
\end{equation*}
704
705
706
Beachte:
\begin{itemize}
 \item Die Rechnungen für die verschiedenen $i$ sind voneinander unabhängig $\implies$ leicht zu parallelisieren.
707

708
709
 \item In der Praxis gilt die Summe natürlich nur über die Einträge der $i$-ten Zeile von $A$, die $\neq 0$ sind.
\end{itemize}
710

711
712
713
714
\begin{exercise}
 Rechnen Sie nach dass dies wirklich ein lineares Verfahren mit Vorkonditionierer
 $C = D^{-1}$ ist.
\end{exercise}
715

716
717
718
719
720
Häufig erweitert man das Jacobi-Verfahren zum \emph{gedämpften} Jacobi-Verfahren.
Dabei wählt man einen Parameter $\eta > 0$ und definiert
\begin{equation*}
 x_i^{k+1}
 =
721
 x^k_i + \frac{\eta}{A_{ii}} \bigg(b_i-\sum_{j=1}^n A_{ij}x_j^k \bigg)
722
723
724
 \qquad \text{bzw.} \qquad
 x^{k+1} = x^k + \eta D^{-1} ( b - Ax^k).
\end{equation*}
725
726
727
728
729

\subsubsection{Konvergenz}

Wir wenden das bekannte Konvergenzkriterium an:
\begin{theorem}
730
Das Jacobi Verfahren konvergiert genau dann, wenn $\rho (I- \eta D^{-1}A ) <1$.
731
732
\end{theorem}

733
734
735
736
737
738
Das Verfahren konvergiert also immer, wenn man nur $\eta$ klein genug wählt.
Allerdings wird das Verfahren für kleine $\eta$ auch immer langsamer.

\medskip

Ein paar direkte Charakterisierungen für den Fall $\eta = 1$ findet man z.B.\ bei \citet{dahmen_reusken:2008}.
739

740
741
Wir greifen die Konvergenztheorie des Jacobi-Verfahrens in Kapitel~\ref{sec:jacobi_smoother}
wieder auf.
742

743
\subsubsection{Konvergenzgeschwindigkeit}
744
745
746
747
748

Die folgende Rechnung stammt aus~\citet[Beispiel~13.10]{dahmen_reusken:2008}.

\bigskip

749
750
Sei $A$ die Matrix des Poisson-Problems auf einem uniformen Gitter
mit Dreiecks- oder Viereckselementen.  Eine Zeile von $A$ ist dann
751
\begin{equation*}
752
    \frac{1}{h^2} (4u_{i,j}-u_{i-1,j}-u_{i+1,j}-u_{i,j-1}-u_{i,j+1} )=f_i.
753
754
\end{equation*}

755
756
757
758
\begin{exercise}
 Rechnen Sie nach dass man tatsächlich diese Matrix erhält.
\end{exercise}

759

760
761
762

Jacobi-Verfahren: $D=4h^{-2}I$.

763
Wir versuchen jetzt, den Spektralradius von $I-CA = I - \eta D^{-1}A$ abzuschätzen:
764
\begin{align*}
765
766
767
768
    \rho (I - \eta D^{-1}A)
    & = \sup_{x \neq 0} \abs[\bigg]{\frac{x^T (I- \eta D^{-1}A) x}{\norm{x}^2}} \\
    & = \sup_{x \neq 0} \abs[\bigg]{1- \eta \frac{x^T D^{-1}Ax}{\norm{x}^2}} \\
    & = \sup_{x \neq 0} \abs[\bigg]{1-\frac{1}{4}\eta h^2 \frac{x^TAx}{\norm{x}^2}}.
769
770
\end{align*}
Der Bruch ist gerade der \emph{Rayleigh-Quotient}
771
772
773
774
775
776
\begin{equation*}
  R(A,x) \colonequals \frac{x^TAx}{\Vert x \Vert^2}.
\end{equation*}
Für symmetrische $A$ gilt $\lambda_\text{min}(A) \le R(A,x) \le \lambda_\text{max}(A)$, und diese Schranken werden
angenommen, wenn $x$ entsprechende Eigenvektoren sind.

777
778
\medskip

779
Daraus folgt dass
780
781
782
\begin{align*}
 \rho (I - \eta D^{-1}A )
 & =
783
 \sup \lbrace \abs[\Big]{1- \frac{1}{4} \eta h^2 \lambda} \; : \; \text{$\lambda$ Eigenwert von $A$} \rbrace \\
784
785
 %
 & =
786
 1-2 \eta \sin^2 ( \frac{1}{2} \pi h ).
787
788
789
790
791
\end{align*}
Mit $\eta = 1$ erhält man:
\todo[inline]{Zeige den allgemeinen Fall!}
\begin{align*}
 \rho (I - \eta D^{-1}A )
792
 =
793
794
 \cos (\pi h).
\end{align*}
795
796
797

Taylor-Entwicklung:
\begin{equation*}
798
    \cos ( \pi h )=1-\frac{1}{2} \pi^2 h^2+\ldots
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
\end{equation*}
Also ist
\begin{equation*}
 \frac{\norm{e^{k+1}}}{\norm{e^k}} \approx \rho (I-D^{-1}A) \approx 1-\frac{1}{2} \pi^2h^2
\end{equation*}
Dieser Ausdruck geht "`quadratisch"' (also ziemlich schnell) gegen $1$ wenn $h \to 0$.

\bigskip

Wir wollen ausrechnen, wie viele Iterationen man ungefähr braucht, um den Anfangsfehler
um einen Faktor $R$ zu reduzieren.  D.h., welches $k$ soll man wählen, um ungefähr
\begin{equation*}
 \frac{\norm{e^k}}{\norm{e^0}} \le \frac{1}{R}
\end{equation*}
zu erhalten?

Da $\frac{\norm{e^k}}{\norm{e^0}} \approx \rho^k$ erhält man
\begin{equation*}
 k
 =
 \log_\rho \frac{1}{R}
 =
821
 \frac{-\ln R}{\ln \rho (I-D^{-1}A )} \approx \frac{-\ln R}{\ln (1-\frac{1}{2}\pi^2 h^2 )} \approx \frac{2}{\pi^2 h^2} \ln R,
822
823
824
\end{equation*}
da $\ln x \approx (x-1)-\frac{1}{2} (x-1)^2 + \ldots$ ist.

825
Für ein doppelt so feines Gitter braucht man also viermal so viele Iterationen
826
827
(und diese sind natürlich auch noch teurer.).

828
829
\medskip

830
Das Jacobi-Verfahren scheint also kein sonderlich gutes Verfahren zu sein.
831
832
833
834
835
836
837

\subsection{Das Gauß-Seidel-Verfahren}

Betrachte noch einmal die Jacobi-Rechenvorschrift:
\begin{equation*}
 x_i^{k+1}
 =
838
839
 \frac{1}{A_{ii}} (b_i-A_{i1}x_1^{k} - A_{i2}x_2^{k}-
   \ldots -A_{i,i-1}x_{i-1}^k - A_{i,i+1}x_{i+1}^k - \ldots - A_{in}x_n^k ).
840
841
842
\end{equation*}
Eigentlich haben wir für $x_1,\ldots,x_{i-1}$ schon bessere Werte als $x_1^k,\ldots,x_{i-1}^k$, nämlich $x_1^{k+1},\ldots,x_{i-1}^{k+1}$.

Sander, Oliver's avatar
Sander, Oliver committed
843
844
\medskip

845
846
Gauß-Seidel-Verfahren:
\begin{align*}
Sander, Oliver's avatar
Sander, Oliver committed
847
848
 x_i^{k+1}
 & =
849
 \frac{1}{A_{ii}} (b_i-A_{i1}x_1^{k+1} - \ldots -A_{i,i-1}x_{i-1}^{k+1}-A_{i,i+1}x_{i+1}^k-\ldots-A_{in}x_n^k ) \\
Sander, Oliver's avatar
Sander, Oliver committed
850
 & =
851
 \frac{1}{A_{ii}} (b_i-\sum_{j=1}^{i-1} A_{ij}x_j^{k+1} -\sum_{j=i+1}^n A_{ij} x_j^k ).
852
853
854
\end{align*}
$x_i^{k+1}$ hängt von $x_{i-1}^{k+1}$ ab $\implies$ keine Parallelisierung möglich.

855
856
857
858
859
\begin{exercise}
 Seien $D$, $L$ Diagonalteil bzw.\ linker Dreiecksteil von $A$.
 Rechnen Sie nach dass das Gauß-Seidel-Verfahren ein lineares Verfahren
 mit $C = (D+L)^{-1}$ ist.
\end{exercise}
860
861
862
863
864
865
866
867
868

\subsubsection{Konvergenz}

Wir erwarten bessere Konvergenzeigenschaften als für das Jacobi-Verfahren.

\medskip

Und in der Tat:
\begin{theorem}
869
870
871
872
Das Gauß-Seidel-Verfahren konvergiert, wenn
\begin{itemize}
 \item $A$ symmetrisch und positiv definit ist und/oder
 \item $A$ irreduzibel und diagonaldominant ist.
873
874
875
876
877
878
\end{itemize}
\end{theorem}
\subsubsection{Konvergenzgeschwindigkeit}
\begin{example}
Sei $A$ die Matrix des Poisson-Problems für ein quadratisches Gebiet.
\begin{itemize}
Sander, Oliver's avatar
Sander, Oliver committed
879
 \item Es gilt $\rho(I-(D+L)^{-1}A) = ( \rho (I-D^{-1}A ) )^2$.
880
881
882
883
  Das heißt der Spektralradius der Jacobi-Methode ist das Quadrat
  des Spektralradius der Gauß-Seidel Methode.  Laut \citet{dahmen_reusken:2008}
  steht das bei \cite{hackbusch:1994}.

884
885
  \todo[inline]{I overlooked Remark 5.5.2. in [3] (On the top of page 131). There is the explicit statement.}

886
887
888
889
890
  \begin{exercise}
   Finden Sie die genaue Stelle, und sagen Sie mir bescheid.
  \end{exercise}


891
892
893
 %
 \item Deshalb
  \begin{equation*}
Sander, Oliver's avatar
Sander, Oliver committed
894
   \rho(I-(D+L)^{-1}A) = \cos^2(\pi h)
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
  \end{equation*}

 \item Taylor-Entwicklung für $\cos^2$
  \begin{equation*}
   \cos^2 \pi h
   \approx
   \Big(1-\frac{1}{2} \pi^2h^2+\ldots \Big)^2
   =
   1-2 \frac{1}{2} \pi^2h^2+\frac{1}{4} \pi^4h^4+\ldots \approx 1-\pi^2h^2
  \end{equation*}
 %
 \item Fehlerreduktion
 \begin{equation*}
  \frac{\norm{e^{k+1}}}{\norm{e^k}} \approx 1-\pi^2h^2
 \end{equation*}
 %
 \item Um den Startfehler um den Faktor $R$ zu reduzieren, braucht man etwa
 \begin{equation*}
  \frac{- \ln R}{\ln \rho (I-(D-L)^{-1}A)}
   \approx
915
  \frac{- \ln R}{\ln (1-\pi^2h^2 )} \approx \frac{\ln R}{\pi^2h^2}
916
917
918
919
920
921
 \end{equation*}
 Iterationen.
\end{itemize}
\end{example}
\emph{Faustregel:} Das Gauß-Seidel-Verfahren braucht nur etwa halb so viele Iterationen wie das Jacobi-Verfahren.

922
923
924
925
926
927
928
929
930
\subsection{Nichtlineare Verfahren}

Zum Vergleich sind hier noch Konvergenzresultate für zwei bekannte
nichtlineare Verfahren.  (Beide ausschließlich für symmetrische Matrizen).

Die Konvergenzresultate sind als Funktion der Konditionszahl $\kappa$
formuliert.

\subsubsection{Das Gradientenverfahren}
931
932
933
934
935
936
937
938
939

\begin{theorem}[{\cite[Kapitel~6]{shewchuk:1994}}]
Sei $\kappa$ die Kondition der Matrix $A$. Dann gilt nach $k$ Schritten des Gradientenverfahrens
\begin{equation*}
 \norm{e^k}_A \leq \Big(\frac{\kappa-1}{\kappa+1}\Big)^k \norm{e^0}_A,
\end{equation*}
wobei $\norm{\cdot}_A$ die Energienorm ist.
\end{theorem}

940
\subsubsection{Das Verfahren der konjugierten Gradienten (CG-Verfahren)}
941
942

\begin{itemize}
943
944
945
946
947
948
 \item Zuerst vorgeschlagen von \citet{hestenes_stiefel:1952} im Jahr 1952.
   Eine gute Beschreibung findet sich bei~\citet{shewchuk:1994}.
 \item Toll: ein direkter Löser für dünnbesetzte lineare GS mit Komplexität $O(n \cdot$ Anzahl der Einträge von $A)$.
 \item Funktioniert in der Praxis aber schlecht: Wegen Rundungsfehlern erhält man
    häufig nach $n$ Iterationen \emph{nicht} die Lösung.

949
950
951
952
 \item Geriet zwischenzeitlich in Vergessenheit
 \item Erlebte Revival als iterative Methode
\end{itemize}

953
954
\begin{theorem}[{\cite{shewchuk:1994}}]
Sei $\kappa$ die Kondition der Matrix $A$. Dann gilt nach $k$ Schritten des CG-Verfahrens
955
\begin{align*}
956
 \norm{e^k}_A \leq 2 ( \frac{\sqrt{\kappa}-1}{\sqrt{\kappa}+1} )^k \cdot \norm{e^0}_A.
957
\end{align*}
958
959
\end{theorem}

960
961
962
\begin{itemize}
 \item Diese Abschätzung für CG ist aber sehr schwach. In vielen Fällen konvergiert das Verfahren deutlich besser!

963
 \item Das CG-Verfahren wird tatsächlich sehr häufig für FE-Probleme verwendet.
964

965
966
 \item Aber man benötigt wieder einen Vorkonditionierer, sonst ist auch
   das CG-Verfahren zu langsam.
967

968
 \item Dabei wird wieder die Gleichung $0 = b - Ax$ durch $0 = C(b-Ax)$ ersetzt.
969

970
971
 \item $C$ muss so gewählt werden dass $\kappa(CA)$ klein ist, dass aber gleichzeitig
   $Cv$ für $v \in \R^n$ billig ausgerechnet werden kann.
972

973
974
 \item Die Konstruktion passender $C$ bildet ein großes Gebiet innerhalb
   der Numerik.
975

976
 \item Eine Möglichkeit: Mehrgitter!
977
978
\end{itemize}

979
\section{Mehrgitterverfahren: Der anschauliche Zugang}
980

981
982
983
Das Ziel ist es, einen iterativen Löser zu konstruieren, dessen Konvergenzgeschwindigkeit
nicht von der Gitterweite $h$ abhängt.

984
\subsection{Glättung}
985

986
Gauß--Seidel- und Jacobi-Verfahren sind für typische Finite-Elemente-Matrizen
987
extrem langsam.
988

989
990
991
992
\medskip

Sie haben allerdings eine wichtige Eigenschaft, die wir uns
für Mehrgitterverfahren zunutze machen:
993

994
995
996
\smallskip

$\longrightarrow$ Die zwei Verfahren \emph{glätten}.
997

998
999
1000
\bigskip

\emph{Beispiel:} Wir betrachten wieder das Poisson-Problem
For faster browsing, not all history is shown. View entire blame