9 Numerische Integration

Ziel der numerischen Integration (auch: numerische Quadratur) ist das Berechnen von bestimmten Integralen

für gegebene Funktionen $f:[a,b]\to\bbR$. Dabei soll nicht die Stammfunktion von $f$ genutzt werden, sondern nur die Funktionswerte an endlichen vielen Stellen $x_1,\ldots,x_n\in[a,b]$.

9.1Idee¶

Verschiedene Ansätze führen zum Ziel. Die wichtigste Gruppe von numerischen Integrationsverfahren sind die sogenannten interpolatorischen Verfahren. Hier wird die Funktion $f$ an zweckmäßig (siehe unten) gewählten Stützstellen interpoliert, um dann das Integral der Interpolante zu berechnen. Sofern der Interpolationsfehler klein genug ist, wird der so berechnete Wert in der Nähe des tatsächlichen Wertes des Integrals liegen.

Neben einfacher Polynominterpolation kommen vorallem stückweise Polynome und Splines zum Einsatz. Integrale von Polynomen lassen sich mit den Computer exakt berechnen. Neben der Frage der Kondition ist hier, analog zu den Lösungsverfahren für nichtlineare Gleichungen, insbesondere die Frage der Konvergenz zu klären, da auch beim exakten Rechnen ohne Rundungsfehler nur eine Näherungslösung zu erwarten ist.

Numerische Integrationsverfahren können auch danach beurteilt werden, für welche Funktionsklassen sie exakte Lösungen liefern. Meist sind diese Polynome bis zu einem gewissen Grad.

9.2Kondition¶

Die Kondition ist somit gut, außer für Funktion $f$ mit $\int_a^b f(x)\diff x\approx 0$. In diesem Fall führen auch kleine absolute Ausgabefehler zu großen relativen Fehlern aufgrund von Auslöschungseffekten.

9.3Konvergenz¶

Bezeichnen wir mit $I_n$ den auf der Grundlage von Funktionsauswertungen an $n$ Stützstellen $x_1,\ldots,x_n$ erhaltenen Näherungswert für das exkate Integral $I$, so sollte ein Verfahren für die numerische Integration die Bedingung

erfüllen. Man sagt dann: Das verfahren konvergiert.

Im Gegensatz zu Lösungsverfahren für nichtlineare Gleichungen spielt der Begriff der Konvergenzordnung bei der numerischen Integration nur eine untergeordnete Rolle.

9.4Exaktheit¶

Werden $n$ Stützstellen verwendet, so ist klar, dass man stets ein Verfahren findet, welches exakt vom Grad $n-1$ ist (ersetze $f$ durch das entsprechende Interpolationspolynom). Weiter unten werden wir sehen, dass auch Verfahren mit Exaktheitsgrad $2\,n-1$ möglich sind.

9.5Allgemeiner Ansatz¶

Da sinnvollerweise $I_n\geq 0$ für alle nichtnegativen Funktionen $f$ gelten soll, ist klar, dass nur positive Gewichte $w_k$ in Frage kommen.

Soll das Verfahren wenigstens exakt vom Grad 0 sein (also konstante Funktionen exakt integrieren), so muss

gelten (IDVID 920).

9.6Newton-Cotes-Formeln¶

9.6.1Idee¶

Wählt man die Stützstellen äquidistant, also

und ersetzt die zu integrierende Funktion $f$ durch das entsprechende Interpolationspolynom vom Grad $n-1$, so erhält man die so genannten Newton-Cotes-Quadraturformeln. Sind $L_{n-1,1},\ldots,L_{n-1,n}$ die Langrange’schen Basispolynome, so folgt

(IDVID 930).

Newton-Cotes-Formeln liefern Verfahren mit Exaktheitsgrad $n-1$. Man kann sogar Exaktheitsgrad $n$ zeigen bei ungeradem $n$ (tun wir hier nicht), wenn man Symmetrieeigenschaften der Lagrange-Basispolynome ausnutzt. Allerdings ist auch bekannt, dass ab $n=9$ negative Gewichte auftreten. Somit ist der sinnvolle Einsatz nur bei sehr geringer Stützstellenanzahl möglich.

Ein Ausweg sind sogenannte summierte Formeln: Für je $m$ benachbarte Stützstellen wende die Newton-Cotes-Formel an und addiere alle so erhaltenen Teilresultate. Ist $n=s\,m$, so wird die Newton-Cotes-Formel also $s$-mal angewendet. Dieses Verfahren ist dann exakt vom Grad $m-1$.

9.6.2Approximationsfehler¶

Ausgehend von der Abschätzung des Approximationsfehlers bei der Polynominterpolation kann man folgendes Resultat zeigen:

Je mehr Stützstellen genutzt werden, desto kleiner ist $h$ und desto kleiner wird der Approximationsfehler werden. Je größer $m$, desto schneller sinkt der Fehler bei Erhöhung der Stützstellenanzahl.

9.6.3Beispiele¶

9.6.3.1Rechteckregel¶

Für $m=1$ (Interpolationspolynome vom Grad 0) wird pro Interpolationspolynom nur eine Stützstelle benötigt, sodass nicht klar ist, auf welchem Intervall die einzelnen Interpolationspolynome zum Einsatz kommen. Zwei Varianten sind üblich ($h:=\frac{1}{n-1}$):

• $[x_1,x_1+\frac{h}{2}]$, $[x_k-\frac{h}{2},x_k+\frac{h}{2}]$, $[x_n-\frac{h}{2},x_n]$, $k=2,\ldots,n-1$,

• $[x_k,x_{k+1}]$, $k=1,\ldots,n-1$.

Bei der ersten Variante ist die Quadraturformel etwas komplizierter:

Bei der zweiten Variante ist die Formel einfacher, aber die letzte Stützstelle (Intervallende) hat keinerlei Einfluss:

Verdoppelt man die Stützstellenanzahl, so wird sich der Approximationsfehler etwa halbieren.

9.6.3.2Trapezregel¶

Für $m=2$ werden pro Interpolationspolynom zwei Stützstellen benötigt. Jedes Interpolationspolynom deckt also den Bereich zwischen zwei benachbarten Stützstellen ab. Statt der Funktion $f$ wird somit das Integral einer stückweise linearen Approximation berechnet.

Bei Verdoppelung der Stützstellenanzahl wird der Approximationsfehler etwa auf ein Viertel reduziert.

9.6.3.3Simpson-Regel¶

Für $m=3$ werden pro Interpolationspolynom drei Stützstellen benötigt. Jedes Interpolationspolynom deckt also den Bereich zwischen drei benachbarten Stützstellen ab. Statt der Funktion $f$ wird somit das Integral einer stückweise quadratischen Approximation berechnet. Achtung: Die Anzahl der Stützstellen muss hier ungerade sein.

9.7Gauss-Quadratur¶

Wählt man die Stützstellen für die Polynominterpolation nicht äquidistant, kann man unter Umständen numerische Integrationsverfahren finden, deren Exaktheitsgrad bei $n$ Stützstellen größer als $n-1$ ist. Der maximal mögliche Exaktheitsgrad liegt bei $2\,n-1$ (IDVID 970).

Man kann zeigen (zun wir hier nicht), dass der maximale Exaktheitsgrad erreicht wird, wenn die Stützstellen gerade die Nullstellen des Legendre-Polynoms vom Grad $n$ sind und die Gewichte wie bei den Newton-Cotes-Formeln als Integrale der Lagrange-Basispolynome gewählt werden. Das so konstruierte, aus Sicht des Exaktheitsgrades optimale numerische Integrationsverfahren wird als Gauß-Quadratur bezeichnet (manchmal auch Gauß-Legendre-Quadratur).

Auch kann man zeigen, dass die Gewichte auch bei größerem $n$ im Gegensatz zu den Newton-Cotes-Formeln stets positiv sind, sodass die Gauss-Quadratur weniger anfällig für Auslöschungseffekte ist als Newton-Cotes-Formeln.

Beachte: Die Gauß-Quadraturformeln werden meist nur für das Intervall $[-1,1]$ angegeben. Soll über ein anderes Intervall $[a,b]$ integriert werden, so sind als Stützstellen und Gewichte

zu verwenden, wobei $\tilde{x}_k$ und $\tilde{w}_k$ die Stützstellen und Gewichte für das Intervall $[-1,1]$ sind.

9.8Monte-Carlo-Verfahren¶

Das Integral

kann als Erwartungswert einer Zufallsgröße interpretiert werden. Sei dazu $X$ eine auf $[a,b]$ gleichverteilte Zufallsgröße und sei $Y:=f(X)$ die durch Anwenden von $f$ daraus entstehende transformierte Zufallsgröße. Mit der Dichtefunktion $p_X(x)=\frac{1}{b-a}$ zu $X$ gilt

Somit

Den Erwartungswert $\bbE X$ der Zufallsgröße $X$ kann man über den Mittelwert aus $n$ Realisierungen $x_1,\ldots,x_n$ annähern. Setzen wir

so gilt nach dem Gesetz der großen Zahlen in gewissem Sinne also $\lim\limits_{n\to\infty}I_n=I$.

Die Konvergenzgeschwindigkeit ist jedoch relativ gering. Der Fehler $|I-I_n|$ verhält sich wie $\frac{1}{\sqrt{n}}$ für große $n$ (siehe Konvergenzgeschwindigkeit im Gesetz der großen Zahlen. Das Gesetz vom iterierten Logarithmus (Josef Steinebach, Uni Köln) und Gesetz des iterierten Logarithmus für Details). Die Trapezregel liefert hingegen beispielsweise einen Fehler in der Größenordnung $\frac{1}{n^2}$. Aber: Möchte man mehrdimensionale Integrale berechnen, so bleibt beim Monte-Carlo-Verfahren der Fehler in der Größenordnung $\frac{1}{\sqrt{n}}$, während er bei der Trapezregel auf $n^{-\frac{2}{d}}$ steigt ($d$ ist die Dimension des Integrationsbereichs).

Neben der niedrigen Konvergenzgeschwindigkeit für eindimensionale Integrale ist ein weiter Nachteil des Monte-Carlo-Verfahrens, dass bei jeder Auswertung des Integrals ein etwas anderer Wert ermittelt wird.