Numerische Integration

Bei der numerischen Quadratur versucht man, den Wert eines Riemann-Integrals näherungsweise zu bestimmen.

Anschauliche Darstellung des Integrals als Flächeninhalt unter einer Kurve

Das Riemann-Intergral

$J(f)=\int _{a}^{b}f(x)dx$

ist der Flächeninhalt unterhalb der Kurve der Funktion f(x), wobei wir zunächst annehmen, dass f(x) > 0 im Intervall [a,b] ist.

Oft kann man das Integral nicht geschlossen lösen, d.h. man kann keine Stammfunktion zu f(x) angeben. Deshalb versucht man, Näherungswerte zu ermitteln.

Dazu unterteilt man die gesuchte Fläche in senkrechte Streifen und nähert jede dieser so erhaltenen Teilflächen durch einfache geometrische Figuren (z.B. Trapez) oder einfache Funktionen (z.B. Polynome) an.

Für die Flächenberechnung dieser einfachen Figuren benötigt man den Wert der Funktion f(x) an den so genannten Stützstellen x₀, ... x_m. Die Summe über diese Teilflächen ergibt eine Näherung des Integrals. Je schmaler man die einzelnen Teilflächen wählt desto genauer wird die Näherung. Von Interesse ist dann noch die Frage, wie groß der Fehler ist, der sich durch die Näherung ergibt. Dieser Fehler wird durch das Restglied beschrieben. Um die Anzahl der Funktionsauswertungen zu minimieren, bei gleichzeitiger Möglichkeit den Fehler zu kontrollieren, verwendet man oft das Rombergsche Extrapolationsverfahren. Hierbei werden die Integralwerte von immer kleiner werdenden 'Streifen' zu einer verschwindenden Breite hin extrapoliert.

Allgemeine Quadraturformel

Mit Hilfe von Interpolationspolynomen und deren Lagrange-Darstellung kann man die folgende allgemeine Quadraturformel und das zugehörige Restglied herleiten.

Die allgemeine Quadraturformel für eine Teilfläche lautet

$Q(f)=\sum _{j=0}^{m}\beta _{j}(b-a)^{j+1}f(x_{0},...,x_{j})$

mit den Koeffizienten

$\beta _{j}=\left\{{\begin{matrix}1&{\mbox{für }}j{\mbox{ = 0}}\\\int _{0}^{1}(z-z_{0})...(z-z_{j-1})dz&{\mbox{für }}j{\mbox{ =1,2,...,m}}\end{matrix}}\right.$

$z_{j}={\frac {x_{j}-a}{b-a}}$

Das Restglied beträgt

$E(f)=\int _{a}^{b}(x-x_{0})...(x-x_{m})f(x_{0},..,x_{m},x)dx$

Ist die Funktion f(x) im Intervall [a,b] (m+1)-mal stetig differenzierbar ("reellwertig" wird nicht gefordert), dann läßt sich das Restglied nach oben abschätzen durch

$\left|E(f)\right|\leq {(b-a)^{m+2} \over (m+1)!}\int _{0}^{1}{\left|(z-z_{0})...(z-z_{m})\right|dz}\ \max _{a\leq x\leq b}{\left|f^{(m+1)}(x)\right|}$

Wenn noch zusätzlich für alle Stützstellen im Intervall [a,b] gilt (x - x^j) >=0 oder alternativ (x - x^j) <=0, dann hat der Integrand keinen Vorzeichenwechsel in [a,b] und man kann zeigen:

$\int _{0}^{1}{\left|(z-z_{0})...(z-z_{m})\right|dz}={\left|\beta _{m+1}\right|}$

Daraus folgt dann die Restgliedabschätzung

$\left|E(f)\right|\leq {(b-a)^{m+2} \over (m+1)!}{\left|\beta _{m+1}\right|}\ \max _{a\leq x\leq b}{\left|f^{(m+1)}(x)\right|}$

Ist die Funktion f(x) zusätzlich noch reellwertig in [a,b], dann kann man mit Hilfe des Mittelwertsatzes der Integralrechnung folgende Darstellung für das Restglied herleiten:

$E(f)={(b-a)^{m+2} \over (m+1)!}{\beta _{m+1}}{f^{(m+1)}(\zeta )}$

mit einer Zwischenstelle ζ im Intervall [a,b].

Summierte Quadraturformeln

Um das Integral noch besser annähern zu können unterteilt man das Intervall [a,b] in nebeneinanderliegende kleinere Teilintervalle. In jedem Teilintervall wendet man die gleiche Näherung für die einzelnen Flächen an und addiert danach die entstandenen Näherungen.

Unterteilen wir also das Intervall [a,b] in N nebeneinanderliegende Teilntervalle

[a₁,b₁],[a₂,b₂], ... [a_N,b_N]     mit a₁=a;    a_k+1=b_k  k=1,...,N-1;   b_N = b

Die Teilintervalle müssen zunächst nicht die gleiche Länge haben.

Dann gilt für jede Teilfläche

$\int _{a_{k}}^{b_{k}}f(x)dx=Q_{k}(f)+E_{k}(f)\quad k=1,...,N$

Daraus folgt für das gesamte Integral

$J(f)=\int _{a}^{b}f(x)dx=\sum _{k=1}^{N}\int _{a_{k}}^{b_{k}}f(x)dx=Q(f)+E(f)$

mit

$Q(f)=\sum _{k=1}^{N}Q_{k}(f)$

$E(f)=\sum _{k=1}^{N}E_{k}(f)$

Sei f(x) nun (m+1) mal stetig differenzierbar im Gesamtintervall [a,b]. Ferner sollen ab jetzt alle Teilintervalle die gleiche Länge h haben, also

$h=b_{k}-a_{k}={\frac {b-a}{N}}\quad k=1,...,N$

Dann gilt für die einzelnen Restglieder (s.oben)

$\left|E_{k}(f)\right|\leq {h^{m+2} \over (m+1)!}{\left|\beta _{m+1}\right|}\ \max _{a_{k}\leq x\leq b_{k}}{\left|f^{(m+1)}(x)\right|}$

Summierung über die einzelnen Restglieder ergibt die Abschätzung für das gesamte Restglied

$\left|E(f)\right|\leq \sum _{k=1}^{N}\left|E_{k}(f)\right|\leq (b-a){h^{m+1} \over (m+1)!}{\left|\beta _{m+1}\right|}\max _{a_{k}\leq x\leq b_{k}}{\left|f^{(m+1)}(x)\right|}$

mit Nh = b - a.

Ist die Funktion f(x) zudem auf [a,b] reellwertig, dann kann man für das Restglied analog herleiten:

$E(f)=(b-a){h^{m+1} \over (m+1)!}{\beta _{m+1}}{f^{(m+1)}(\zeta )}$

Spezielle Quadraturformeln

Man hat nun verschiedene Möglichkeiten, die einzelnen Teilfächen durch spezielle einfachere Flächen anzunähern. Die Anwendung der obigen allgemeinen Quadraturformel auf diese speziellen Flächen liefert einige bekannte und wichtige spezielle Quadraturformeln.

Sehnentrapezformel

Man ersetzt die Kurve f(x) durch die Verbindungsgerade zwischen den Punkten (a,f(a)) und (b,f(b)) - also durch die Sehne - und erhält somit ein Trapez.

Mit

• m = 1
• x₀ = a
• x₁ = b

erhält man die Koeffizienten β_j

• z₀ = 0
• z₁ = 1
• β₀ = 1
• β₁ = 0,5
• β₂ = -1/6

und daraus schließlich die Sehnentrapezformel

$Q(f)={\frac {b-a}{2}}(f(a)+f(b))$

Tangententrapezformel

Man legt an die Kurve f(x) im Punkt c in der Mitte des Intervalls [a,b] die Tangente und erhält so wieder ein Trapez.

Mit

• m = 1
• x₀ = c
• x₁ = c

erhält man die Koeffizienten β_j

• z₀ = 0,5
• z₁ = 0,5
• β₀ = 1
• β₁ = 0
• β₂ = 1/12

und daraus schließlich die Tangententrapezformel

$Q(f)=(b-a)\ f({a+b \over 2})$

Simpsonsche Formel oder Keplersche Fassregel

Interpoliert man die Funktion f(x) mittels eines quadratischen Polynoms in den Punkten a,b,c (c liegt in der Mitte von [a,b]), dann erhält man die Simpsonsche Formel.

Mit

• m = 3
• x₀ = a
• x₁ = b
• x₂ = c
• x₃ = c

erhält man die Koeffizienten β_j

• z₀ = 0
• z₁ = 1
• z₂ = 0,5
• z₃ = 0,5
• β₀ = 1
• β₁ = 1/2
• β₂ = -1/6
• β₃ = 0
• β₄ = -1/120

und damit die Simpsonsche Formel

$Q(f)={\frac {b-a}{6}}\ (f(a)+4f({\frac {a+b}{2}})+f(b))$

Hermitsche Formel

Mit

• m = 3
• x₀ = a
• x₁ = b
• x₂ = a
• x₃ = b

ergeben sich

• z₀ = 0
• z₁ = 1
• z₂ = 0
• z₃ = 1
• β₀ = 1
• β₁ = 1/2
• β₂ = -1/6
• β₃ = -1/12
• β₄ = -1/30