Höhere Mathematik sehen und verstehen

Springer Spektrum ISBN 978-3-662-62576-7
www.mathematik-sehen-und-verstehen.de
https://masuv.web.leuphana.de
Dörte Haftendorn

Dieter Riebesehl

Hubert Dammer
Kapitel 05 Numerik
5.1 Fehler 5.2 Nullstellensuche 5.3 Interpolation, Gestaltung 5.4 B-Splines, NURBS 5.5 Numerische Integration 5.6 Numerische DGLn
ALLE BILDER (außer den Icons) erscheinen beim Anklicken in neuem Tab in dreifacher Größe.

5.1 Numerische Fehler
Abschnitt
5.1.1.1
Seite 351
Iteration für \(\pi\)
Seitenzahl: \(n=12\)
Segmenthöhe: \(k_n=0.5\)
Vieleckfläche: \(\displaystyle A_n=\frac{n}{2}k_n=3\)
Verdopplung:
\(k_{2n} = \) \(\displaystyle\sqrt\frac{1-\sqrt{1-k_n^2}}{2}\), \(A_{12\cdot2^m} = 3\cdot2^m\stackrel{\mbox{ \(m+1\) Wurzeln}}{\sqrt{2 - \sqrt{2+\sqrt{2+\dots\sqrt{3}}}}}\) \(\xrightarrow{m\to\infty}\pi\) für \(m\ge1\)

Zwölfeck.ggb
Abschnitt
5.1.1.1
Seite 351
Abb. 5.1 Iterative Berechnung von \(\pi\). In der Excel-Tabelle sind ein paar Zeilen ausgelassen. In der rechten Grafik ist aber für alle Zeilen der erreichte relative Fehler der Näherung für \(\pi\) aufgezeichnet. An der y-Achse stehen die dekadischen Logarithmen des relativen Fehlers, also z.B. \(-4\) für \(10^{-4}\).


ApproximationPi.xlsx
Abschnitt
5.1.1.1
Seite 352
Abb. 5.2 Stabile iterative Berechnung von \(\pi\).
\(n\) ist die Eckenzahl der regelmäßigen \(n\)-Ecke aus der 12-Eck-Familie, \(k_n\) sind die Höhen der Segmente, \(A_n\) ist die Gesamtfläche des \(n\)-Ecks, die im Grenzwert für \(n\to\infty\) gerade \(\pi\) sein muss. Nur der Anfang und der Schluss der Rechnung sind gezeigt.


ApproximationPi.xlsx

5.2 Nullstellensuche
Abschnitt
5.2.1
Seite 353
Abb. 5.3 Intervallschachtelungen. Bei allen Verfahren ist unten eine Intervallschachtelung angegeben, In a) und b) wird sie durch fortgesetzte Halbierung, Bisektionen, gefunden. c) zeigt die Lösung des Problems mit einer rekursiven Folge zur Fixpunktgleichung \(\cos x = x\), siehe Abschnitt 5.2.4.


bisektion.ggb
Abschnitt
5.2.2
Seite 354
Abb. 5.4 Sekantenverfahren. Neue Näherung \(x_2\) ist die Nullstelle der Sekante durch die Graphenpunkte an den Enden des Intervalls \([x_0,x_1]\). Mit einfacher Rechnung und \(y_i=f(x_i)\), \(i=0,1\), ergibt sich \[ x_2 =\frac{y_0 x_1-y_1 x_0}{y_0-y_1}. \] Das neue Intervall ist \([x_0,x_2]\) oder \([x_2,x_1]\), je nachdem, welches Vorzeichen \(f(x_2)\) hat. Man sieht hier, dass man für die Funktion \(f(x) = \cos x - x\) und Start bei \(x_0=0.2,\ x_1=1\), mit \(x_3=0.74\) schon 2-stellige Genauigkeit der Nullstelle erreicht hat.


sekanten-verfahren.ggb
Abschnitt
5.2.3
Seite 355
Abb. 5.5 Nullstellensuche mit dem Newtonverfahren für \(f(x) = {\rm e}^x-4x\). Es sind je zwei gekoppelte Grafikfenster nebeneinander dargestellt. Links ist jeweils der Funktionsgraph und die definierende Tangentenfolge, rechts die Trägerfunktion der Newton-Iteration mit der Folge der Näherungspunkte zu sehen. a) und b) zeigen ein "übliches" Beispiel, bei dem die Näherungsfolge monoton fallend auf die gesuchte Nullstelle zuläuft. c) und d) zeigen einen anderen Ausschnitt derselben Funktion. Damit erklären sich einige Überraschungen, die man mit der Newton-Iteration erleben kann. Der Text geht darauf in Abschnitt 5.2.4.2 ein.


newton_it-DH1.ggb
newton_it-DH2.ggb
newton_it-DH3.ggb
Abschnitt
5.2.4.1
Seite 357
Abb. 5.6 Wie schnell konvergiert das Newtonverfahren? Man misst die Konvergenzgeschwindigkeit daran, mit welchem Faktor die Anzahl der richtig berechneten Stellen von Schritt zu Schritt wächst. Hier werden die gelben Balken, die dies anzeigen, bei jedem Schritt etwa doppelt so lang.


Abschnitt
5.2.4.2
Seite 358
Abb. 5.7 CAS mit GeoGebra und Überraschung beim Newtonverfahren. a) zeigt die Berechnung der Trägerfunktion des Newtonverfahrens und die Berechnung von Folgenwerten. In b) und c) sind wieder zwei gekoppelte Grafikfenster nebeneinander dargestellt, aber diesmal für die Nullstelle von \({\rm e}^x-2.7x\). Betrachtet man nur berechnete Folgenglieder, ohne auf einen Graphen zu achten, könnte man meinen, die Folge konvergiere gegen einen Wert zwischen \(x_3\) und \(x_6\). In c) erkennt man aber, dass die rekursive "Treppchenfolge" überhaupt nicht konvergieren wird. Ein Zoom im Fenster b) hätte auch gezeigt: Es gibt gar keine Nullstelle.


newton-it-cas1.ggb
newton-it-cas2.ggb
newton-it-cas4.ggb
Abschnitt
5.2.5.1
Seite 359
Abb. 5.8 Das Heronverfahren als Rekursion. a) zeigt eine Parabel mit der Nullstelle \(\sqrt{r}\), die mit dem Newtonverfahren approximiert wird. b) verfolgt mit dem Startwert \(x_0\) die Heronfolge im gekoppelten Fenster in der iterativen Treppendarstellung. Wenn Sie in GeoGebra \(x_0\) oder den Radikanden \(r\) verändern, reagieren beide Darstellungen gemeinsam. In Bild b) ist die Trägerfunktion eine Hyperbel.


heron-newton_it.ggb

5.3 Interpolation und Kurven gestalten
Abschnitt
5.3.1.2
Seite 362
Abb. 5.9 Interpolation nach Lagrange. a) zeigt Lagrange-Basis-Polynome, die an genau 3 Stützstellen Nullstellen haben. In b) wird z.B. der Faktor \(c_b\) für \(L_B(x)\) gefunden, der \(L_B\) so staucht, dass \(B\) auf dem Graphen liegt. Das rote Interpolationspolynom entsteht auf diese Weise.


LagrangePolynome.ggb
Abschnitt
5.3.1.2
Seite 363
Abb. 5.10 Kostenfunktion als Lagrange-Polynom, aus ihr wird alles, was noch dargestellt ist, berechnet. Erfahren Sie durch Ziehen an den gegebenen Punkten, wie überraschend deutlich die Wirtschaftsaussagen, die sich daraus ergeben, auf kleine Änderungen der Lage reagieren.


lagrange_wirtschaft3.ggb
Abschnitt
5.3.1.3
Seite 364
Abb. 5.11 Das Newton'sche-Interpolations-Polynom ist eine Linearkombination der Newton'schen Basispolynome, die aus immer mehr Linearfaktoren \((x-\mbox{Stützstelle})\) aufgebaut sind und mit den passenden Koeffizienten immer mehr Stützpunkte genau treffen.


newton_interpolation.ggb
Abschnitt
5.3.2
Seite 365
Abb. 5.12 Gauß'sche Glockenkurve und Interpolationspolynom. a) Die Punkte \(A,\,\dots,\,G\) liegen auf der Glockenkurve, das blaue Interpolationspolynom passt nach Sicht brauchbar in deren Bereich. b) Wie gut die Passung wirklich ist, sieht man durch Visualisierung der Differenzen zur Glockenkurve. Für verschiedene Stellungen von \(F\) zwischen \(A\) und \(G\) ist je eine solche "Fehlerkurve" gezeichnet.


interpol-gauss-deutlich.ggb
Abschnitt
5.3.2.1
Seite 366
Abb. 5.13 Poeler Kogge und Splines. Einsatz von Splines (drei durchgezogene bunte Stücke) für die Balkenform der Kogge, die Nägel sind in den Punkten \(A,B,C,D\) zu denken.
Leider kann man die Punkte (noch) nicht ziehen, weil wir die Koeffizienten aus einer Berechnung mit Mathematica übertragen haben. Später kann man dann selbst ein anderes Bild hinterlegen.

poelerkogge3.ggb
Abschnitt
5.3.2.2
Seite 367
Abb. 5.14 Basisfunktionen für kubische Splines. Links: \(b_0,\ b_3\), und 6-fach überhöht \(b_1,\ b_2\). Rechts: Ein damit erstellter natürlicher Spline durch 5 Nägel, in vier Farben durchgezogen dargestellt. Die Nägel sind frei verschiebbar. Zum Vergleich ist das Interpolationspolynom schwarz gestrichelt dazu gezeichnet. Man sieht, dass der Spline wesentlich weniger "ausschwingt". Orange gestrichelt ist das Ergebnis des GeoGebra-Befehls Spline(A,B,C,D,E), siehe dazu den Text.


BasisKubSplinesNatuerl.ggb
SplineKub5Punkte.ggb
Abschnitt 5.3.2.2 S. 368 Hier finden Sie die Ableitungen der Spline-Basisfunktionen und die Herleitung der Gleichungen zur Splineberechnung.
Zunächst tabellarisch die Basisfunktionen und deren Ableitungen an den Intervallenden 0 und 1:

\[ \begin{array}{|l|l|l|l|} \hline b_0(x) = 1 - x & b'_0(x) = -1 & b'_0(0) = -1 & b'_0(1) = -1 \\ \hline b_1(x) = \frac{1}{6}(1-x)(x-2)x & b'_1(x) = \frac{1}{6}\left((2-x)x + (1-x)x + (1-x)(2-x)\right) & b'_1(0) = -\frac{1}{3} & b'_1(1) = \frac{1}{6} \\ \hline b_2(x) = \frac{1}{6}(x-1)(x+1)x & b'_2(x) = \frac{1}{6}\left((x+1)x + (x-1)x + (x-1)(x+1)\right) & b'_2(0) = -\frac{1}{6} & b'_2(1) = \frac{1}{3} \\ \hline b_3(x) = x & b'_3(x) = 1 & b'_3(0) = 1 & b'_3(1) = 1 \\ \hline \end{array} \] Die im Punkt \(B\) zusammentreffenden Splinepolynome sind mit den Bezeichnungen aus dem Buch

\( p_{AB} = A_y b_0^{AB} + c_A b_1^{AB} + c_B b_2^{AB} + B_y b_3^{AB} \quad \) und \( \quad p_{BC} = B_y b_0^{BC} + c_B b_1^{BC} + c_C b_2^{BC} + C_y b_3^{BC} \).

Die Bedingung einer stetigen 1. Ableitung im Punkt \(B\) führt zu \( p_{AB}'(B) = p_{BC}'(B).\) Hier setzt man die oben angegebenen Ableitungen ein und nutzt aus, dass \( b_i^{BC} \) einfach das nach rechts verschobene \( b_i^{AB} \) ist:

\( A_y b_0'(1) + c_A b_1'(1) + c_B b_2'(1) + B_y b_3'(1) = B_y b_0'(0) + c_B b_1'(0) + c_A b_2'(0) + C_y b_3'(0)\)

\( -A_y + \frac{1}{6}c_A + \frac{1}{3}c_B + B_y = -B_y - \frac{1}{3}c_B - \frac{1}{6}c_A + C_y \)

\( \frac{1}{6}c_A + \frac{2}{3}c_B + \frac{1}{6}c_C = C_y - 2B_y + A_y\)

\( \frac{1}{6}(c_A + 4c_B + c_C) = C_y - 2B_y + A_y\)

Das ist die erste der Gleichungen, die anderen ergeben sich daraus durch Verschiebung zu den nächsten Punkten.

Abschnitt
5.3.3
Seite 369
Abb. 5.15 Bézierkurve a) zusammen mit ihren Steuerpunkten \(A,B,C,D\), in b) erzeugt aus ihrem "Gerüst": Man wählt einen Prozentsatz \(t\) zwischen 0 und 1, hier \(t=0.4=40\%\), und markiert auf jedem der schwarzen Gerüstvektoren den \(40\%\)-Punkt von \(A\), von \(B\) und von \(C\) aus. Auf den zwei sich daraus ergebenden Vektoren (orange) werden wieder die \(40\%\)-Punkte markiert. Im nächsten Schritt (rosa) hat man so den Punkt \(P\) der Bézierkurve erzeugt. Diese ist der geometrische Ort von \(P\) bei jeder Wahl von \(t\). Siehe auch Abb. 5.16 und Seite 372.


bezier_Geruest.ggb
Abschnitt
5.3.3.1
Seite 370
Abb. 5.16 Bézierkurve als Parameterkurve mit Bernstein-Polynomen. a) zeigt den Bézierspline, mit dem man einen Notenbogen anpassen kann, siehe Abb. 5.15 und Seite 372. b) Die vier Bernstein-Polynome sind eine Basis für den Polynomraum \(\Pi_3\), wir verwenden sie im Intervall [0,1]. c) zeigt, wie die Koordinaten von \(B\) als Koeffizienten von \(b_1\) in der Parameterkurve wirken.


bezier_parameterkurve-summanden-buch.ggb
Abschnitt
5.3.3.2
Seite 371
Abb. 5.17 Bézierkurve als Parameterkurve.
Allgemeine Erklärung für Parameterkurven.
Die Koordinaten eines Kurvenpunktes \(P\) sind Funktionen \(x(t),\ y(t)\) eines Parameters \(t\). Diese beiden Funktionen sind hier zusätzlich dargestellt, und zwar so, dass die Parameterkurve durch die waagerechten und senkrechten gestrichelten Linien geometrisch entstehen kann. Sie können an \(t\) ziehen, damit \(P\) wandern lassen, aber auch das "Gerüst" verändern.


Bezier anschaulich.ggb
Abschnitt
5.3.3.3
Seite 372
Abb. 5.18 Krümmung des Béziersplines in \(A\), dargestellt durch den blauen Krümmungskreis, mit seiner Spur in GeoGebra. a) zeigt, dass sich bei Verlängerung von \(\overrightarrow{AB}\) die Krümmung ändert. In b) ändert sich die Höhe \(k\) des Dreiecks \(ABC\) nicht, wenn \(C\) auf der Parallelen zu \(AB\) wandert, also bleibt die Krümmung gleich. c) zeigt beispielhaft, dass eine beliebige Veränderung von \(D\) keinen Einfluss auf die Krümmung in \(A\) hat.


Bezier-kruemmung.ggb
Abschnitt 5.3.3.3
Seite 372
Die Krümmungsformel für den Beziérspline im Punkt \(A\) lautet \[ \kappa = \frac{2k}{3\Big|\Big|\overrightarrow{AB}\Big|\Big|} \] mit der Höhe \(k\) auf der Seite \(AB\) im Dreieck \(ABC\). Die Herleitung benutzt die Formel \[ \kappa = \frac{\ddot y\dot x - \ddot x \dot y}{(\dot x^2 + \dot y^2)^{\frac{3}{2}}}. \] Sie ist anzuwenden auf die Parameterdarstellung \[ \vec P = \begin{pmatrix} P_x \cr P_y\end{pmatrix},\quad \begin{array}{cccccc} &P_x(t)=&A_x b_0(t)&+\;B_x b_1(t)&+\;C_x b_2(t)&+\;D_xb_3(t) = x(t),\\ &P_y(t)=&A_y b_0(t)&+\;B_y b_1(t)&+\;C_y b_2(t)&+\;D_yb_3(t) = y(t), \end{array} \] oder ausgeschrieben \[ \vec P = \begin{pmatrix} x(t) \cr y(t)\end{pmatrix} = \vec A (1-t)^3 + 3\vec B (1-t)^2 t + 3\vec C (1-t)t^2 + \vec D t^3. \] Wir rechnen die nötigen Ableitungen im Punkt \(A\), also für \(t=0\) aus. Dabei nehmen wir von vornherein nur Terme mit, die keinen Faktor \(t\) enthalten, da sie ohnehin wegfallen. Bei der Anwendung der Produktregel schreiben wir also nur die Terme auf, in denen alle \(t\)-Faktoren "wegdifferenziert" wurden: \[ \begin{pmatrix} \dot x(t) \cr \dot y(t)\end{pmatrix}\Bigg|_{t=0} = \Big(-3\vec A (1-t)^2 + 3\vec B(1-t)^2\Big)\Big|_{t=0} = 3\vec B - 3\vec A. \] Vorsicht: für die zweiten Ableitungen müssen wir wieder auf \(\vec P\) schauen, weil nun auch \(\vec C\) einen Beitrag leisten wird: \[ \begin{pmatrix} \ddot x(t) \cr \ddot y(t)\end{pmatrix}\Bigg|_{t=0} = \Big(6\vec A (1-t) - 12\vec B(1-t) + 6\vec C(1-t)\Big)\Big|_{t=0} = 6\vec C - 12\vec B + 6\vec A. \] \(\vec B\) hat einen Faktor 12 statt 6, weil wir bei diesem Summanden bei Anwendung der Produktregel offenbar einmal den \(t\)-Faktor und einmal den \((1-t)\)-Faktor ableiten müssen, aber dafür gibt es zwei Reihenfolgen, die bei der zweimaligen Anwendung der Produktregel auch beide vorkommen. Wenn Sie diesem Argument nicht folgen können, dann müssen Sie alles "zu Fuß" nachrechnen.

Nun ist alles in die Krümmungsformel einzusetzen: \[ \kappa = \frac{3(B_x-A_x)\cdot6(A_y-2B_y+C_y) - 3(B_y-A_y)\cdot6(A_x-2B_x+C_x)} {\big( 9(B_x-A_x)^2 + 9(B_y-A_y)^2 \big)^{\frac{3}{2}}} = \frac{18}{27}\frac{B_xC_y + A_xB_y - A_xC_y - B_yC_x - A_yB_x + A_yC_x} {\Big|\Big|\overrightarrow{AB}\Big|\Big|^3}. \] Jetzt kommt die Dreieckshöhe \(k\) ins Spiel. Die Fläche des Dreiecks \(ABC\) errechnet sich zu \(\frac{1}{2}k\Big|\Big|\overrightarrow{AB}\Big|\Big|\), das Doppelte davon ist aber gerade die Fläche des von \(\overrightarrow{AB}\) und \(\overrightarrow{BC}\) aufgespannten Parallelogramms. Letztere lässt sich aber auch als \(\Big|\overrightarrow{AB}\times\overrightarrow{BC}\Big|\) schreiben. Dieses Kreuzprodukt hat einen Wert, der sich als Determinante mal einen Einheitsvektor \( \vec e_{\bot} \)schreiben lässt, welcher auf der Ebene, die von \(ABC\) aufgespannt wird, senkrecht steht: \[ \overrightarrow{AB}\times\overrightarrow{BC} = \begin{vmatrix} B_x-A_x & C_x-B_x \cr B_y-A_y & C_y-B_y \end{vmatrix}\cdot\vec e_{\bot}, \] was ausmultipliziert genau den obigen Zähler für \(\kappa\) ergibt! Damit hat man \[ \kappa = \frac{2}{3}\frac{k\Big|\Big|\overrightarrow{AB}\Big|\Big|} {\Big|\Big|\overrightarrow{AB}\Big|\Big|^3} = \frac{2}{3}\frac{k} {\Big|\Big|\overrightarrow{AB}\Big|\Big|^2} \]

Abschnitt
5.3.4
Seite 373
Abb. 5.19 Bézierfläche, erzeugt aus ihrem "Gerüst":
a)
Realisierung in Mathematica.
b) Für die Realisierung in GeoGebra nennt man im Tabellenfenster die z-Ordinaten für die Stützpunkte. Daraus wird eine Matrix \(ordi\) erzeugt, auf deren Komponenten man mit Indizes zugreifen kann. In griffiger Mathematica-Syntax ist in b) unten gezeigt, wie daraus mit den Bernsteinpolynomen die z-Komponente der Parameterdarstellung erzeugt wird. Sie ist ein Polynom 6. Grades in s und t.
c) Realisierung im 3D-Fenster von GeoGebra.
zu a) BezierFlaechen.nb, Mathematica Quelltext
Dies zum Lesen und Verstehen
 
zu b) und c) bezier-flaeche+geruest.ggb
Anmerkung: Die Funktion Oberfläche(s,t,z(s,t),s,0,3,t,0,3) geht stets verloren beim Öffnen. Laden Sie die ggb herunter und öffnen Sie sie dann von Hand im geöffneten GeoGebra. Kopieren Sie von hier    1+s-s^2/3-s^3/27+t+t^2/3-(s t^2)/3+(2 s^3 t^2)/81-(4 t^3)/27+(2 s t^3)/27-(s^2 t^3)/81-(2 s^3 t^3)/729
Tragen Sie diesen Term als z(s,t) in den Oberfläche-Befehl in die Eingabezeile ein. Dann entsteht Abb. 5.19 c).
Leider ist dies nur statisch. Wenn Sie aber stpoly (steht im Algebafenster) anstelle von z eintragen, können Sie mit unserer Datei können Sie dynamisch experimentieren.
Leider bleibt dies nicht beim Speichern erhalten. Die lange Summe f von oben, aber auch die für stpoly gemäß Abb. 5.19 b), lässt sich in GeoGebra schlecht eintippen. Man erzeugt sie besser in einem anderen Editor (wir nahmen Mathematica) und überträgt das mit Copy und Paste.
Wir haben dieses merkwürdige Verhalten dem GeoGebra-Team berichtet.

5.4 B-Splines und NURBS
Abschnitt
5.4
Seite 374
Abb. 5.20 Legende Bernsteinpolynome 3. Grades, wie sie für Béziersplines verwendet werden. a) als Graphen, b) gestapelt, ihre Summe ist \(=1\) für jedes \(t\) im Intervall \([0,\,1]\), c) Erweiterung mit drei weiteren Steuerpunkten und einer weiteren Bézierkurve: die Gesamtkurve ist i.A. an der Nahtstelle nicht differenzierbar.


bernsteinSumme.ggb
Abschnitt
5.4.2.1
Seite 376
Abb. 5.21 Aufbau der Basen für B-Splines. a) oben zwei Basissplines zu \(p=0\), darunter die Gewichtungsfunktionen in passender Farbe: orange zu rot und hellblau zu blau, aus diesen errechnet sich in b) die erste Basisfunktion zu \(p=1\), dazu die zweite in blau/blau gestrichelt, darunter wieder die Gewichtsfunktionen, aus diesen errechnet sich in c) die erste Basisfunktion zu \(p=2\) usw. bis man in d) Basisfunktionen zu \(p=3\) sieht.


B-Spline-Basis2.ggb
Abschnitt
5.4.2.2
Seite 378
Abb. 5.22 Elemente eines langen B-Splines für 18 Punkte (\(A\) bis \(R\)).
a) Basiskurven vom Grad 3, die für 18 Punkte benötigt werden. Sie sind jeweils in einem Intervall der Breite 4 definiert.
b) 18 Punkte, die das Gerüst für den B-Spline bilden, und der B-Spline für \(3\le t\le 18\).
c) Die Graphen von \(P_x(t)=A_x B_{0,3}(t)+\dots+ R_x B_{17,3}(t)\) in Grün und \(P_y(t)\) entsprechend in Blau.


B-Splines3Grad.ggb
B-Splines3Grad-PxPy-lang.ggb
Abschnitt
5.4.3.1
Seite 379
Abb. 5.23 Hinführung zu NURBS.
a) Nicht-uniforme Basiskurven vom Grad 3, je fünf benachbarte Knoten aus der Liste \(T\) gehören zu einer Basis-Funktion 3. Grades für B-Splines.
b) Mit Gewichten und rationalen Basen kann man u.a. auch exakte Standardkurven mit NURBS modellieren, es ist hier gezeigt für den rot-grün-gestrichelten Viertelkreis, der mittlere der fünf Bögen. Variation eines Gewichtes erzeugt die anderen Kurven.
  c) Basiskurven und d) ein geschlossener NURBS daraus mit 7 Steuerpunkten und nicht uniformen, rationalen B-Splines zweiten Grades, der aus vier Ellipsenstücken besteht. Hier in GeoGebra frei ziehbar.

a) und b) NURBS1 zum Viertelkreis in vielfätiger Datei, Splines-und-NURBS.nb Mathematica Quelltext    Dies zum Lesen und Verstehen      c) und d) NURBS2.ggb
Zusatz zu
Abschnitt
5.4.3.1
Seite 379
Trisektrix als rationaler Bézierspline
Die Grundidee ist es, eine Kurve dritten Grades mit rationalen Béziersplines darzustellen, denn wir wollen die Bernsteinpolyome \(b_i\) verwenden, die ja Grad 3 haben. Wir haben uns für die Trisektrix entschieden, von der wir aus dem Kurvenbuch Seite 64 die implizite Gleichung \((a+x)y^2=(3a-x)x^2\) kennen.
Nun benötigen wir eine rationale Parameterdarstellung.
Eine solche bekommt man wenn man durch eine Singularität, also hier den Doppelpunkt im Ursprung, eine Gerade \(y=tx\) legt, sie mit der Kurve zum Schnitt bringt und t als Parameter wählt. Sofort folgt \((a+x)t^2=(3a-x)\) und damit
\(x(t)=\frac{a(3-t^2)}{1+t^2}\). Aus der Geradengleichung folgt dann \(y(t)=\frac{at(3-t^2)}{1+t^2}\).
Die rationale Basis, die wir verwenden wollen, ist nach Seite 379 unten \(R_i=\frac{w_i b_i}{\Sigma_{j=0}^3{w_j b_j}}\), für \(i=0\dots3\).
In der Linearkombination der \(R_i\) mit den Steuerpunktkoordinaten als Koeffizienten gilt der gemeinsame Nenner \(\Sigma_{j=0}^3{w_j b_j}=1+t^2\) ohne Einfluss der Steuerpunkte. Durch Koeffizientenvergleich der Potenzen von \(t\) bestimmen wir \(w_0=1,\;w_1=1,\;w_2=\frac{4}{3},\;w_3=2\). Beim Kreis (s.o.) ist es (zufällig) derselbe Nenner.
Mit diesen Gewichten im Zähler ergibt sich, wie im Buch Seite 380, die rationale Basis:
\(R_0(t)=\frac{(1-t)^3}{t^2+1},\;R_1(t)=\frac{3 (1-t)^2 t}{t^2+1},\;R_2(t)=\frac{4 (1-t) t^2}{t^2+1},\;R_3(t)=\frac{2 t^3}{t^2+1}\)
Die Steuerpunkte müssen nun so gewählt werden, dass der Gesamtzähler der Linearkomination \(A\,R_0+B\,R_1+C\,R_2+D\,R_3\) für \(x(t)\) und \(y(t)\) der Zähler der Parameterdarstellung der Trisektrix ist:
\( a(3-t^2)=A_x(1-t)^3+B_x 3 (1-t)^2 t+C_x4 (1-t) t^2+D_x 2 t^3\)
\( at(3-t^2)=A_y(1-t)^3+B_y 3 (1-t)^2 t+C_y4 (1-t) t^2+D_y 2 t^3\)
Durch Koeffizientenvergleich erhält man aus einem Gleichungssystem
\(A=(3 a,0),B=(3 a,a),C=\left(2 a,\frac{3 a}{2}\right),D=(a,a)\) als Steuerpunkte einer Trisektix mit dem Doppelpunkt im Ursprung und der Schlaufenbreite 3\(a\).
Trisektrix-rat-Bezierspline.nb, Mathematica Quelltext    
Dies zum Lesen und Verstehen
Zusatz zu
Abschnitt
5.4.3.1
Seite 380

Durchlauf des blauen Kreises für t von 0 bis 1
Bild a) mathem. negativ, Uhrzeigersinn,
Bild b) mathem. positiv, entgegen dem Uhrzeigersinn,
Bild c) mathem. negativ, Uhrzeigersinn, aber unpassend zugeordnet.

Im vorigen Absatz ist a) animiert zu sehen. Die "Umwendung" hat ihren Grund in der gegenläufigen Parametrisierung.
Metamorphose einer Trisektrix in einen Kreis mit Hilfe von rationalen Béziersplines
Zwei Parametrisierungen des Kreises werden hier verwirklicht und zwar genau nach dem eben beschriebenen Verfahren. Daher sind die beiden folgenden Mathematica-Realsierungen lediglich winzige Modifizierungen der oben genannten Datei.
Kreis-Buch-als-rat-Bezierspline.nb, Mathematica Quelltext
Dies zum Lesen und Verstehen
Kreis-alternativ-als-rat-Bezierspline.nb, Mathematica Quelltext
Dies zum Lesen und Verstehen
Parameterdarstellung, Kreis aus dem Buch Seite 379
\(x(t)=\frac{2rt}{1+t^2}\).     \(y(t)=\frac{r(1-t^2)}{1+t^2}\).
Parameterdarstellung, Kreis anders herum durchlaufen
\(x(t)=\frac{r(1-t^2)}{1+t^2}\).     \(y(t)=\frac{2rt}{1+t^2}\).
Die beiden NURBS vom Kreis und die Trisektrix unterscheiden sich also nur durch die Stellung der Steuerpunkte. Darum verwirklicht die GeoGebra-Datei eine gleichzeitige lineare Bewegung von einem Steuerpunkt zum entsprechenden anderen, wie in den Bildern angedeutet. Wir erhalten bei a) und b) eine Metamorphose von der Trisektrix zum Kreis. Bei c) ist die Zuordnung wie bei b) aber der Kreis wie bei a). Damit gelingt die Verwandlung nicht.
a) echte-Trisektrix-Kreis.ggb \(A\rightarrow P_0,\;B\rightarrow P_1,\;C\rightarrow P_2,\;D\rightarrow P_3,\;\), die \(P_i\) gemäß dem Buch
b) echte-Trisektrix-anderer-Kreis.ggb
c) Trisektrix-Kr-Buch-o-Wenden.ggb
Zusatz zu
Abschnitt
5.4.3.3
Seite 380
NURBSKnotenDoppelt28.jpg
NURBSKnotenDoppelt28.jpg
Doppelte Knoten in NURBS, Erläuterungen
Im Buch steht auf Seite 380 zu den geschlossenen NURBS, dass als nötiger Trick Elemente in der Knotenliste verdoppelt werden, und dass man die dadurch in der Rekursion entstehenden Intervalle der Länge 0 einfach weglassen soll. Als Erklärung wird nur geboten, dass das "tatsächlich funktioniert". Wir wollen hier aber begründen, dass dieses Weglassen der richtige Weg ist.

NurbsKnotendoppelt.ggb     

In dieser GeoGebra-Datei werden zur Knotenliste {0,1,2,3,4} die Basispolynome \(B_{00},\ B_{10},\ B_{20},\ B_{30},\ B_{01},\ B_{11},\ B_{21},\ B_{02},\ B_{12}\) und \(B_{03}\) eines NURBS berechnet. Sie können in der GeoGebra-Datei links im Algebrafenster die Funktionen auswählen, die sie sehen wollen. Dort können Sie auch rekursiven Definitionen verfolgen und die Funktionsterme sowie die Graphen ansehen. Links sind durchgezogen\( B_{02},\ B_{03}\) und \(B_{12}\) zu sehen, die gestrichelten Graphen werden im Folgenden erklärt.



Der Knoten 2 ist über einen Schieberegler \(a\) veränderbar, er kann im Bereich 2 bis 3 verschoben werden. Das heißt, dass Sie den Knoten 3 verdoppeln können. Dabei können Sie sich anschauen, wie sich der Verlauf der Funktionen ändert. Bei \(B_{11}\) werden Sie den Bruch \(\frac{3-x}{3-a}\) als Teilterm entdecken. Dieser ist für \(a=3\), also den verdoppelten Knoten, nicht definiert. Ziehen Sie \(a\) auf 3, wird \(B_{11}\) nicht mehr angezeigt, und alle vom Intervall [a,3] abhängigen Funktionen auch nicht. Deshalb gibt es die Funktionsversionen \(B_{10x},\ B_{01x},\ B_{11x},\ B_{21x},\ B_{02x},\ B_{03x}\) und \(B_{12x}\) (alle gestrichelt), die undefinierte Teilterme weglassen und \(a=3\) setzen. \(B_{20}\) wird nicht mehr gebraucht. Im ersten Bild links oben sind \(B_{02},\ B_{03},\ B_{12},\ B_{02x},\ B_{03x}\) und \(B_{12x}\) gezeigt für \(a=2.8\), im Bild links Mitte nur \(B_{02x},\ B_{03x}\) und \(B_{12x}\) für \(a=3\). Sie können sehen, wie sich die durchgezogenen Funktionen an die gestrichelten annähern, wenn \(a\to3\) geht.
NURBSKnotenDoppelt28.jpg Wir wollen zum Verstehen durch Anschauen aber auch noch einen rechnerischen Beweis liefern, und zwar für \(B_{11}\), welche Funktion ja sozusagen der "Anfang allen Übels" ist. Dazu können Sie sich auch die Funktion \(\Delta B_{11} = B_{11} - B_{11x}\) anzeigen lassen, wie es im unteren Bild für \(a=2.8\) geschehen ist. Es ist \[ \Delta B_{11} = \left(\frac{x-1}{a-1} - \frac{x-1}{3-1}\right)B_{10} + \left(\frac{3-x}{3-a} - \frac{x-1}{3-1}\right)B_{20}= \] \[ = (3-a)\frac{x-1}{(a-1)(3-1)}B_{10} + \text{ irgendwas }\cdot B_{20}. \] Der linke Summand geht gegen 0, wenn \(a\to3\) geht, der rechte ist nur auf dem Intervall [a,3] von null verschieden. Daraus folgt, dass \(\Delta B_{11}\to0\) punktweise geht für \(a\to3\), oder anders gesagt \(B_{11}\stackrel{a\to3}\longrightarrow B_{11x}\) punktweise.
Eine analoge Rechnung gelingt in allen anderen Fällen ebenso. Schauen Sie sich auch an, wie die Grenzfunktionen für \(a=3\) im doppelten Knoten einmal weniger differenzierbar sind als ohne doppelten Knoten (0-mal differenzierbar heißt dabei einfach "stetig") und überlegen Sie, dass das daran liegt, dass die Gewichtsfunktionen am Intervallrand am doppelten Knoten keine Nullstelle mehr haben: der dazu nötige Term wird gerade weggelassen!
Abschnitt
5.4.3.4
Seite 381
Abb. 5.24 Torus als NURBS. In der Praxis gibt es ausgereifte Software für NURBS, die dem Anwender die Mühen abnimmt. Das gilt insbesondere für den Einsatz im 3D-Bereich. Im Gegensatz zu einfacheren Konzepten für Splines können NURBS sowohl Freiformen als auch Standardformen wie Zylinder, Kegel, Ringe, Kugel u.s.w. gestalten.

NURBS-Torus.nb, Mathematica Quelltext
Dies zum Lesen und Verstehen

5.5 Numerische Integration
Abschnitt
5.5.1.1
Seite 383
Abb. 5.25 Beweise der Kepler'schen Regel:
a) beweist mit der linken Parabel eine Sondersituation, die dann rechts daneben verwendet wird.
b) zeigt den historischen, geometrischen Weg von Archimedes und Kepler. Ein Parabelsegment (grün) nimmt immer \(\frac{2}{3}\) des umbeschriebenen (blau gestreiften) Parallelogramms ein.

kepler-beweise.ggb
Abschnitt
5.5.1.1
Seite 384
Abb. 5.26 Allgemeine Aussage über Parabelsegmente, bewiesen mit Integralrechnung und Scherung. Näheres dazu im Text.


parabelsegment.ggb
Abschnitt
5.5.1.2
Seite 385
Abb. 5.27 Kepler'sche Regel für eine Parabel und ein Polynom $p_3$. a) Die beiden Integralflächen haben zwar nicht dieselbe Form, aber dieselbe Fläche. b) Der Grund ist, dass das Integral über ihre Differenzfunktion null ist. c) Kleines Lehrstück zur Bestimmung der Interpolationspolynome mit in GeoGebra definierten Matrizen (siehe Text und den folgenden Zusatz XXX).


kepler1.ggb
kepler2.ggb
Abschnitt
5.5.2
Seite 386
Abb. 5.28 Simpson'sche Regel an Beispielen.
a) zeigt die Simpson-Regel mit 4 Streifen der Breite \(h=0.25\) für die Funktion \(\sin(\pi x^2)\).
In b) wird die Viertelkreisfläche mit ebenfalls 4 Streifen der Breite 1 angenähert, aber mit bescheidenem Erfolg.
Genaueres im folgenden Zusatz
c) zeigt eine Näherung für die wichtige Gauß'sche Glockenkurve mit 8 Streifen, der Fehler liegt unter 3 Millionstel.
Weiteres zu den gezeigten Näherungen steht im Abschnitt 5.5.2.2 nach der Simpson-Formel in Satz 5.5 und der Fehlerbetrachtung 5.5.2.1.
Zusatz
zum
Teil b)
simpson-a.ggb Sinus      simpson-b.ggb Viertelkreis     simpson-c.ggb Glockenkurve
Die Breite zweier benachbarter Balken sollte so gewählt werden, dass man das zugehörige Kurvenstück annähernd als Parabel (oder Poymom 3. Grades) auffassen kann.

Bei dem Viertelkreis in b) sind die zugehörigen Parabeln eingezeichnet. Es ist daher klar, dass der linke Doppelbalken einen deutlichen Fehler erzeugt.
In der rechten Bildhälfte ist GeoGebra-CAS zu sehen. Es werden damit die Parabeln berechnet, die die Grundlage der Rechnung sind.
Achten Sie in simpson-b.ggb Viertelkreis auf das CAS-Fenster. Man macht es bei Ansicht mit Strg+Umschalt+K (zweimal diesen Hotkey klicken) auf.
Abschnitt
5.5.3.1
Seite 392
Abb. 5.29 Fehler einfacher Quadraturformeln. Es ist in jedem Bild die approximierte Integralfläche im Vergleich mit der wirklichen Integralfläche zu sehen. Die Fehlerflächen sind nach links zusammengeschoben und liegen auf einem Streifen der Breite \(h\) bzw. bei d) \(\frac{h}{2}\). Einzelheiten stehen im Text. In der GeoGebra-Datei können Sie mit \(n\) die Anzahl der Streifen verändern und das Fehlerverhalten verstehen.


FehlerObersumme.ggb FehlerUntersumme.ggb
FehlerTrapezsumme.ggb FehlerMittenregel.ggb
Abschnitt
5.5.3.2
Seite 392
Abb. 5.30 Abschätzung des Fehlers bei der Mittenregel. In der linken Hälfte überschätzt die Mittenregel das Integral um die grüne Fläche, in der rechten unterschätzt sie ebenfalls um die grüne Fläche. Durch Punktspiegelung dieser an der Mitte kann man links die Differenzfläche (rot schraffiert) als Gesamtfehler erkennen. Sie hat jeweils Ordinatenstücke der Länge \(d_1-d_2\), (gelb). Diese roten Differenzflächen sind in Abb. 5.29 d) bei der Mittenregel abgebildet.


RechnungMittenregel.ggb
Abschnitt
5.5.3.4
Seite 395
Abb. 5.31 Die erste Gauß'sche Quadraturformel integriert Polynome bis zum Grad 3 exakt: von links nach rechts sind die Funktionen \(f(x)=1,\ x,\ x^2,\ x^3\) gezeigt. In blau sieht man jedes Mal die Ordinaten an den Stützstellen \(\pm\frac{1}{\sqrt3}\), diese Werte sind aber nur ganz links eingetragen. Außerdem ist für \(x^2\) das (blau schraffierte) Rechteck gezeigt, das die Fläche nach der Gauß'schen Quadraturformel darstellt.


GaussQuadratur.ggb
Abschnitt
5.5.4.1
Seite 396
Abb. 5.32 Raumintegral numerisch über ein Rechteck. \(int_{[a,b]\times[c,d]} f(x,y)\, {\rm d} x\,{\rm d} y\) für \(f(x,y)=(x-2)^2+\frac{1}{2}(y-3)^2-7\). Bezüglich \(x\) wird die Drei-Achtel-Regel, bzgl. \(y\) die Kepler'sche Regel verwendet. Die Gewichte für die Punkte im Rechteck ergeben sich als die Produkte der Gewichte für die einzelnen Variablen. Z.B. ist das Gewicht zu \((x_1,y_1)\) gegeben als \(\frac{3}{8}\cdot\frac{4}{6} = \frac{1}{4}\). Gleiche Gewichte haben gleiche Farben. Die Realisierung und die Rechnungen in diesem Fall, erstellt mit GeoGebra, finden Sie im folgenden Zusatz.


RaumintegralRechteck.ggb
Zusatz
Abschnitt 5.5.4.1
Seite 396
Der wahre Wert des Integrals errechnet sich zu \[ \int_1^3\int_1^4 f(x,y)\,{\rm d} x\,{\rm d} y = -32. \] Für die numerische Integration brauchen wir zunächst eine Wertetabelle: \[ \begin{array}{r|c|c|c|c|} f(x,y) & x=1 & 2 & 3 & 4 \\ \hline y=1 & -4.0 & -5.0 & -4.0 & -1.0 \\ \hline 2 & -5.5 & -6.5 & -5.5 & -2.5 \\ \hline 3 & -6.0 & -7.0 & -6.0 & -3.0 \\ \hline \end{array} \] Damit erhält man für die Näherung \[ \int_1^3\int_1^4 f(x,y)\,{\rm d} x\,{\rm d} y \approx (4-1)\cdot(3-1)\cdot \begin{pmatrix} -4.0 \cdot\frac{1}{48} -5.0 \cdot\frac{1}{16} -4.0 \cdot\frac{1}{16} -1.0\cdot\frac{1}{48} \\ -5.5 \cdot\frac{1}{12} -6.5 \cdot\frac{1}{4}\ -5.5 \cdot\frac{1}{4}\ -2.5\cdot\frac{1}{12} \\ -6.0 \cdot\frac{1}{48} -7.0 \cdot\frac{1}{16} -6.0 \cdot\frac{1}{16} -3.0\cdot\frac{1}{48} \\ \end{pmatrix} = -32 \] die natürlich exakt ist, weil der Polynomgrad kleiner als 3 ist.
Abschnitt
5.5.4.2
Seite 397
Abb. 5.33 Raumintegral numerisch über ein achsenparalleles Dreieck. a) Stützpunkt ist der Dreiecksschwerpunkt, d.h. der Schnitt der Seitenhalbierenden. b) Die Stützpunkte liegen auf den Seitenmitten. c) Der blaue Stützpunkt ist der Schwerpunkt. Die grünen Punkte sind drei Punkte eines Rechtecks mit linker unterer Ecke \(\frac{9-2\sqrt{15}}{21}(h,k)\) und rechter oberer Ecke \(\frac{6+\sqrt{15}}{21}(h,k)\), die roten Punkte sind drei Punkte eines Rechtecks mit linker unterer Ecke \(\frac{6-\sqrt{15}}{21}(h,k)\) und rechter oberer Ecke \(\frac{9+2\sqrt{15}}{21}(h,k)\); beide Male gehört eine Rechteckecke nicht zu den Stützstellen. Die zugehörigen Gewichte sind eingetragen.


RaumintegralDreieck.ggb
Abschnitt
5.5.4.2
Seite 398
Abb. 5.34 Beispiel für ein Raumintegral. Das Integral \(\int_{-1}^1\int_{-1}^1 \frac{1}{2\pi}{\rm e}^{-\frac{1}{2}(x^2+y^2)}\,{\rm d} x\,{\rm d} y\) wird auf vier Weisen numerisch berechnet. Der richtige Wert ist \(0.466065\), als Näherungen ergeben sich
a) \(0.569672\), b) \(0.485574\),
c) \(0.469626\), d) \(0.466559\).
Hier sieht man beim Vergleich von b) und c), dass mit weniger Punkten und verbesserter Fehlerordnung eine höhere Genauigkeit erzielt wird.


RaumintegralBsp.ggb

5.6 DGLn numerisch lösen
 In den Abb. 5.35 bis 5.44 sind alle als "Lösung" bezeichneten Kurven mit GeoGebra numerisch erzeugt, und zwar mit dem Runge-Kutta-Verfahren.
Abschnitt 5.6.5 zeigt den Aufwand für die exakten Lösungen in Abb. 5.45.
In den folgenden Dateien ist das Richtungsfeld ein Bild. Das in GeoGebra verfügbare Richtungsfeld besteht nur aus Strichlein. Aber bei Abb. 4.4 haben wie ein für alle DGLn 1. Ordnung nutzbares schönes Richtungsfeld in GeoGebra verwirklicht.
Abschnitt
5.6
Seite 399
Abb. 5.35 Untersuchungen zum Richtungsfeld von \(y'=y^2-x=f(x,y)\).
a) Richtungsfeld mit zwei auffälligen Kurven. b) Richtungsfeld mit Isoklinen für \(m=\) \(2,\,0,\,-2,\,-4,\,-6\). Berechnung der Isoklinen bei Abb. 4.6
c) Mit dem Spurmodus von GeoGebra eingezeichnete Scharen von Lösungskurven, die Sie hier selbst erkunden können. In der ggb-Datei sind einige davon als Folge von Lösungskruven fest realisiert.
Hervorzuheben ist die schwarze Lösung, die die blauen und grünen Lösungen voneinander trennt. Diese ist in Abschnit 5.65 exakt berechnet und in Abb. 5.45 violett dargestellt.
a) y2minusx-Richtungsfeld+eine+Loesung.ggb
b) y2minusx-isoklinen.ggb
c) y2minusx-isoklinen+Loesung.ggb
Abschnitt
5.6
Seite 400
Abb. 5.36 Richtungsfeld von \(y'=y^2-x=f(x,y)\). a) Nullisokline, sowie farblich hervorgehobene Richtungsvektoren aus anderen Isoklinen, b) die grünen und die blauen Lösungskurven, werden exakt von der Grenzlinie getrennt, die in Abb. 5.46 violett eingetragen ist. In deren Nähe werden die Lösungen (von links nach rechts gesehen) gestreut.
Die blauen Lösungskurven haben alle einen Scheitel auf der Nullisokline, von dort an ist ihre Farbe hellblau und sie haben diese Nullisokline auch als als Asymptote. Deren Nähe werden die die hellblauen Lösungen gesammelt.
a) siehe oben
b) y2minusx-isoklinen+Loesung.ggb
Wichtiger Hinweis für die Abbildungen 5.37 bis 5.41:
Beim Ändern der Schrittweite h werden oft falsche Graphen angezeigt.
Mit Strg r muss man eine Neuberechnung erzwingen.
Abschnitt
5.6.1
Seite 401
Abb. 5.37 Eulerverfahren für \(y'=y^2-x\).
Schrittweite ist \(h=0.5\), Anfangswerte sind \(y(0) = 0.7\), also ist \(A=P_0=(0,\,0.7)\). Es folgen \(P_1=(0.5,\,0.945), P_2=(1.0,\,1.142)\dots\)
Ab etwa \(x=4\) ist die blaue Lösungskurve recht dicht an der (rot gepunkteten) Nullisokline. Wenn ein Eulerschritt diese nach unten überschreitet, wird der nächste Schritt nach oben gehen. Im Bild wechseln nun positive und negative Steigungen. Man beachte aber, dass die blaue Lösung rechts niemals wieder Steigung null haben wird, sie liegt stets über der Nullisokline. Das Eulervefahren "entgleist" hier.
y2minusx-Euler.ggb
Kleiner Trick: Es wird mit dem Schieberegler eigentlich \(n\) gesteuert und daraus dann \(h=2^{-n}\) bestimmt.
Abschnitt
5.6.1
Seite 402
Abb. 5.38 Eulerverfahren für \(y'=y^2-x\). Schrittweiten sind \(h=0.5,\ 0.25,\ 0.125\), Anfangswerte wieder \(y(0) = 0.7\). Die schwarzen Balken dokumentieren die Halbierung des globalen Fehlers an der Stelle 2 bei Halbierung der Schrittweite.

In der Datei (=vorige) ist der Fehlerbalken einzuschalten.
y2minusx-Euler.ggb
Abschnitt
5.6.2
Seite 405
Abb. 5.39 Heunverfahren für \(y'=y^2-x\) schrittweise und im ganzen Verlauf.
a) ist ein Ausschnitt aus c), und b) zeigt den in a) sichtbaren 4. Verfahrensschritt vergrößert im zweiten Grafikfenster von GeoGebra. Daraus sind \(P_3\), \(P_4\) und die Richtung nach \(P_4\) (orange) aus der Konstruktion in b) übertragen worden und demonstrieren die Korrektheit des geometrischen Verfahrens.
Speziell b): Mit der in \(P_3\) gültigen (lila) Pfeilrichtung erreicht man im Abzissenabstand \(h=0.5\) den Zwischenpunkt \(Z\), in dem die durch den hellgrünen Vektor bzw. die grüne Gerade gezeigte Richtung gilt. Der Mittelwert der beiden Geradensteigungen ist hier geometrisch durch die durch \(M\) gemittelte Strecke und die blau-gestrichelte Gerade \(MZ\) visualisiert. Auf einer Parallelen zu ihr durch \(P_3\) erreicht man an der Abszisse von \(Z\) Punkt \(P_4\).
y2minusx-Heun.ggb
  In c) ist ein Bild einer "multifunktionalen" GeoGebra-Datei gezeigt. Hier können Sie nicht nur die Schrittweite \(h\) und den Start \(A\) steuern, sondern mit \(\alpha\) mehrere der im Folgenden dargestellten Verfahren auswählen. Hier gilt \(\alpha=1\) und \(A:\, y(0) = 0.7\). Vergleichen Sie mit dem Eulerverfahren Abb. 5.37, das man mit \(\alpha=0\) erhält, und mit dem modifizierten Eulerverfahren in Abb. 5.41, das \(\alpha=0.5\) erfordert. y2minusx-Heun-mit-alpha.ggb
Abschnitt
5.6.2
Seite 405
Abb. 5.40 Heunverfahren für \(y'=y^2-x\). Schrittweiten sind \(h=1,\ 0.5,\ 0.25\). Drei Fehlerbalken sind hervorgehoben und zeigen, dass sich der globale Fehler viertelt, wenn die Schrittweite halbiert wird. Vergleichen Sie mit Abb. 5.38.

In der Datei (=vorige) ist die Fehlergerade einzuschalten.
y2minusx-Heun-mit-alpha.ggb
Abschnitt
5.6.2.2
Seite 408

y2minusx-Euler-modi.ggb
y2minusx-Euler-modi-mit-alpha.ggb
Abb. 5.41 Modifiziertes Eulerverfahren \(y'=y^2-x\) schrittweise und im gesamten Verlauf.
a) und b) Der Schritt von \(P_3=(1.5,0.964)\) nach \(P_4\) wird genau betrachtet. Zum Zwischenpunkt \(Z\) (violett) gelangt man mit einem Eulerschritt mit \(m_3=f(x_3,y_3)\). Der Hilfspunkt \(H\) auf der Hälfte der Strecke \(\overline{P_3 Z}\) ist \(H=(x_3+\frac{h}{2}, y_3+m_3\cdot\frac{h}{2})\). Dort hat das Richtungfeld die Steigung \(m_H=f(x_H,y_H)\) (hellblau). Mit dieser Steigung vollzieht man einen ganzen Eulerschritt und gelangt zu \(P_4=(x_3+h,y_3+m_H\cdot h)\).
b) zeigt dieses Vorgehen in größerem Maßstab im zweiten Grafikfenster von GeoGebra. Die gelbe Gerade und \(P_4\) sind in das Hauptfenster a) übertragen worden. Man erkennt die Übereinstimmung.
c) Gesamtverlauf mit Startpunkt \(A\colon\ y(0)=0.7\). Vergleichen Sie mit Abb. 5.37, Abb. 5.39 c) und Abb. 5.43 b).
Abschnitt
5.6.3
Seite 410
Abb. 5.42 Einsatz von numerischer Integration für DGLn. a) Riemann'sche Summe \(\to\) Eulerverfahren, b) Mittenregel \(\to\) modifiziertes Eulerverfahren, c) Trapezregel \(\to\) Heunverfahren, d) Kepler'sche Regel \(\to\) ? Lesen Sie den Text im Buch!
Abschnitt
5.6.4
Seite 411
Abb. 5.43 Runge-Kutta-Verfahren.
a) Skizze eines Schrittes von \(P_3\) aus:
Die Hilfspunkte werden in folgender Reihenfolge aufgesucht:
Mit \(m_3\) (lila) einen halben Eulerschritt zu \(H\), dort findet man \(m_H\) (hellblau).
Mit \(m_H\) einen halben Eulerschritt zu \(R\), dort findet man \(m_r\) (grün).
Mit \(m_r\) einen ganzen Eulerschritt zu \(S\), dort findet man \(m_s\) (violett).
Der Zielpunkt \(P_4{=}(x_4,y_4)\) wird mit einem ganzen Eulerschritt erreicht (rot), bei dem man \(m = \frac{1}{6}(m_3+2m_H+2m_R+m_S)\) verwendet.
b) Die Punkte des Runge-Kutta-Verfahrens liegen praktisch perfekt auf der Lösungskurve.
y2minusx-RungeKutta.ggb
Abschnitt
5.6.4
Seite 411
Abb. 5.44 Numerische DGL-Verfahren für \(y'=y^2-x\) im Vergleich mit \(h=0.5\). Der exakte Wert (blau) an der Stelle \(x=3\) ist \(y(3)=-1.47393\). Das Runge-Kutta-Verfahren (rot) bringt: \(y(3)\approx -1.50044\), das ist ein Fehler von nur \(1.4\%\). Modifiziertes Eulerverfahren (ocker): \(y(3)\approx -1.3515\), Heun (lila): \(y(3)\approx -1.1338\), Euler (grün): \(y(3)\approx 0.9641\). \(A\) ist beweglich.
y2minusx-alle-Vergleich-Pfeile.ggb
y2minusx-alle-Vergleich.ggb
Abschnitt
5.6.5
Seite 412
Abb. 5.45 \(y'=y^2-x\) exakt gelöst. eine Vielzahl von numerisch bestimmen Lösungen ist in Abb.5.35 und 5.36 dargestellt.


Lösung-y'=y2-x, Mathematica Quelltext

Dies zum Lesen und Verstehen
Abschnitt
5.6.5
Seite 413
Abb. 5.46 \(y'=y^2-x\) exakt gelöst, aber numerisch ungenau. Wenn man statt der Airy-Funktionen die ursprünglichen Besselfunktionen verwendet, gerät selbst Mathematica hier offenbar in numerische Schwierigkeiten, wie man am rechten Bildrand deutlich sieht. Auch ist diese "Lösungskurve" bei \(x=0\) plötzlich unstetig. Tatsächlich wurden hier die Realteile der scheinbar komplexen Lösungen gezeichnet, was aber am Ergebnis nichts ändern dürfte, da die Lösungen ja eigentlich reell sind.


Lösung-y'=y2-x, Mathematica Quelltext
Dies zum Lesen und Verstehen


Abschnitt
5.6.6.2
Seite 415

Bei allen ist nur des Bild anzusehen.
Abb. 5.47 Mehrschrittverfahren veranschaulicht. a) explizites Mehrschrittverfahren, der rote Punkt ergibt sich aus Formel 5.1, b) implizites Mehrschrittverfahren, \(p_{n+1}\) wird mit unbekanntem \(y_1\) aufgestellt, c) leapfrog-Verfahren. Bei letzterem ist der Integrationsbereich ein anderer, nämlich von \(x_{-1}\) bis \(x_1\). Das Integral wird durch die Mittenregel approximiert, weshalb nur \(f(x_0,y_0)\) in die Rechnung eingeht.
DGLnumMehrschritt1.ggb
DGLnumMehrschritt2.ggb
DGLnumleapfrog3.ggb
Abschnitt
5.6.6.4
Seite 417
Abb. 5.48 Das Eulerverfahren ist instabil. Die DGL \(y' = -\lambda y\) zu den Anfangswerten \(y(0) =1\) hat die exakte Lösung \(y(x) = {\rm e}^{-\lambda x}\). Sie wird mit dem Eulerverfahren approximiert. Im Bild ist \(\lambda = 5\). In a), b), c) sind die Schrittweiten \(h=0.45,\ 0.4,\ 0.35\). Bei \(h>0.4\) wird das Verfahren instabil!
DGLnumEulerinstabil.ggb
Website zum Buch: Höhere Mathematik sehen und verstehen 05 Numerik   START und TIPPS     Liste der Zusätze Erstellt: Dez.2020, Update 06. Dezember 2021
 Prof. Dr. Dörte Haftendorn  Prof. Dr. Dieter Riebesehl  Dipl.-Math. Dipl.-Ing. Hubert Dammer