[WS] Numerische Integration - Beispiele

Neue Frage »

01.04.2007, 23:51

tigerbine

Auf diesen Beitrag antworten »

[WS] Numerische Integration - Beispiele

Gliederung

Abgeschlossene Newton-Cotes Formeln
1a. Trapezregel
1b. Simpson-Regel oder Keplersche Faßregel
1c. Pulcherrima oder 3/8 Regel
Offene Newton-Cotes Formeln
2a. Mittelpunktsregel
2b. offene NC für n=1
2c. offene NC für n=2
Bezug zum Riemann-Integral
Vorgegebene Genauigkeit erreichen
4a. Summierte Trapezregel
4b. Ordnung von Quadraturformeln
Adaptive Formeln
Romberg-Integration

02.04.2007, 00:20

tigerbine

Auf diesen Beitrag antworten »

1. Abgeschlossene Newton-Cotes Formeln
Berechnung der Gewichte und Quadraturfehler für n = 1, 2, 3

Dafür verwenden wir die im Theorieteil erhaltene Gestalt:

$\begin{eqnarray*} \alpha_j = (b-a)\cdot \int_{0}^{1}~\prod_{k=0 \atop k \neq j}^{n}~\frac{nt-k}{j-k}~dt \end{eqnarray*}$

02.04.2007, 00:50

tigerbine

Auf diesen Beitrag antworten »

1a. Trapezregel
Berechnung der Gewichte

$\begin{eqnarray*} n=1 \end{eqnarray*}$

$\begin{eqnarray*} \alpha_0 &&= (b-a)\cdot \int_{0}^{1}~\prod_{k=0 \atop k \neq j}^{1}~\frac{t-k}{j-k}~dt =(b-a)\cdot \int_{0}^{1}~\frac{t-1}{0-1}~dt =(b-a) \cdot [-\frac{1}{2}t^2+t]_{0}^{1} \\&& = (b-a) \cdot \frac{1}{2} \end{eqnarray*}$

$\begin{eqnarray*} \alpha_1 &&= (b-a)\cdot \int_{0}^{1}~\prod_{k=0 \atop k \neq j}^{1}~\frac{t-k}{j-k}~dt =(b-a)\cdot \int_{0}^{1}~\frac{t-0}{1-0}~dt =(b-a) \cdot [\frac{1}{2}t^2]_{0}^{1} \\&& = (b-a) \cdot \frac{1}{2} \end{eqnarray*}$

Quadraturformel

Als Näherung für das Integral $\begin{eqnarray*} \int_a^b~f(x)~dx = I(f) \end{eqnarray*}$ erhält man hier die Quadraturformel:

$\begin{eqnarray*} I_1(f)&&=\alpha_0 \cdot f(x_0) + \alpha_1\cdot f(x_1) = h\cdot \frac{1}{2} \cdot f(a) + h \cdot \frac{1}{2}\cdot f(b) = \\ &&=(b-a) \cdot \frac{f(a)+f(b)}{2} \end{eqnarray*}$

Geometrische Interpretation

http://upload.wikimedia.org/wikipedia/commons/thumb/b/bb/Sehnentrapezformel.svg/373px-Sehnentrapezformel.svg.png

Anschaulich hätte man die Formel auch aus der Berechnung des Flächeninhalts des roten Trapezes erhalten können. Desweiteren erkennt man, dass die interpolierende Polynomfunktion linear ist.

Berechnung des Quadraturfehlers

Hier ist n=1 ungerade, daher nehmen wir die Standardfehlerformel.

$\begin{eqnarray*} I(f)-I_1(f) &&= \int_a^b~\frac{f^{(2)}(\xi_x)}{2!} \cdot (x-a)(x-b)~dx \\&& \stackrel{(1),(2)} = \frac{f^{(2)}(\xi)}{2} \cdot (b-a)^3\cdot \int_0^1~(t-0)(t-1)~dt \\&& = -\frac{f^{(2)}(\xi)}{12} \cdot (b-a)^3 \end{eqnarray*}$

Dabei ist die Funktion (t-0)(t-1) auf [0,1] nicht positiv.

02.04.2007, 02:52

tigerbine

Auf diesen Beitrag antworten »

1b. Simpson-Regel oder Keplersche Faßregel
Berechnung der Gewichte

$\begin{eqnarray*} n=2 \end{eqnarray*}$

$\begin{eqnarray*} \alpha_0 &&= (b-a)\cdot \int_{0}^{1}~\prod_{k=0 \atop k \neq j}^{2}~\frac{2t-k}{j-k}~dt =(b-a)\cdot \int_{0}^{1}~\frac{2t-1}{0-1} \cdot \frac{2t-2}{0-2}~dt =(b-a) \cdot [\frac{2}{3}t^3-\frac{3}{2}t^2+t]_{0}^{1} \\&& = (b-a) \cdot \frac{1}{6} \end{eqnarray*}$

$\begin{eqnarray*} \alpha_1 &&= (b-a)\cdot \int_{0}^{1}~\prod_{k=0 \atop k \neq j}^{2}~\frac{2t-k}{j-k}~dt =(b-a)\cdot \int_{0}^{1}~\frac{2t-0}{1-0} \cdot \frac{2t-2}{1-2}~dt =(b-a) \cdot [-\frac{4}{3}t^3+2t^2]_{0}^{1} \\&& = (b-a) \cdot \frac{4}{6} \end{eqnarray*}$

$\begin{eqnarray*} \alpha_2 &&= (b-a)\cdot \int_{0}^{1}~\prod_{k=0 \atop k \neq j}^{2}~\frac{2t-k}{j-k}~dt =(b-a)\cdot \int_{0}^{1}~\frac{2t-0}{2-0} \cdot \frac{2t-1}{2-1}~dt =(b-a) \cdot [\frac{2}{3}t^3-\frac{1}{2}t^2]_{0}^{1} \\&& = (b-a) \cdot \frac{1}{6} \end{eqnarray*}$

Quadraturformel

Als Näherung für das Integral $\begin{eqnarray*} \int_a^b~f(x)~dx = I(f) \end{eqnarray*}$ erhält man hier die Quadraturformel:

$\begin{eqnarray*} I_2(f)&&=\alpha_0 \cdot f(x_0) + \alpha_1\cdot f(x_1) + \alpha_2\cdot f(x_2)= \\ && = \frac{b-a}{2}\cdot \frac{1}{3}\cdot f(a) +\frac{b-a}{2}\cdot \frac{4}{3}\cdot f(\frac{a+b}{2})+\frac{b-a}{2}\cdot \frac{1}{3}\cdot f(b) = \\ &&= \frac{b-a}{6}\cdot \Bigl\{ f(a) + 4\cdot f(\frac{a+b}{2}) + f(b)\Bigr\} \end{eqnarray*}$

Geometrische Interpretation

http://upload.wikimedia.org/wikipedia/commons/thumb/0/08/Simpson_rule.png/800px-Simpson_rule.png

Warum die Regel auch Keplersche Faßregel heißt, kann hier zum Einstieg nachgelesen werden. Ansonsten erkennt man in der Grafik schön, dass die interpolierende Polynomfunktion den Grad 2 hat.

Berechnung des Quadraturfehlers

Da n=2 gerade ist, wählen wir einen weiteren Knoten $\begin{eqnarray*} x_0=a,~x_1=x_2=\frac{b-a}{2},~x_3=b \end{eqnarray*}$ . Mit der Fehlerformel für gerades n ergibt sich dann:

$\begin{eqnarray*} I(f)-I_2(f) && = \int_a^b~\frac{f^{(4)}(\xi_x)}{4!} \cdot (x-a)(x-\frac{(b-a)}{2})^2(x-b)~dx \\&& \stackrel{(1),(2)} = \frac{f^{(4)}(\xi)}{24} \cdot (b-a)^5 \cdot \int_0^1~(t-0)(t-0.5)^2(t-1)~dt \\&& = -\frac{f^{(4)}(\xi)}{2880} \cdot (b-a)^5 \end{eqnarray*}$

Dabei ist die Funktion nichtpositiv auf [0,1]:

12.04.2007, 15:52

tigerbine

Auf diesen Beitrag antworten »

1c. Pulcherrima oder 3/8 Regel
Berechnung der Gewichte

$\begin{eqnarray*} n=3 \end{eqnarray*}$

Quadraturformel

Als Näherung für das Integral $\begin{eqnarray*} \int_a^b~f(x)~dx = I(f) \end{eqnarray*}$ erhält man hier die Quadraturformel:

$\begin{eqnarray*} I_3(f)&&=\alpha_0 \cdot f(x_0) + \alpha_1\cdot f(x_1) + \alpha_2\cdot f(x_2) +\alpha_3\cdot f(x_3)= \\ && = \frac{b-a}{3}\cdot \frac{3}{8}\cdot f(a) +\frac{b-a}{3}\cdot \frac{9}{8}\cdot f(x_1)+\frac{b-a}{3}\cdot \frac{9}{8}\cdot f(x_2)+\frac{b-a}{3}\cdot \frac{3}{8}\cdot f(b) = \\ &&= \frac{b-a}{8}\cdot \Bigl\{ f(a) + 3\cdot f(x_1) +3\cdot f(x_2)+ f(b)\Bigr\} \end{eqnarray*}$

12.04.2007, 17:01

tigerbine

Auf diesen Beitrag antworten »

2. Offene Newton-Cotes Formeln
Berechnung der Gewichte und Quadraturfehler für n = 0, 1, 2

Dafür verwenden wir die im Theorieteil erhaltene Gestalt:

$\begin{eqnarray*} \alpha_j=(b-a) \cdot \int_{0}^{1}~\prod_{k=0 \atop k \neq j}^n \frac{(n+2)t-(k+1)}{j-k}~dt \end{eqnarray*}$

12.04.2007, 19:30

tigerbine

Auf diesen Beitrag antworten »

2a. Mittelpunktsregel
Berechnung der Gewichte

$\begin{eqnarray*} n=0 \end{eqnarray*}$

$\begin{eqnarray*} \alpha_0 = (b-a) \cdot \int_{0}^{1}~1~ dt = (b-a) \cdot [t]_0^1=(b-a)\cdot 1 \end{eqnarray*}$

Quadraturformel

Als Näherung für das Integral $\begin{eqnarray*} \int_a^b~f(x)~dx = I(f) \end{eqnarray*}$ erhält man hier die Quadraturformel:

$\begin{eqnarray*} I_0(f)&&=\alpha_0 \cdot f(x_0) = h\cdot 2 \cdot f(\frac{b-a}{2})= \\ &&=(b-a) \cdot f(\frac{b-a}{2}) \end{eqnarray*}$

Geometrische Interpretation

http://upload.wikimedia.org/wikipedia/commons/thumb/b/bc/Mittelpunktsregel.svg/386px-Mittelpunktsregel.svg.png

Quadraturfehler

Da n gerade ist, wählen wir einen weiteren Knoten $\begin{eqnarray*} x_1:=x_0 = \frac{b-a}{2} \end{eqnarray*}$ . Mit der Fehlerformel für gerades n ergibt sich dann:

$\begin{eqnarray*} I(f)-I_0(f) && = \int_a^b~\frac{f^{(2)}(\xi_x)}{2!} \cdot (x-\frac{(b-a)}{2})^2~dx \\&& \stackrel{(1),(2)} = \frac{f^{(2)}(\xi)}{2} \cdot (b-a)^3 \cdot \int_0^1~(t-0.5)^2~dt \\&& = \frac{f^{(2)}(\xi)}{24} \cdot (b-a)^3 \end{eqnarray*}$

12.04.2007, 19:49

tigerbine

Auf diesen Beitrag antworten »

2b. offene NC für n=1
Berechnung der Gewichte

$\begin{eqnarray*} n=1 \end{eqnarray*}$

$\begin{eqnarray*} \alpha_0 &&=(b-a)\cdot \int_0^1~\frac{3t-2}{0-1}~dt =(b-a)\cdot \int_0^1~-3t+2~dt = (b-a) \cdot [-\frac{3}{2}t^2+2t]_0^1 =\\&&=(b-a) \cdot \frac{1}{2} \end{eqnarray*}$

$\begin{eqnarray*} \alpha_1 &&=(b-a)\cdot \int_0^1~\frac{3t-1}{1-0}~dt =(b-a)\cdot \int_0^1~3t-1~dt = (b-a) \cdot [\frac{3}{2}t^2-t]_0^1 =\\&&=(b-a) \cdot \frac{1}{2} \end{eqnarray*}$

Quadraturformel

Als Näherung für das Integral $\begin{eqnarray*} \int_a^b~f(x)~dx = I(f) \end{eqnarray*}$ erhält man hier die Quadraturformel:

$\begin{eqnarray*} I_1(f)&&=\alpha_0 \cdot f(x_0) + \alpha_1\cdot f(x_1) =\\&&=(b-a)\cdot \frac{f(x_0)+f(x_1)}{2} \end{eqnarray*}$

12.04.2007, 20:40

tigerbine

Auf diesen Beitrag antworten »

2c. offene NC für n=2
Berechnung der Gewichte

$\begin{eqnarray*} n=2 \end{eqnarray*}$

$\begin{eqnarray*} \alpha_0 &&= (b-a) \cdot \int_0^1~\frac{4t-2}{0-1} \cdot \frac{4t-3}{0-2}~dt =(b-a) \cdot \int_0^1~8t^2-10t+3~dt = (b-a)\cdot [\frac{8}{3}t^3-\frac{10}{2}t^2 + 3t]_0^1 \\&&= (b-a) \cdot \frac{2}{3} \end{eqnarray*}$

$\begin{eqnarray*} \alpha_1 &&= (b-a) \cdot \int_0^1~\frac{4t-1}{1-0} \cdot \frac{4t-3}{1-2}~dt =(b-a) \cdot \int_0^1~-16t^2+16t-3~dt = (b-a)\cdot [-\frac{16}{3}t^2+\frac{16}{2}t^2-3t ]_0^1 \\&&= (b-a) \cdot (-\frac{1}{3}) \end{eqnarray*}$

$\begin{eqnarray*} \alpha_2 &&= (b-a) \cdot \int_0^1~\frac{4t-1}{2-0} \cdot \frac{4t-2}{2-1}~dt =(b-a) \cdot \int_0^1~8t^2-6t+1~dt = (b-a)\cdot [ \frac{8}{3}t^2-\frac{6}{2}t^2+t]_0^1 \\&&= (b-a) \cdot \frac{2}{3} \end{eqnarray*}$

Quadraturformel

Als Näherung für das Integral $\begin{eqnarray*} \int_a^b~f(x)~dx = I(f) \end{eqnarray*}$ erhält man hier die Quadraturformel:

$\begin{eqnarray*} I_1(f)&&=\alpha_0 \cdot f(x_0) + \alpha_1\cdot f(x_1) + \alpha_2 \cdot f(x_2) =\\&&=(b-a)\cdot \Bigl\{\frac{2}{3} \cdot f(x_0) - \frac{1}{3} \cdot f(x_1) + \frac{2}{3}\cdot f(x_2) \Bigr\} \end{eqnarray*}$

16.04.2007, 00:44

tigerbine

Auf diesen Beitrag antworten »

3. Bezug zum Riemann-Integral
Hier soll das Integral für die Funktion $\begin{eqnarray*} f(x)=\cos(x) \end{eqnarray*}$ durch die oben aufgeführten Verfahren visualisiert werden. Hierbei soll der im Theorie-Teil erwähnte Ansatz der Treppenfunktion verdeutlicht werden.

Geschlossene NC-Formeln

[attach]8484[/attach][attach]8485[/attach][attach]8486[/attach]

Offene NC-Formeln

[attach]8487[/attach][attach]8488[/attach]

18.04.2007, 22:43

tigerbine

Auf diesen Beitrag antworten »

4. Vorgegebene Genauigkeit erreichen
Nun wollen wir uns der praktischen Aufgabe stellen, eine gegebene Funktion mittels der Bisher genannten Verfahren auf eine bestimmte Genauigkeit zu integrieren. Das erste i.A. konvergente Verfahren waren die summierten Newton-Cotes Formeln. Generell gilt es auch zu bedenken, dass das Auswerten einer Funktion (was ja immer nötig ist) auch mit Rechenaufwand verbunden ist und i.A. auch schon mit Fehlern belastet ist.

Beispielaufgabe:

Bestimme das Integral $\begin{eqnarray*} \int_0^1~f(x)~dx = \int_0^1~e^{-\frac{1}{2}x^2}~dx \end{eqnarray*}$ mit einer Fehlergenauigkeit von $\begin{eqnarray*} |e(x)| \leq 0.01 \end{eqnarray*}$

Sammeln wir zunächst ein paar Informationen über die Funktion.

Grob ergibt sich durch Abschätzung der Faktoren:

$\begin{eqnarray*} f(x) = e^{-\frac{1}{2}x^2}\Rightarrow \max_{x \in [0,1]}~|f(x)| = 1 \end{eqnarray*}$

$\begin{eqnarray*} f'(x) = e^{-\frac{1}{2}x^2} \cdot (-2x) \Rightarrow \max_{x \in [0,1]}~|f'(x)| < 2 \end{eqnarray*}$

$\begin{eqnarray*} f''(x) = e^{-\frac{1}{2}x^2} \cdot (4x^2-2) \Rightarrow \max_{x \in [0,1]}~|f''(x)| = 2 \end{eqnarray*}$

25.04.2007, 09:32

tigerbine

Auf diesen Beitrag antworten »

4a. Trapezregel
Aus dem Theorieteil wissen wir, dass es sich bei der Trapezregel (abgeschlossene NC-Regel mit gerader Anzahl von Stützstellen) um eine Quadraturformel der Ordnung 2 handelt, d.h. alle $\begin{eqnarray*} p \in \Pi_1 \end{eqnarray*}$ werden exakt integriert. Es sind2 Funktionsauswertungen nötig. Wie gut wird das Ergebnis sein? Werfen wir einen Blick in die Fehlerabschätzung:

$\begin{eqnarray*} |I(f)-I_1(f)|\leq \frac{(1-0)^3}{12}\cdot \max_{x \in[0,1]}~|f^{(2)}(x)| = \frac{1}{12} \cdot 2 = \frac{1}{6} \approx 0.167 \end{eqnarray*}$

Wir erreichen so (mit 2 Auswertungen) noch nicht die geforderte Genauigkeit. Daher wollen wir es mit der Summierten Trapezregel versuchen. An wie vielen Stellen muss dabei die Funktion (in Abhängigkeit von der Teilintervallanzahl N) ausgewertet werden?

$\begin{eqnarray*} I^T(f)=\frac{h}{2} \Bigl\{f(a) +2\cdot \sum_{p=1}^{N-1}~f(a+p\cdot h)+f(b) \Bigr\} \end{eqnarray*}$

Ab N>1 ist die Funktion also an (N+1)-Stellen auswerten. Wie können wir dann das Integral approximieren?

$\begin{eqnarray*} |I(f)-I^T(f)|\leq \frac{(1-0)}{12}\cdot h^2 \cdot \max_{x \in[0,1]}~|f^{(2)}(x)| \leq \frac{1}{6} \cdot h^2 \end{eqnarray*}$

Damit ergibt sich für die geforderte Genauigkeit:

$\begin{eqnarray*} 0.01 \leq \frac{1}{6} \cdot h^2 \end{eqnarray*}$

$\begin{eqnarray*} 0.06 \leq h^2 \end{eqnarray*}$

$\begin{eqnarray*} \sqrt{0.06} \leq h \end{eqnarray*}$

Für die Anzahl der Teilintervalle bedeutet dies dann:

$\begin{eqnarray*} N = \frac{1-0}{h} = 0.06^{-0.5} >4 \end{eqnarray*}$

D.h. mit N=5 sollte die gewünschte Genauigkeit erreicht worden sein.

[attach]8506[/attach]

25.04.2007, 09:32

tigerbine

Auf diesen Beitrag antworten »

4b. Ordnung der Quadraturformel
Die Potenz mit der h ist der Fehlerabschätzung auftaucht spielt also eine entscheidende Rolle. Was wäre, wenn wir nur eine Formel der Ordnung 1 hätten? Diese müssen wir uns nun erst noch basteln, da aufgrund der symmetrischen Lage der Knoten bei den NC-Formeln man immer Ordnugnen > 1 erhält. Wir wählen einfach die "Links-Regel" und setzten

$\begin{eqnarray*} I_L(f) = (b-a)\cdot f(a) \end{eqnarray*}$

Es ist dann:

$\begin{eqnarray*} \int_{a}^{b}~(x-a)~dx = [\frac{1}{2}x^2-ax]_a^b = \frac{1}{2}(b-a)^2 \end{eqnarray*}$

Somit erhalten wir

$\begin{eqnarray*} |I(f)-I_L(f)| \leq \frac{(b-a)^2}{2} \cdot \max_{x \in [a,b]} |f'(x)| \end{eqnarray*}$

Für die Summierte Formel ergibt sich dann:

$\begin{eqnarray*} |I(f)-I^{L}(f)| \leq N \cdot \frac{h^2}{2} \cdot \max_{x \in [a,b]} |f'(x)| = \frac{(b-a)}{2} \cdot h \cdot \max_{x \in [a,b]} |f'(x)| \end{eqnarray*}$

Analoge Rechnungen wie in 4a ergeben dann, dass N>16 sein um das Integral der Funktion mit der gewünschten Genauigkeit zu berechnen.

27.08.2008, 05:03

tigerbine

Auf diesen Beitrag antworten »

5. Adaptive Formeln
Ziel ist die adaptive Berechnung des Integrals. Beginnen wir mit der Bestimmung des ersten Teilintervalls, also [-1,?]. Zunächst einmal ein paar Vergleichswerte:

$\begin{eqnarray*} \int_{-1}^{1}~\frac{1}{10^{-4}+x^2}~dx = [100 \cdot \arctan(100x)]_{-1}^{+1} \approx 312.1593 \end{eqnarray*}$

$\begin{eqnarray*} \int_{-1}^{0}~\frac{1}{10^{-4}+x^2}~dx = [100 \cdot \arctan(100x)]_{-1}^{0} \approx 156.0797 \end{eqnarray*}$

$\begin{eqnarray*} \int_{-1}^{-0.5}~\frac{1}{10^{-4}+x^2}~dx = [100 \cdot \arctan(100x)]_{-1}^{-0.5} \approx 0.9998 \end{eqnarray*}$

1. Durchlauf

$\begin{eqnarray*} h_0 = 1 - (-1) = 2 \end{eqnarray*}$

$\begin{eqnarray*} P_0(f) = \frac{2}{6} \cdot \{f(-1) + 4 f(0) + f(1)\} \approx 13333.9999 \end{eqnarray*}$

$\begin{eqnarray*} Q_0(f) = \frac{2}{12} \cdot \{f(-1) + 4 f(-0.5) + 2 f(0) + 4f(0.5)+ f(1)\} \approx 3338.9978 \end{eqnarray*}$

$\begin{eqnarray*} \Rightarrow |Q_0(f)-P_0(f)| = 9995.0021 > \frac{15 \cdot 2 \cdot 0.1}{2} = 1.5 \end{eqnarray*}$

2. Durchlauf

$\begin{eqnarray*} P_1(f) = \frac{1}{6} \cdot \{f(-1) + 4 f(-0.5) + f(0)\} \approx 1669.4989 \end{eqnarray*}$

$\begin{eqnarray*} Q_1(f) = \frac{1}{12} \cdot \{f(-1) + 4 f(-0.75) +2 f(-0.5) + 4f(-0.25)+ f(0)\} \approx 840.0004 \end{eqnarray*}$

$\begin{eqnarray*} \Rightarrow |Q_1(f)-P_1(f)| = 829.4986 > \frac{15 \cdot 1 \cdot 0.1}{2} = 0.75 \end{eqnarray*}$

3. Durchlauf

$\begin{eqnarray*} P_2(f) = \frac{0.5}{6} \cdot \{f(-1) + 4 f(-0.75) + f(-0.5)\} \approx 1.0090 \end{eqnarray*}$

$\begin{eqnarray*} Q_2(f) = \frac{0.5}{12} \cdot \{f(-1) + 4 f(-0.75) +2 f(-0.5) + 4f(-0.25)+ f(0)\} \approx 1.0006 \end{eqnarray*}$

$\begin{eqnarray*} \Rightarrow |Q_1(f)-P_1(f)| = 0.0084 < \frac{15 \cdot 0.5 \cdot 0.1}{2} = 0.375 \end{eqnarray*}$

D.h. nach 3 Durchläufen ist für das Intervall [-1,-0.5] die geforderte Genauigkeit erreicht.

Mit einen kleinen Matlabfile, welches die Länge des letzten "erfolgreichen" Integrals für den neuen Durchlauf verdoppelt, erhält man:

code:

1:
2:
3:
4:
5:
6:
7:
8:
9:
10:
11:
12:
13:
14:
15:
16:
17:
18:
19:
20:
21:
22:
23:
24:
25:
26:
27:
28:
29:
30:

I_Matlab=3.121593e+002 
 
 i      [x_i,x_i+1]                  Q_i                  QS              m_i 
--------------------------------------------------------------------------------------
 0      [-1.000000,-0.500000]        1.000600e+000        1.000600e+000      2
 1      [-0.500000,-0.250000]        1.999794e+000        3.000394e+000      2
 2      [-0.250000,-0.125000]        3.988376e+000        6.988771e+000      2
 3      [-0.125000,-0.093750]        2.643535e+000        9.632306e+000      3
 4      [-0.093750,-0.062500]        5.239414e+000        1.487172e+001      1
 5      [-0.062500,-0.046875]        5.152812e+000        2.002453e+001      2
 6      [-0.046875,-0.031250]        9.952541e+000        2.997707e+001      1
 7      [-0.031250,-0.023438]        9.357928e+000        3.933500e+001      2
 8      [-0.023438,-0.015625]        1.660337e+001        5.593837e+001      1
 9      [-0.015625,-0.011719]        1.371128e+001        6.964966e+001      2
10      [-0.011719,-0.007813]        2.011666e+001        8.976631e+001      1
11      [-0.007813,-0.003906]        2.908028e+001        1.188466e+002      1
12      [-0.003906,-0.001953]        1.795142e+001        1.367980e+002      2
13      [-0.001953,+0.000000]        1.928844e+001        1.560864e+002      1
14      [+0.000000,+0.001953]        1.928844e+001        1.753749e+002      1
15      [+0.001953,+0.005859]        3.371303e+001        2.090879e+002      0
16      [+0.005859,+0.009766]        2.435244e+001        2.334404e+002      1
17      [+0.009766,+0.017578]        2.800233e+001        2.614427e+002      0
18      [+0.017578,+0.025391]        1.420411e+001        2.756468e+002      1
19      [+0.025391,+0.033203]        8.265489e+000        2.839123e+002      1
20      [+0.033203,+0.048828]        9.053302e+000        2.929656e+002      0
21      [+0.048828,+0.080078]        7.778394e+000        3.007440e+002      0
22      [+0.080078,+0.142578]        5.423128e+000        3.061671e+002      0
23      [+0.142578,+0.267578]        3.268479e+000        3.094356e+002      0
24      [+0.267578,+0.517578]        1.804851e+000        3.112404e+002      0
25      [+0.517578,+1.000000]        9.324854e-001        3.121729e+002      0

26.09.2008, 00:28

tigerbine

Auf diesen Beitrag antworten »

6. Romberg-Integration
Es soll das Integral $\begin{eqnarray*} \int_{1}^{2}~\frac{e^x}{x}~dx \end{eqnarray*}$ approximiert werden. Dazu bestimmt man zu verschiedenen Schrittweiten $\begin{eqnarray*} h_i:=\frac{(b-a)}{2^{i}} \end{eqnarray*}$ den Wert der summierten Trapezregel $\begin{eqnarray*} I^T(f)=\frac{h}{2} \cdot \Bigl\{f(a) +2\cdot \sum_{p=1}^{N-1}~f(a+p\cdot h)+f(b) \Bigr\} \end{eqnarray*}$

code:

1:
2:
3:
4:
5:
6:
7:
8:
9:
10:
11:
12:
13:
14:
15:
16:
17:
18:
19:
20:
21:
22:
23:
24:
25:
26:
27:
28:
29:
30:
31:
32:
33:
34:
35:
36:
37:
38:
39:
40:
41:
42:
43:
44:
45:
46:
47:
48:
49:
50:
51:
52:
53:
54:
55:
56:
57:
58:
59:
60:
61:
62:
63:
64:
65:
66:
67:
68:
69:
70:
71:
72:
73:
74:
75:
76:
77:
78:
79:
80:
81:
82:
83:
84:
85:
86:
87:
88:
89:
90:
91:
92:
93:
94:
95:
96:
97:
98:
99:
100:
101:
102:
103:
104:
105:

Es wird ein Integral durch Romberg-Extrapolation approximiert.
 
Beachte: Der Datensatz hat die Form
         Knoten:            h_0,...,h_n
         Integralwerte:     i_0,...,i_n
         Iterationen:       n
         Funktion in rombergf.m anlegen
Intervallanfang eingeben:         1
Intervallende eingeben:           2
Genauigkeit eps = 0.0001
 
--------------------------------------------------------------------------------------------
 
I_Matlab=3.059117e+000 
 
Neville-Romberg Schema  - Funktionswerte bei h=0
===============================================
 
IR =
    1.0000    3.2064    3.0607    3.0591    3.0591    3.0591
    0.5000    3.0971    3.0592    3.0591    3.0591         0
    0.2500    3.0687    3.0591    3.0591         0         0
    0.1250    3.0615    3.0591         0         0         0
    0.0625    3.0597         0         0         0         0
 
Fehler = 0.000000 
Neville-Romberg Schema - Polynomeinträge 
=========================================
 
Spalte  0 
-----------
 
p_ 0, 0(x)= 
 
     + 3.206405 * x^0 
 
p_ 1, 0(x)= 
 
     + 3.097099 * x^0 
 
p_ 2, 0(x)= 
 
     + 3.068704 * x^0 
 
p_ 3, 0(x)= 
 
     + 3.061520 * x^0 
 
p_ 4, 0(x)= 
 
     + 3.059718 * x^0 
 
Spalte  1 
-----------
 
p_ 0, 1(x)= 
 
     + 3.060663 * x^0     + 0.145741 * x^2 
 
p_ 1, 1(x)= 
 
     + 3.059239 * x^0     + 0.151439 * x^2 
 
p_ 2, 1(x)= 
 
     + 3.059125 * x^0     + 0.153267 * x^2 
 
p_ 3, 1(x)= 
 
     + 3.059117 * x^0     + 0.153767 * x^2 
 
Spalte  2 
-----------
 
p_ 0, 2(x)= 
 
     + 3.059144 * x^0     + 0.153338 * x^2     - 0.006077 * x^4 
 
p_ 1, 2(x)= 
 
     + 3.059117 * x^0     + 0.153877 * x^2     - 0.007803 * x^4 
 
p_ 2, 2(x)= 
 
     + 3.059117 * x^0     + 0.153934 * x^2     - 0.008532 * x^4 
 
Spalte  3 
-----------
 
p_ 0, 3(x)= 
 
     + 3.059117 * x^0     + 0.153913 * x^2     - 0.008379 * x^4     + 0.001754 * x^6 
 
p_ 1, 3(x)= 
 
     + 3.059117 * x^0     + 0.153938 * x^2     - 0.008775 * x^4     + 0.002961 * x^6 

Spalte  4 
-----------
 
p_ 0, 4(x)= 
 
     + 3.059117 * x^0     + 0.153938 * x^2     - 0.008801 * x^4     + 0.003363 * x^6     - 0.001212 * x^8

[attach]8709[/attach]

Neue Frage »

Antworten »

[WS] Numerische Integration - Beispiele

Verwandte Themen