Numerische Integration - Wissenschaftliches Rechnen (i788)

Führen ein konkretes Verfahren zur numerischen Integration ein, welches uns gelegentlich gute Dienste leisten wird, insbesondere auch für das Lösen partieller Differentialgleichungen. Für einen anderen Ansatz und zusätzliche theoretische Untermauerung sei auf die Veranstaltung zur numerischen Mathematik verwiesen.

4.1Idee¶

Um das Integral

I:=\int_a^b f(x)\diff x

(1)

einer stetigen Funktion $f:[a,b]\to\bbR$ näherungsweise zu berechnen, stellen wir $f$ als Linearkombination einfach integrierbarer Basisfunktionen dar:

f(x)\approx\sum_{i=1}^n c_i\,b_i(x).

(2)

Dabei sind $c_1,\ldots,c_n\in\bbR$ die von $f$ abhängenden Koeffizienten und $b_1,\ldots,b_n:[a,b]\to\bbR$ von $f$ unabhängige Funktionen. Aus den Eigenschaften des Integrals (Linearität!) folgt sofort

I=\sum_{i=1}^n c_i\,\int_a^b b_i(x)\diff x.

(3)

Die Idee, allgemeine Funktionen als Linearkombinationen von zweckmäßig gewählten Basisfunktionen darzustellen, ist im wissenschaftlichen Rechnen weit verbreitet, sodass im Folgenden auch über die Bedürfnisse der numerischen Integration hinausgehende Beispiele für Basisfunktionen genannt werden.

4.2Übliche Basisfunktionen¶

4.2.1Stückweise konstante Funktionen¶

Für $a=:x_0<x_1<\cdots<x_n:=b$ setze

b_i(x):=\begin{cases} 1,&\text{falls }x\in[x_{i-1},x_i),\\ 0,&\text{sonst},\end{cases}

(4)

für $i=1,\ldots,n$ .

Dann gilt

f(x)\approx\sum_{i=1}^n f\left(\frac{x_i+x_{i-1}}{2}\right)\,b_i(x)

(5)

und die Integrale der Basisfunktionen sind

\int_a^b b_i(x)\diff x=\frac{1}{x_i-x_{i-1}}.

(6)

4.2.2Stückweise lineare Funktionen¶

Für $a=:x_1<x_2<\cdots<x_n:=b$ setze

{\begin{align} b_1(x)&:=\begin{cases} \frac{x_2-x}{x_2-x_1},&\text{falls }x\in[x_1,x_2],\\ 0,&\text{sonst}, \end{cases}\\ b_n(x)&:=\begin{cases} \frac{x-x_{n-1}}{x_n-x_{n-1}},&\text{falls }x\in[x_{n-1},x_n],\\ 0,&\text{sonst},\end{cases} \end{align}}

(7)

und

b_i(x):=\begin{cases} \frac{x-x_{i-1}}{x_i-x_{i-1}},&\text{falls }x\in[x_{i-1},x_i],\\ \frac{x_{i+1}-x}{x_{i+1}-x_i},&\text{falls }x\in[x_i,x_{i+1}],\\ 0,&\text{sonst},\end{cases}

(8)

für $i=2,\ldots,n-1$ .

Dann gilt

f(x)\approx\sum_{i=1}^n f(x_i)\,b_i(x)

(9)

und die Integrale der Basisfunktionen sind

\int_a^b b_i(x)\diff x=\begin{cases} \frac{x_2-x_1}{2},&\text{für }i=1,\\ \frac{x_{i+1}-x_{i-1}}{2},&\text{für }i=2,\ldots,n-1,\\ \frac{x_n-x_{n-1}}{2},&\text{für }i=n. \end{cases}

(10)

4.2.3Radiale Basisfunktionen¶

Seien $x_1,\ldots,x_n\in[a,b]$ paarweise verschieden und sei $g:\bbR\to\bbR$ eine Funktion, die nur von $|x|$ und nicht direkt von $x$ abhängt. Dann heißen die Funktionen

b_i(x):=g(x-x_i)

(11)

radiale Basisfunktionen. Diese entstehen durch Verschiebung der Funktion $g$ und sind symmetrisch. Üblicherweise hat $g$ Hütchen- oder Glockenform.

Beispiele:

Für $x_i:=(i-1)\,\frac{b-a}{n-1}$ liefert
$g(x):=\begin{cases}1-\frac{n-1}{b-a}\,|x|,&\text{falls }|x|\leq\frac{b-a}{n-1},\\0,&\text{sonst}\end{cases}$
(12)
gerade die Basisfunktionen zur Darstellung stückweise linearer Funktionen.
Mit der Glockenkurve
$g(x)=\rme^{-\frac{x^2}{\sigma^2}}$
(13)
erhält man sehr glatte Approximationen von $f$ , wobei $\sigma>0$ ein Parameter zur Steuerung der Glockenbreite ist. Da sich die Träger (Intervalle mit von Null verschiedenen Funktionswerten) der Basisfunktionen überlappen, können die Koefficienten $c_i$ zur Darstellung von $f$ nicht direkt aus $f$ abgelesen werden. Stattdessen ist ein lineres Gleichungssystem zu lösen:
$f(x_j)=\sum_{i=1}^n c_i\,g(x_j-x_i),\quad j=1,\ldots,n.$
(14)
Die Systemmatrix hat Einsen auf der Diagonale und (überlicherweise) kleine Werte außerhalb der Diagonale.

4.2.4Fourier-Basis¶

Für glatte Funktionen $f:[a,b]\to\bbR$ , die einen Ausschnitt einer periodischen Funktion darstellen, also $f(a)=f(b)$ erfüllen, bietet sich die sogenannte Fourier-Basis an. Zur Vereinfachung der Formeln wird diese gern anders indiziert:

{\begin{align} b_0(x)&:=\frac{1}{b-a}\\ b^{\mathrm{sin}}_i(x)&:=\sqrt{\frac{b-a}{2}}\,\sin\left(2\,\pi\,i\,\frac{x-a}{b-a}\right),\\ b^{\mathrm{cos}}_i(x)&:=\sqrt{\frac{b-a}{2}}\,\cos\left(2\,\pi\,i\,\frac{x-a}{b-a}\right). \end{align}}

(15)

Man erhält also die Folge $b_0,b^{\mathrm{sin}}_1,b^{\mathrm{cos}}_1,b^{\mathrm{sin}}_2,b^{\mathrm{cos}}_2,\ldots$ von Basisfunktionen. Die Koeffizienten zur Darstellung von $f$ sind dann

{\begin{align} c_0&:=\int_a^b f(x)\,b_0(x)\diff x,\\ c^{\mathrm{sin}}_i&:=\int_a^b f(x)\,b^{\mathrm{sin}}_i(x)\diff x,\\ c^{\mathrm{cos}}_i&:=\int_a^b f(x)\,b^{\mathrm{cos}}_i(x)\diff x \end{align}}

(16)

Bezeichnen wir mit

c:=(c_0,c^{\mathrm{sin}}_1,c^{\mathrm{cos}}_1,\ldots,c^{\mathrm{sin}}_n,c^{\mathrm{cos}}_n)

(17)

einen Koeffizientenvektor und mit

g_c(x):=c_0\,b_0(x)+\sum_{i=1}^n c^{\mathrm{sin}}_i\,b^{\mathrm{sin}}_i(x)+c^{\mathrm{cos}}_i\,b^{\mathrm{cos}}_i(x)

(18)

die zugehörige Linearkombination der Basisfunktionen, so erhält man die zu $f$ passenden Koeffizienten als Lösung des Minimierungsproblems

\int_a^b\bigl(f(x)-g_c(x)\bigr)^2\diff x\to\min_{c\in\bbR^{2n+1}}

(19)

(ohne Beweis).

Für die Berechnung der Koeffizienten aus endlichen vielen Funktionswerten von $f$ existieren sehr effiziente Algorithmen (Schnelle Fourier-Transformation).

4.2.5Wavelet-Basen¶

Wavelets sind Funktionen, aus denen man durch systematisches Skalieren und Verschieben Funktionensysteme mit ähnlichen Eigeschaften wie der Fourier-Basis gewinnen kann:

beste Approximation im Sinne des integrierten quadratischen Fehlers,
Berechnung der Koeffizienten mittels sehr effizienter Algorithmen,
Kodierung grober Strukturen in den ersten Koeffizienten, Details in den späteren.

Jedoch bieten Wavelet-Basen viel mehr Freiheiten, da sie nicht zwingend im Frequenzbereich arbeiten (und damit stets überall im Ortsbereich wirken würden), sondern auch lokal im Ortsbereich Manipulationen erlauben.

Das einfachste Beispiel für Wavelets sind die Haar-Wavelets, welche ein Spezialfall der häufig verwendeten Daubechies-Wavelets sind.

4.3Diskretisierte Stammfunktion¶

Ist zu $f:[a,b]\to\bbR$ die Stammfunktion

F:[a,b]\to\bbR,\quad F(z)=f(a)+\int_a^z f(x)\diff x

(20)

gesucht, so muss (theoretisch) für jedes $z$ eine numerische Integration von $f$ über dem Intervall $[a,z]$ durchgeführt werden. Praktisch sind wir nur an $F(z)$ für endlich viele $z=z_i$ mit $a=:z_0<z_1<\cdots<z_n:=b$ interessiert. Wegen $F(z_0)=0$ sollen also $n$ Werte bestimmt werden.

Wurde $f$ auf $[a,b]$ mittels $n$ Basisfunktionen diskretisiert, so besteht zwischen den Koeffizienten $c_1,\ldots,c_n$ und den gesuchten Funktionswerten $F(z_1),\ldots,F(z_n)$ ein linearer Zusammenhang

F(z_j)=\sum_{i=1}^n c_i\,\int_a^{z_j}b_i(x)\diff x,\quad j=1,\ldots,n,

(21)

der als Matrix-Vektor-Produkt

Z=A\,c

(22)

mit

c:={\begin{bmatrix}c_1&\cdots&c_n\end{bmatrix}}^\rmT,

(23)

Z:={\begin{bmatrix}F(z_1)&\cdots&F(z_n)\end{bmatrix}}^\rmT

(24)

und

A:={\begin{bmatrix} \int_a^{z_1}b_1(x)\diff x&\cdots&\int_a^{z_1}b_n(x)\diff x\\ \vdots&&\vdots\\ \int_a^{z_n}b_1(x)\diff x&\cdots&\int_a^{z_n}b_n(x)\diff x \end{bmatrix}}

(25)

geschrieben werden kann. Üblicherweise wählt man die Basisfunktionen so, dass diese zu den Integrations Grenzen $z_1,\ldots,z_n$ passen, also z.B. stückweise lineare Hütchenfunktionen zu den Punkten $z_1,\ldots,z_n$ oder auf den Intervallen $[z_i,z_{i+1}$ konstante Basisfunktionen.

4 Numerische Integration