• Einführung in die Programmierung für Studierende der Physik

• von Dr.phil.nat.Dr.rer.pol. Matthias Hanauske

Theorie: Numerische Integration

Wir betrachten in diesem Unterpunkt die Methode der numerischen Integration mittels der "Geschlossenen Newton-Cotes Gleichungen" (closed Newton-Cotes formulars). Die Vorgehensweise der Herleitung dieser Gleichungen erfolgt, indem man die zu integrierende Funktion $f(x)$ in ein Lagrange Polynom vom Grade N ($P_n(x)$) entwickelt (siehe Anwendungsbeispiel: Interpolation und Polynomapproximation) und dann durch analytische Integration zur Approximation gelangt, die für $N=1$ in die Trapez-Regel, für $N=2$ in die Simpson's-Regel und für $N=3$ in die Simpson's-3/8-Regel übergehen.
Wir möchten das folgende bestimmte Integral numerisch approximieren: $$ \begin{equation} \int_a^b f(x) \, dx \approx \sum_{i=0}^N c_i \,f(x_i) \quad, \end{equation} $$ wobei man das Integrationsintervall in $N$ äquidistante Punkte unterteilt: $x_0=a$, $x_N=b$ und $x_i=a+i\cdot(b-a)/N=a+i\cdot h$ mit ($i=0,1,2,..,N$). Im Folgenden wollen wir die jeweiligen Gleichungen mittels der Methode der Lagrange Polynome ($P_N(x)$) herleiten. Das N-te Lagrange Polynom $P_N(x)$ einer Funktion ist wie folgt definiert: $$ \begin{equation} P_N(x) = \sum_{k=0}^N f(x_k) \cdot L_{N,k}(x) \quad, \end{equation} $$ mit $$ \begin{equation} L_{N,k}(x) = \prod_{i=0 , \, i \neq k}^N \frac{(x-x_i)}{(x_k-x_i)} \quad. \end{equation} $$ Setzt man das Lagrange Polynom $P_N(x)$ in das bestimmte Integral ein, erhält man: $$ \begin{equation} \int_a^b f(x) \, dx \approx \int_a^b \sum_{k=0}^N f(x_k) \cdot L_{N,k}(x) \, dx = \sum_{k=0}^N f(x_k) \cdot \underbrace{\int_a^b L_{N,k}(x) \, dx}_{c_k} := \sum_{k=0}^N f(x_k) \, c_k \quad . \end{equation} $$

Die Trapez-Regel (N=1)

Die einfachste Approximation erhält man mittels eines linearen Lagrange Polynoms $P_1(x)$. Man benutzt dabei lediglich zwei Punkte, nämlich die Grenzen des Intervalls ($x_0=a$ und $x_1=b=x_0+h$ mit $h=(b-a)$). Mittels $P_1(x)$ berechnen wir das folgende bestimmte Integral $$ \begin{equation} \int_a^b f(x) \, dx \approx {\cal I}_1 := \int_{x_0}^{x_1} \sum_{k=0}^1 f(x_k) \cdot L_{1,k}(x) \, dx = \int_{x_0}^{x_1} \sum_{k=0}^1 f(x_k) \cdot \prod_{i=0 , \, i \neq k}^1 \frac{(x-x_i)}{(x_k-x_i)} \, dx \quad, \end{equation} $$ und durch analytische Integration erhält man die folgende Trapez-Regel: $$ \begin{equation} {\cal I}_1 = \frac{(x_1 - x_0)}{2} \cdot \left( f(x_0) + f(x_1) \right) = \frac{h}{2} \cdot \left( f(x_0) + f(x_1) \right) \quad. \end{equation} $$

Die Simpson's-Regel (N=2)

Die Simpson's Approximation erhält man mittels eines quadratischen Lagrange Polynoms $P_2(x)$. Man benutzt dabei drei Punkte, nämlich die Grenzen des Intervalls ($x_0=a$ und $x_2=b =a + 2h$) und einen Wert in der Mitte des Intervalls ($x_1=a + h$) mit $h=(b-a)/2$. Wieder erhält man durch analytische Integration die folgende Simpson's-Regel: $$ \begin{equation} {\cal I}_2 = \frac{h}{3} \cdot \left( f(x_0) + 4 f(x_1) + f(x_2) \right) \end{equation} $$

Die Simpson's-3/8-Regel (N=3)

Die Simpson's-3/8 Approximation erhält man mittels des Lagrange Polynoms $P_3(x)$. Man benutzt dabei vier Punkte, nämlich die Grenzen des Intervalls ($x_0=a$ und $x_3=b =a + 3h$) und zwei zusätzliche Werte ($x_1=a + h$ und $x_2=a + 2h$) mit $h=(b-a)/3$. Man erhält: $$ \begin{equation} {\cal I}_3 = \frac{3 \, h}{8} \cdot \left( f(x_0) + 3 f(x_1) + 3 f(x_2) + f(x_3) \right) \end{equation} $$

Die N=4 Regel

Die $N=4$ Approximation erhält man mittels des Lagrange Polynoms $P_4(x)$. Man benutzt dabei fünf Punkte, nämlich die Grenzen des Intervalls ($x_0=a$ und $x_4=b =a + 4h$) und zwei zusätzliche Werte ($x_1=a + h$, $x_2=a + 2h$ und $x_3=a + 3h$) mit $h=(b-a)/4$. Man erhält: $$ \begin{equation} {\cal I}_4 = \frac{2 \, h}{45} \cdot \left( 7 f(x_0) + 32 f(x_1) + 12 f(x_2) + 32 f(x_3) + 7 f(x_4) \right) \end{equation} $$ Wir betrachten nun z.B. die Funktion $f(x)=10 \cdot {\rm sin}(x) \cdot e^{-\frac{x}{10}}$ und approximieren den Wert für das bestimmte Integral $\int_2^4 f(x) \, dx$. Die nebenstehende Abbildung veranschaulicht für dieses Beispiel die unterschiedlichen Integrationsregeln grafisch.

Stückweise numerische Integration

Möchte man die Funktion jedoch über ein großes Intervall integrieren, so liefern auch die höheren N-Regeln keine genaue Approximation. Man sollte dann eine iterative, stückweise numerische Integration verwenden, wobei man innerhalb der einzelnen Teilstücke eine der oberen Integrationsregeln verwendet. Wir unterteilen dazu das Integrations-Intervall $[a,b]$ in $n/N$-Teilstücke, wobei $N$ die Ordnung der Integrationsregel und $n$ das kleinste gemeinsame Vielfache von $N$ ist. $$ \begin{eqnarray} \int_a^b f(x) \, dx &=& \sum_{j=1}^{n/N} \int_{x_{N \cdot j-N}}^{x_{N \cdot j}} f(x) \, dx \\ &=& \int_{a=x_0}^{x_N} f(x) \, dx + \int_{x_N}^{x_{N \cdot 2}} f(x) \, dx \, + ... \, + \int_{x_{N \cdot (n/N-1)-N}}^{x_{N \cdot (n/N-1)}} f(x) \, dx + \int_{x_{N \cdot (n/N)-N}}^{x_{N \cdot (n/N)}} f(x) \, dx \,\, =\\ &=& \underbrace{\int_{a=x_0}^{x_N} f(x) \, dx}_{\hbox{Teilintervall 1}} + \underbrace{\int_{x_N}^{x_{N \cdot 2}} f(x) \, dx}_{\hbox{Teilintervall 2}} \, + ... \, + \underbrace{\int_{x_{n-2 \cdot N}}^{x_{n - N}} f(x) \, dx}_{\hbox{Teilintervall } n/N-1} + \underbrace{\int_{x_{n - N}}^{x_{n}} f(x) \, dx}_{\hbox{Teilintervall } n/N} \quad, \end{eqnarray} $$ und benutzen in den Teilintervallen eine der oberen Integrationsregeln, so geht man zu einer stückweisen numerischen Integration über.

Unterteilen wir z.B. das Integrations-Intervall in drei Teilstücke und benutzen die Simpson's-Regel ($N=2$, $n=6$ $\rightarrow$ $n/N=3$ Teilstücke), so erhalten wir z.B. die folgenden Integrationsterme: $$ \begin{eqnarray} \int_a^b f(x) \, dx &=& \underbrace{\int_{a=x_0}^{x_2} f(x) \, dx}_{\hbox{Teilintervall 1}} + \underbrace{\int_{x_2}^{x_{4}} f(x) \, dx}_{\hbox{Teilintervall 2}} \, + \underbrace{\int_{x_{4}}^{x_{6}} f(x) \, dx}_{\hbox{Teilintervall 3}} \,\, =\\ &=& \int_{a=x_0}^{x_2} \sum_{k=0}^2 f(x_k) \cdot L_{2,k}(x) \, dx + \int_{x_2}^{x_{4}} \sum_{k=2}^4 f(x_k) \cdot L_{2,k}(x) \, dx \, + \int_{x_{4}}^{x_{6}} \sum_{k=4}^6 f(x_k) \cdot L_{2,k}(x) \, dx \,\, =\\ &=& \underbrace{\frac{h}{3} \cdot \left( f(x_0) + 4 f(x_1) + f(x_2) \right)}_{\hbox{Teilintervall 1}} + \underbrace{\frac{h}{3} \cdot \left( f(x_2) + 4 f(x_3) + f(x_4) \right)}_{\hbox{Teilintervall 2}} \, + \underbrace{\frac{h}{3} \cdot \left( f(x_4) + 4 f(x_5) + f(x_6) \right)}_{\hbox{Teilintervall 3}} \,\, =\\ &=& \frac{h}{3} \cdot \left( f(x_0) + 4 f(x_1) + 2 f(x_2) + 4 f(x_3) + 2 f(x_4) + 4 f(x_5) + f(x_6) \right)\quad. \end{eqnarray} $$

Wir betrachten wieder die Funktion $f(x)=10 \cdot {\rm sin}(x) \cdot e^{-\frac{x}{10}}$ und approximieren den Wert für das bestimmte Integral $\int_2^4 f(x) \, dx$. Wir unterteilen das Integrations-Intervall in drei Teilstücke und benutzen die Simpson's-Regel ($N=2$, $n=6$ $\rightarrow$ $n/N=3$ Teilstücke). Die nebenstehende Abbildung veranschaulicht diese stückweise numerische Integration.
Die gesamte Herleitung der Integrationsregeln ist in dem Jupyter Notebook 'NumericalIntegration.ipynb' einsehbar ( View Notebook, Download Notebook).