Einschrittverfahren

Einschrittverfahren s​ind in d​er numerischen Mathematik n​eben den Mehrschrittverfahren e​ine große Gruppe v​on Rechenverfahren z​ur Lösung v​on Anfangswertproblemen. Diese Aufgabenstellung, b​ei der e​ine gewöhnliche Differentialgleichung zusammen m​it einer Startbedingung gegeben ist, spielt i​n allen Natur- u​nd Ingenieurwissenschaften e​ine zentrale Rolle u​nd gewinnt beispielsweise a​uch in d​en Wirtschafts- u​nd Sozialwissenschaften i​mmer mehr a​n Bedeutung. Anfangswertprobleme werden verwendet, u​m dynamische Vorgänge z​u analysieren, z​u simulieren o​der vorherzusagen.

Einschrittverfahren nähern die Lösung (blau) eines Anfangswert­problems an, indem vom gegebenen Startpunkt ${\displaystyle A_{0}}$ aus nacheinander Punkte ${\displaystyle A_{1}}$, ${\displaystyle A_{2}}$ usw. bestimmt werden

Die namensgebende Grundidee d​er Einschrittverfahren ist, d​ass sie ausgehend v​on dem gegebenen Anfangspunkt Schritt für Schritt entlang d​er gesuchten Lösung Näherungspunkte berechnen. Dabei verwenden s​ie jeweils n​ur die zuletzt bestimmte Näherung für d​en nächsten Schritt, i​m Gegensatz z​u den Mehrschrittverfahren, d​ie auch weiter zurückliegende Punkte i​n die Rechnung miteinbeziehen. Die Einschrittverfahren lassen s​ich grob i​n zwei Gruppen einteilen: i​n die expliziten Verfahren, d​ie die n​eue Näherung direkt a​us der a​lten berechnen, u​nd in d​ie impliziten Verfahren, b​ei denen d​azu eine Gleichung gelöst werden muss. Letztere eignen s​ich auch für sogenannte steife Anfangswertprobleme.

Das einfachste u​nd älteste Einschrittverfahren, d​as explizite Euler-Verfahren, w​urde 1768 v​on Leonhard Euler veröffentlicht. Nachdem 1883 e​ine Gruppe v​on Mehrschrittverfahren vorgestellt worden war, entwickelten u​m 1900 Carl Runge, Karl Heun u​nd Wilhelm Kutta deutliche Verbesserungen d​es eulerschen Verfahrens. Aus diesen g​ing die große Gruppe d​er Runge-Kutta-Verfahren hervor, d​ie die wichtigste Klasse v​on Einschrittverfahren bildet. Weitere Entwicklungen d​es 20. Jahrhunderts s​ind beispielsweise d​ie Idee d​er Extrapolation, v​or allem a​ber Überlegungen z​ur Schrittweitensteuerung, a​lso zur Wahl geeigneter Längen d​er einzelnen Schritte e​ines Verfahrens. Diese Konzepte bilden d​ie Grundlage, u​m schwierige Anfangswertprobleme, w​ie sie i​n modernen Anwendungen auftreten, effizient u​nd mit d​er benötigten Genauigkeit d​urch Computerprogramme lösen z​u können.

Einführung

Gewöhnliche Differentialgleichungen

Die Entwicklung d​er Differential- u​nd Integralrechnung d​urch den englischen Physiker u​nd Mathematiker Isaac Newton u​nd unabhängig d​avon durch d​en deutschen Universalgelehrten Gottfried Wilhelm Leibniz i​m letzten Drittel d​es 17. Jahrhunderts w​ar ein wesentlicher Impuls für d​ie Mathematisierung d​er Wissenschaft i​n der frühen Neuzeit. Diese Methoden bildeten d​en Startpunkt d​es mathematischen Teilgebiets d​er Analysis u​nd sind i​n allen Natur- u​nd Ingenieurwissenschaften v​on zentraler Bedeutung. Während Leibniz v​on dem geometrischen Problem, Tangenten a​n gegebene Kurven z​u bestimmen, z​ur Differentialrechnung geführt wurde, g​ing Newton v​on der Fragestellung aus, w​ie sich Änderungen e​iner physikalischen Größe z​u einem bestimmten Zeitpunkt bestimmen lassen.[1]

Zum Beispiel ergibt sich bei der Bewegung eines Körpers dessen Durchschnittsgeschwindigkeit einfach als die zurückgelegte Strecke dividiert durch die dafür benötigte Zeit. Um jedoch die Momentangeschwindigkeit ${\displaystyle v(t)}$ des Körpers zu einem bestimmten Zeitpunkt ${\displaystyle t}$ mathematisch zu formulieren, ist ein Grenzübergang notwendig: Man betrachtet kurze Zeitspannen der Länge ${\displaystyle \Delta t}$, die dabei zurückgelegten Wegstrecken ${\displaystyle \Delta x}$ und die zugehörigen Durchschnittsgeschwindigkeiten ${\displaystyle {\tfrac {\Delta x}{\Delta t}}}$. Wenn man nun die Zeitspanne ${\displaystyle \Delta t}$ gegen null konvergieren lässt und wenn sich dabei die Durchschnittsgeschwindigkeiten ebenfalls einem festen Wert annähern, dann wird dieser Wert die (Momentan-)Geschwindigkeit ${\displaystyle v(t)}$ zu dem gegebenen Zeitpunkt ${\displaystyle t}$ genannt. Bezeichnet ${\displaystyle x(t)}$ die Position des Körpers zur Zeit ${\displaystyle t}$, dann schreibt man ${\displaystyle v(t)=x'(t)}$ und nennt ${\displaystyle v}$ die Ableitung von ${\displaystyle x}$.

Der entscheidende Schritt in die Richtung der Differentialgleichungsmodelle ist nun die umgekehrte Fragestellung: Im Beispiel des bewegten Körpers sei also zu jedem Zeitpunkt ${\displaystyle t}$ die Geschwindigkeit ${\displaystyle v(t)}$ bekannt und daraus soll seine Position ${\displaystyle x(t)}$ bestimmt werden. Es ist anschaulich klar, dass zusätzlich die Anfangsposition des Körpers zu einem Zeitpunkt ${\displaystyle t_{0}}$ bekannt sein muss, um dieses Problem eindeutig lösen zu können. Es ist also eine Funktion ${\displaystyle x(t)}$ mit ${\displaystyle x'(t)=v(t)}$ gesucht, die die Anfangsbedingung ${\displaystyle x(t_{0})=x_{0}}$ mit gegebenen Werten ${\displaystyle t_{0}}$ und ${\displaystyle x_{0}}$ erfüllt.

Im Beispiel der Bestimmung der Position ${\displaystyle x}$ eines Körpers aus seiner Geschwindigkeit ist die Ableitung der gesuchten Funktion explizit gegeben. Meist liegt jedoch der wichtige allgemeine Fall gewöhnlicher Differentialgleichungen für eine gesuchte Größe ${\displaystyle y}$ vor: Aufgrund der Naturgesetze oder der Modellannahmen ist ein Funktionszusammenhang bekannt, der angibt, wie die Ableitung ${\displaystyle y'(t)}$ der zu bestimmenden Funktion aus ${\displaystyle t}$ und aus dem (unbekannten) Wert ${\displaystyle y(t)}$ berechnet werden kann. Zusätzlich muss wieder eine Anfangsbedingung gegeben sein, die beispielsweise aus einer Messung der gesuchten Größe zu einem fest gewählten Zeitpunkt erhalten werden kann. Zusammengefasst liegt also der folgende allgemeine Aufgabentyp vor: Man finde die Funktion ${\displaystyle y}$, die die Gleichungen

${\displaystyle y'(t)=f(t,y(t)),\quad y(t_{0})=y_{0}}$

erfüllt, wobei ${\displaystyle f}$ eine gegebene Funktion ist.

Die dargestellte Lösung der Differentialgleichung des Lorenz-Attraktors ist eine sehr komplizierte Kurve im dreidimensionalen Raum

Ein einfaches Beispiel ist eine Größe ${\displaystyle y}$, die exponentiell wächst. Das bedeutet, dass die momentane Änderung, also die Ableitung ${\displaystyle y'(t)}$ proportional zu ${\displaystyle y(t)}$ selbst ist. Es gilt also ${\displaystyle y'(t)=\lambda y(t)}$ mit einer Wachstumsrate ${\displaystyle \lambda }$ und beispielsweise einer Anfangsbedingung ${\displaystyle y(0)=y_{0}}$. Die gesuchte Lösung ${\displaystyle y}$ lässt sich in diesem Fall bereits mit elementarer Differentialrechnung finden und mithilfe der Exponentialfunktion angeben: Es gilt ${\displaystyle y(t)=y_{0}e^{\lambda t}}$.

Die gesuchte Funktion ${\displaystyle y}$ in einer Differentialgleichung kann vektorwertig sein, das heißt, für jedes ${\displaystyle t}$ kann ${\displaystyle y(t)=(y_{1}(t),\dotsc ,y_{d}(t))}$ ein Vektor mit ${\displaystyle d}$ Komponenten sein. Man spricht dann auch von einem ${\displaystyle d}$-dimensionalen Differentialgleichungssystem. Im Anschauungsfall eines bewegten Körpers ist dann ${\displaystyle y(t)}$ seine Position im ${\displaystyle d}$-dimensionalen euklidischen Raum und ${\displaystyle y'(t)}$ seine Geschwindigkeit zur Zeit ${\displaystyle t}$. Durch die Differentialgleichung ist also in jedem Zeit- und Raumpunkt die Geschwindigkeit der gesuchten Bahnkurve mit Richtung und Betrag vorgegeben. Daraus soll die Bahn selbst berechnet werden.

Grundidee der Einschrittverfahren

Bei der oben als Beispiel betrachteten einfachen Differentialgleichung des exponentiellen Wachstums ließ sich die Lösungsfunktion direkt angeben. Das ist bei komplizierteren Problemen im Allgemeinen nicht mehr möglich. Man kann dann zwar unter gewissen Zusatzvoraussetzungen an die Funktion ${\displaystyle f}$ zeigen, dass eine eindeutig bestimmte Lösung des Anfangswertproblems existiert; diese kann aber dann nicht mehr durch Lösungsverfahren der Analysis (wie beispielsweise Trennung der Variablen, einen Exponentialansatz oder Variation der Konstanten) explizit berechnet werden. In diesem Fall können numerische Verfahren verwendet werden, um Näherungen für die gesuchte Lösung zu bestimmen.

Die Verfahren zur numerischen Lösung von Anfangswertproblemen gewöhnlicher Differentialgleichungen lassen sich grob in zwei große Gruppen einteilen: die Einschritt- und die Mehrschrittverfahren. Beiden Gruppen ist gemeinsam, dass sie schrittweise Näherungen ${\displaystyle y_{0},y_{1},y_{2},\dotsc }$ für die gesuchten Funktionswerte ${\displaystyle y(t_{0}),y(t_{1}),y(t_{2}),\dotsc }$ an Stellen ${\displaystyle t_{0}<t_{1}<t_{2}<\ldots }$ berechnen. Die definierende Eigenschaft der Einschrittverfahren ist dabei, dass zur Bestimmung der folgenden Näherung ${\displaystyle y_{j+1}}$ nur die „aktuelle“ Näherung ${\displaystyle y_{j}}$ verwendet wird. Bei Mehrschrittverfahren gehen im Gegensatz dazu zusätzlich bereits zuvor berechnete Näherungen mit ein; ein Dreischrittverfahren würde also beispielsweise außer ${\displaystyle y_{j}}$ auch noch ${\displaystyle y_{j-1}}$ und ${\displaystyle y_{j-2}}$ zur Bestimmung der neuen Näherung ${\displaystyle y_{j+1}}$ verwenden.

Zwei Schritte des expliziten Euler-Verfahrens

Das einfachste und grundlegendste Einschrittverfahren ist das explizite Euler-Verfahren, das der Schweizer Mathematiker und Physiker Leonhard Euler 1768 in seinem Lehrbuch Institutiones Calculi Integralis vorstellte.[2] Die Idee dieser Methode ist es, die gesuchte Lösung durch eine stückweise lineare Funktion anzunähern, bei der in jedem Schritt von der Stelle ${\displaystyle t_{j}}$ zur Stelle ${\displaystyle t_{j+1}}$ die Steigung des Geradenstücks durch ${\displaystyle f(t_{j},y_{j})}$ gegeben ist. Genauer betrachtet: Durch die Problemstellung ist bereits ein Wert der gesuchten Funktion gegeben, nämlich ${\displaystyle y(t_{0})=y_{0}}$. Aber auch die Ableitung an dieser Stelle ist bekannt, denn es gilt ${\displaystyle y'(t_{0})=f(t_{0},y_{0})}$. Damit kann die Tangente an den Graphen der Lösungsfunktion bestimmt und als Näherung verwendet werden. An der Stelle ${\displaystyle t_{1}>t_{0}}$ ergibt sich mit der Schrittweite ${\displaystyle h_{0}:=t_{1}-t_{0}}$

${\displaystyle y(t_{1})\approx y_{0}+h_{0}f(t_{0},y_{0})=:y_{1}}$.

Dieses Vorgehen k​ann nun i​n den folgenden Schritten fortgesetzt werden. Insgesamt ergibt s​ich damit für d​as explizite Euler-Verfahren d​ie Rechenvorschrift

${\displaystyle y_{j+1}=y_{j}+h_{j}f(t_{j},y_{j}),\quad j=0,1,2,\dotsc }$

mit den Schrittweiten ${\displaystyle h_{j}=t_{j+1}-t_{j}}$.[3]

Das explizite Euler-Verfahren ist der Ausgangspunkt für zahlreiche Verallgemeinerungen, bei denen die Steigung ${\displaystyle f(t_{j},y_{j})}$, durch Steigungen ersetzt wird, die das Verhalten der Lösung zwischen den Stellen ${\displaystyle t_{j}}$ und ${\displaystyle t_{j+1}}$ genauer annähern. Eine zusätzliche Idee für Einschrittverfahren bringt das implizite Eulerverfahren, das ${\displaystyle f(t_{j+1},y_{j+1})}$ als Steigung verwendet. Diese Wahl erscheint auf den ersten Blick wenig geeignet, da ja ${\displaystyle y_{j+1}}$ unbekannt ist. Als Verfahrensschritt erhält man aber nun die Gleichung

${\displaystyle y_{j+1}=y_{j}+h_{j}f(t_{j+1},y_{j+1})}$

aus der ${\displaystyle y_{j+1}}$ (gegebenenfalls durch ein numerisches Verfahren) berechnet werden kann. Wählt man beispielsweise als Steigung das arithmetische Mittel aus den Steigungen des expliziten und des impliziten Euler-Verfahrens, so erhält man das implizite Trapez-Verfahren. Aus diesem lässt sich wiederum ein explizites Verfahren gewinnen, wenn man zum Beispiel das unbekannte ${\displaystyle y_{j+1}}$ auf der rechten Seite der Gleichung durch die Näherung des expliziten Euler-Verfahrens nähert, das sogenannte Heun-Verfahren.[4] All diesen Verfahren und allen weiteren Verallgemeinerungen ist die Grundidee der Einschrittverfahren gemeinsam: der Schritt

${\displaystyle y_{j+1}=y_{j}+h_{j}\Phi }$

mit einer Steigung ${\displaystyle \Phi }$, die von ${\displaystyle t_{j}}$, ${\displaystyle y_{j}}$ und ${\displaystyle h_{j}}$ sowie (bei impliziten Verfahren) von ${\displaystyle y_{j+1}}$ abhängen kann.

Definition

Mit den Überlegungen aus dem Einführungsabschnitt dieses Artikels kann der Begriff des Einschrittverfahrens wie folgt definiert werden: Gesucht sei die Lösung ${\displaystyle y}$ des Anfangswertproblems

${\displaystyle y'(t)=f(t,y(t))}$, ${\displaystyle \quad y(t_{0})=y_{0}}$.

Dabei w​erde vorausgesetzt, d​ass die Lösung

${\displaystyle y\colon I\to \mathbb {R} ^{d}}$

auf einem gegebenen Intervall ${\displaystyle I=[t_{0},T]}$ existiert und eindeutig bestimmt ist. Sind

${\displaystyle t_{0}<t_{1}<t_{2}<\ldots <t_{n}=T}$

Zwischenstellen im Intervall ${\displaystyle I}$ und ${\displaystyle h_{j}=t_{j+1}-t_{j}}$ die zugehörigen Schrittweiten, dann heißt das durch

${\displaystyle y_{j+1}=y_{j}+h_{j}\Phi (t_{j},y_{j},y_{j+1},h_{j})}$, ${\displaystyle \quad j=0,\dotsc ,n-1}$

gegebene Verfahren Einschrittverfahren mit Verfahrensfunktion ${\displaystyle \Phi }$. Wenn ${\displaystyle \Phi }$ nicht von ${\displaystyle y_{j+1}}$ abhängt, dann nennt man es explizites Einschrittverfahren. Anderenfalls muss in jedem Schritt ${\displaystyle j}$ eine Gleichung für ${\displaystyle y_{j+1}}$ gelöst werden, und das Verfahren wird implizit genannt.[5]

Konsistenz und Konvergenz

Konvergenzordnung

Globaler Fehler bei verschiedenen Schrittweiten ${\displaystyle h}$ für drei Einschritt­verfahren. In doppelt-logarithmischer Auftragung erscheint der Zusammen­hang jeweils ungefähr linear, die Steigungen entsprechen dabei den Konvergenz­ordnungen 1, 2 und 4.

Für ein praxistaugliches Einschrittverfahren sollen die berechneten ${\displaystyle y_{j}}$ gute Näherungen für die Werte ${\displaystyle y(t_{j})}$ der exakten Lösung ${\displaystyle y}$ an der Stelle ${\displaystyle t_{j}}$ sein. Da im Allgemeinen die Größen ${\displaystyle d}$-dimensionale Vektoren sind, misst man die Güte dieser Näherung mit einer Vektornorm als ${\displaystyle \|y_{j}-y(t_{j})\|}$, dem Fehler an der Stelle ${\displaystyle t_{j}}$. Es ist wünschenswert, dass diese Fehler für alle ${\displaystyle j}$ schnell gegen null konvergieren, falls man die Schrittweiten gegen null konvergieren lässt. Um auch den Fall nicht konstanter Schrittweiten zu erfassen, definiert man dazu genauer ${\displaystyle h}$ als das Maximum der verwendeten Schrittweiten und betrachtet das Verhalten des maximalen Fehlers an allen Stellen ${\displaystyle t_{j}}$ im Vergleich zu Potenzen von ${\displaystyle h}$. Man sagt, das Einschrittverfahren zur Lösung des gegebenen Anfangswertproblems habe die Konvergenzordnung ${\displaystyle p\geq 1}$, wenn die Abschätzung

${\displaystyle \max _{j=0,\dotsc ,n}\|y_{j}-y(t_{j})\|\leq Ch^{p}}$

für alle hinreichend kleinen ${\displaystyle h}$ mit einer von ${\displaystyle h}$ unabhängigen Konstante ${\displaystyle C>0}$ gilt.[6]

Die Konvergenzordnung ist die wichtigste Kenngröße für den Vergleich verschiedener Einschrittverfahren.[7] Ein Verfahren mit höherer Konvergenzordnung ${\displaystyle p}$ liefert im Allgemeinen bei vorgegebener Schrittweite einen kleineren Gesamtfehler bzw. umgekehrt sind weniger Schritte nötig, um eine vorgegebene Genauigkeit zu erreichen. Bei einem Verfahren mit ${\displaystyle p=1}$ ist zu erwarten, dass sich bei einer Halbierung der Schrittweite auch der Fehler nur ungefähr halbiert. Bei einem Verfahren der Konvergenzordnung ${\displaystyle p=4}$ kann man hingegen davon ausgehen, dass sich dabei der Fehler ungefähr um den Faktor ${\displaystyle {\bigl (}{\tfrac {1}{2}}{\bigr )}^{4}={\tfrac {1}{16}}}$ verringert.

Globaler und lokaler Fehler

Die in der Definition der Konvergenzordnung betrachteten Fehler ${\displaystyle \|y_{j}-y(t_{j})\|}$ setzen sich auf eine zunächst kompliziert erscheinende Weise aus zwei Einzelkomponenten zusammen: Natürlich hängen sie zum einen von dem Fehler ab, den das Verfahren in einem einzelnen Schritt macht, indem es die unbekannte Steigung der gesuchten Funktion durch die Verfahrensfunktion annähert. Zum anderen ist aber zusätzlich zu berücksichtigen, dass bereits der Startpunkt ${\displaystyle (t_{j},y_{j})}$ eines Schrittes im Allgemeinen nicht mit dem exakten Startpunkt ${\displaystyle (t_{j},y(t_{j}))}$ übereinstimmt; der Fehler nach diesem Schritt hängt also auch von allen Fehlern ab, die bereits in den vorangegangenen Schritten gemacht wurden. Aufgrund der einheitlichen Definition der Einschrittverfahren, die sich nur in der Wahl der Verfahrensfunktion ${\displaystyle \Phi }$ unterscheiden, lässt sich aber beweisen, dass man (unter gewissen technischen Voraussetzungen an ${\displaystyle \Phi }$) direkt von der Fehlerordnung in einem einzelnen Schritt, der sogenannten Konsistenzordnung, auf die Konvergenzordnung schließen kann.

Der Begriff d​er Konsistenz i​st ein allgemeines u​nd zentrales Konzept d​er modernen numerischen Mathematik. Während m​an bei d​er Konvergenz e​ines Verfahrens untersucht, w​ie gut d​ie numerischen Näherungen z​ur exakten Lösung passen, stellt m​an sich b​ei der Konsistenz vereinfacht gesprochen d​ie „umgekehrte“ Frage: Wie g​ut erfüllt d​ie exakte Lösung d​ie Verfahrensvorschrift? In dieser allgemeinen Theorie gilt, d​ass ein Verfahren g​enau dann konvergent ist, w​enn es konsistent u​nd stabil ist. Um d​ie Notation z​u vereinfachen, s​oll in d​er folgenden Überlegung angenommen werden, d​ass ein explizites Einschrittverfahren

${\displaystyle y_{j+1}=y_{j}+h\Phi (t_{j},y_{j},h)}$

mit konstanter Schrittweite ${\displaystyle h}$ vorliegt. Mit der wahren Lösung ${\displaystyle t\mapsto y(t)}$ definiert man den lokalen Abschneidefehler (auch lokaler Verfahrensfehler genannt) ${\displaystyle \eta }$ als[8]

${\displaystyle \eta (t,h)=y(t)+h\Phi (t,y(t),h)-y(t+h)}$.

Man nimmt also an, dass die exakte Lösung bekannt ist, startet einen Verfahrensschritt an dem Punkt ${\displaystyle (t,y(t))}$ und bildet die Differenz zur exakten Lösung an der Stelle ${\displaystyle t+h}$. Damit definiert man: Ein Einschrittverfahren hat die Konsistenzordnung ${\displaystyle p\geq 1}$, wenn die Abschätzung

${\displaystyle \|\eta (t,h)\|\leq Ch^{p+1}}$

für alle hinreichend kleinen ${\displaystyle h}$ mit einer von ${\displaystyle h}$ unabhängigen Konstante ${\displaystyle C>0}$ gilt.

Der auffällige Unterschied zwischen den Definitionen der Konsistenzordnung und der Konvergenzordnung ist die Potenz ${\displaystyle h^{p+1}}$ anstelle von ${\displaystyle h^{p}}$. Das lässt sich anschaulich so deuten, dass beim Übergang vom lokalen zum globalen Fehler eine Potenz der Schrittweite „verloren geht“. Es gilt nämlich der folgende, für die Theorie der Einschrittverfahren zentrale Satz:[9]

Ist die Verfahrensfunktion ${\displaystyle \Phi }$ Lipschitz-stetig und hat das zugehörige Einschrittverfahren die Konsistenzordnung ${\displaystyle p}$, dann hat es auch die Konvergenzordnung ${\displaystyle p}$.

Die Lipschitz-Stetigkeit der Verfahrensfunktion als Zusatzvoraussetzung für die Stabilität ist im Allgemeinen immer dann erfüllt, wenn die Funktion ${\displaystyle f}$ aus der Differentialgleichung selbst Lipschitz-stetig ist. Diese Forderung muss für die meisten Anwendungen sowieso vorausgesetzt werden, um die eindeutige Lösbarkeit des Anfangswertproblems zu garantieren. Nach dem Satz genügt es also, die Konsistenzordnung eines Einschrittverfahrens zu bestimmen. Das lässt sich prinzipiell durch Taylor-Entwicklung von ${\displaystyle \eta (t,h)}$ nach Potenzen von ${\displaystyle h}$ erreichen. In der Praxis werden die entstehenden Formeln für höhere Ordnungen sehr kompliziert und unübersichtlich, sodass zusätzliche Konzepte und Notationen benötigt werden.[10]

Steifheit und A-Stabilität

→ Hauptartikel: Steifes Anfangswertproblem und A-Stabilität

Die Konvergenzordnung e​ines Verfahrens i​st eine asymptotische Aussage, d​ie das Verhalten d​er Näherungen beschreibt, w​enn die Schrittweite g​egen null konvergiert. Sie s​agt jedoch nichts darüber aus, o​b das Verfahren für e​ine gegebene f​este Schrittweite a​uch tatsächlich e​ine brauchbare Näherung berechnet. Dass d​ies bei bestimmten Typen v​on Anfangswertproblemen tatsächlich e​in großes Problem darstellen kann, beschrieben zuerst Charles Francis Curtiss u​nd Joseph O. Hirschfelder 1952. Sie hatten beobachtet, d​ass bei manchen Differentialgleichungssystemen d​er chemischen Reaktionskinetik d​ie Lösungen n​icht mit expliziten numerischen Verfahren berechnet werden können, u​nd nannten solche Anfangswertprobleme „steif“.[11] Es existieren zahlreiche mathematische Kriterien, u​m für e​in gegebenes Problem festzustellen, w​ie steif e​s ist. Anschaulich handelt e​s sich b​ei steifen Anfangswertproblemen m​eist um Differentialgleichungssysteme, b​ei denen einige Komponenten s​ehr schnell konstant werden, während andere Komponenten s​ich nur langsam ändern. Ein solches Verhalten t​ritt typischerweise b​ei der Modellierung chemischer Reaktionen auf. Die für d​ie praktische Anwendung nützlichste Definition v​on Steifheit i​st dabei jedoch: Ein Anfangswertproblem i​st steif, w​enn man b​ei seiner Lösung m​it expliziten Einschrittverfahren d​ie Schrittweite „zu klein“ wählen müsste, u​m eine brauchbare Lösung z​u erhalten. Solche Probleme können a​lso nur m​it impliziten Verfahren gelöst werden.[12]

Zur Berechnung einer exponentiell fallenden Lösung (blau) ist das explizite Euler-Verfahren (rot) bei zu großer Schrittweite völlig unbrauchbar; das implizite Euler-Verfahren (grün) bestimmt die Lösung für beliebige Schrittweiten qualitativ richtig.

Dieser Effekt lässt s​ich genauer darstellen, w​enn man untersucht, w​ie die einzelnen Verfahren m​it exponentiellem Zerfall zurechtkommen. Dazu betrachtet m​an nach d​em schwedischen Mathematiker Germund Dahlquist d​ie Testgleichung

${\displaystyle y'(t)=\lambda y(t),\quad y(0)=1}$

mit der für ${\displaystyle \lambda <0}$ exponentiell abfallenden Lösung ${\displaystyle y(t)=e^{\lambda t}}$. Die nebenstehende Grafik zeigt – exemplarisch für das explizite und das implizite Euler-Verfahren – das typische Verhalten dieser beiden Verfahrensgruppen bei diesem so einfach erscheinenden Anfangswertproblem: Verwendet man bei einem expliziten Verfahren eine zu große Schrittweite, dann ergeben sich stark oszillierende Werte, die sich im Laufe der Rechnung aufschaukeln und sich immer weiter von der exakten Lösung entfernen. Implizite Verfahren berechnen hingegen typischerweise die Lösung für beliebige Schrittweiten qualitativ richtig, nämlich als exponentiell fallende Folge von Näherungswerten.[13]

Etwas allgemeiner betrachtet man die obige Testgleichung auch für komplexe Werte von ${\displaystyle \lambda }$. In diesem Fall sind die Lösungen Schwingungen, deren Amplitude genau dann beschränkt bleibt, wenn ${\displaystyle \operatorname {Re} (\lambda )\leq 0}$ gilt, also der Realteil von ${\displaystyle \lambda }$ kleiner oder gleich 0 ist. Damit lässt sich eine wünschenswerte Eigenschaft von Einschrittverfahren formulieren, die für steife Anfangswertprobleme eingesetzt werden sollen: die sogenannte A-Stabilität. Ein Verfahren heißt A-stabil, wenn es für beliebige Schrittweiten ${\displaystyle h>0}$ angewendet auf die Testgleichung für alle ${\displaystyle \lambda }$ mit ${\displaystyle \operatorname {Re} (\lambda )\leq 0}$ eine Folge von Näherungen ${\displaystyle y_{0},y_{1},y_{2},\dotsc }$ berechnet, die (wie die wahre Lösung) beschränkt bleibt. Das implizite Euler-Verfahren und das implizite Trapez-Verfahren sind die einfachsten Beispiele A-stabiler Einschrittverfahren. Andererseits lässt sich zeigen, dass ein explizites Verfahren niemals A-stabil sein kann.[14][15]

Spezielle Verfahren und Verfahrensklassen

Einige Einschrittverfahren im Vergleich

Einfache Verfahren der Ordnung 1 und 2

Wie der französische Mathematiker Augustin-Louis Cauchy um 1820 bewies, besitzt das Euler-Verfahren die Konvergenzordnung 1.[16] Wenn man die Steigungen ${\displaystyle f(t_{j},y_{j})}$ des expliziten Euler-Verfahrens und ${\displaystyle f(t_{j+1},y_{j+1})}$ des impliziten Euler-Verfahrens, wie sie an den beiden Endpunkten eines Schritts vorliegen, mittelt, kann man hoffen, eine bessere Näherung über das ganze Intervall zu erhalten. Tatsächlich lässt sich beweisen, dass das so erhaltene implizite Trapez-Verfahren

${\displaystyle y_{j+1}=y_{j}+{\frac {h}{2}}{\Big (}f(t_{j},y_{j})+f(t_{j+1},y_{j+1}){\Big )}}$

die Konvergenzordnung 2 hat. Dieses Verfahren weist sehr gute Stabilitätseigenschaften auf, ist allerdings implizit, sodass in jedem Schritt eine Gleichung für ${\displaystyle y_{j+1}}$ gelöst werden muss. Nähert man diese Größe auf der rechten Seite der Gleichung durch das explizite Euler-Verfahren an, so entsteht das explizite Verfahren von Heun

${\displaystyle y_{j+1}=y_{j}+{\frac {h}{2}}{\Big (}f(t_{j},y_{j})+f{\big (}t_{j+1},y_{j}+hf(t_{j},y_{j}){\big )}{\Big )}}$,

das ebenfalls die Konvergenzordnung 2 besitzt. Ein weiteres einfaches explizites Verfahren der Ordnung 2, das verbesserte Euler-Verfahren, erhält man durch folgende Überlegung: Eine „mittlere“ Steigung im Verfahrensschritt wäre die Steigung der Lösung ${\displaystyle y}$ in der Mitte des Schritts, also an der Stelle ${\displaystyle t_{j}+{\tfrac {h}{2}}}$. Da die Lösung aber unbekannt ist, nähert man sie durch einen expliziten Euler-Schritt mit halber Schrittweite an. Es ergibt sich die Verfahrensvorschrift

${\displaystyle y_{j+1}=y_{j}+hf{\big (}t_{j}+{\tfrac {h}{2}},y_{j}+{\tfrac {h}{2}}f(t_{j},y_{j}){\big )}}$.

Diese Einschrittverfahren d​er Ordnung 2 wurden a​ls Verbesserungen d​es Euler-Verfahrens a​lle 1895 v​on dem deutschen Mathematiker Carl Runge veröffentlicht.[17]

Runge-Kutta-Verfahren

Das klassische Runge-Kutta-Verfahren vierter Ordnung mittelt in jedem Schritt vier Hilfssteigungen (rot)

→ Hauptartikel: Runge-Kutta-Verfahren

Die erwähnten Ideen für einfache Einschrittverfahren führen bei weiterer Verallgemeinerung zur wichtigen Klasse der Runge-Kutta-Verfahren. Zum Beispiel lässt sich das Verfahren von Heun übersichtlicher so präsentieren: Zuerst wird eine Hilfssteigung ${\displaystyle k_{1}=f(t_{j},y_{j})}$ berechnet, nämlich die Steigung des expliziten Euler-Verfahrens. Damit wird eine weitere Hilfssteigung bestimmt, hier ${\displaystyle k_{2}=f(t_{j}+h,y_{j}+hk_{1})}$. Die tatsächlich verwendete Verfahrenssteigung ${\displaystyle \Phi }$ ergibt sich anschließend als ein gewichtetes Mittel der Hilfssteigungen, im Verfahren von Heun also ${\displaystyle {\tfrac {1}{2}}k_{1}+{\tfrac {1}{2}}k_{2}}$. Dieses Vorgehen lässt sich auf mehr als zwei Hilfssteigungen verallgemeinern. Ein ${\displaystyle s}$-stufiges Runge-Kutta-Verfahren berechnet zunächst Hilfssteigungen ${\displaystyle k_{1},\dotsc ,k_{s}}$ durch Auswertung von ${\displaystyle f}$ an geeigneten Stellen und anschließend ${\displaystyle \Phi }$ als gewichtetes Mittel. Bei einem expliziten Runge-Kutta-Verfahren werden die Hilfssteigungen ${\displaystyle k_{1},k_{2},k_{3},\dotsc }$ der Reihe nach direkt berechnet, bei einem impliziten ergeben sie sich als Lösungen eines Gleichungssystems. Als typisches Beispiel sei das explizite klassische Runge-Kutta-Verfahren der Ordnung 4 angeführt, das mitunter einfach als das Runge-Kutta-Verfahren bezeichnet wird: Dabei werden zunächst die vier Hilfssteigungen

${\displaystyle {\begin{aligned}k_{1}&=f(t_{j},y_{j})\\k_{2}&=f(t_{j}+{\tfrac {h}{2}},y_{j}+{\tfrac {h}{2}}k_{1})\\k_{3}&=f(t_{j}+{\tfrac {h}{2}},y_{j}+{\tfrac {h}{2}}k_{2})\\k_{4}&=f(t_{j}+h,y_{j}+hk_{3})\end{aligned}}}$

berechnet u​nd dann a​ls Verfahrenssteigung d​as gewichtete Mittel

${\displaystyle {\tfrac {1}{6}}k_{1}+{\tfrac {1}{3}}k_{2}+{\tfrac {1}{3}}k_{3}+{\tfrac {1}{6}}k_{4}}$

verwendet.[18] Dieses bekannte Verfahren veröffentlichte d​er deutsche Mathematiker Wilhelm Kutta i​m Jahr 1901, nachdem e​in Jahr z​uvor Karl Heun e​in dreistufiges Einschrittverfahren d​er Ordnung 3 gefunden hatte.[19]

Die Konstruktion von expliziten Verfahren noch höherer Ordnung mit möglichst kleiner Stufenzahl ist ein mathematisch recht anspruchsvolles Problem. Wie John C. Butcher 1965 zeigen konnte, gibt es zum Beispiel für die Ordnung 5 nur minimal sechsstufige Verfahren; ein explizites Runge-Kutta-Verfahren der Ordnung 8 benötigt mindestens 11 Stufen. 1978 fand der österreichische Mathematiker Ernst Hairer ein Verfahren der Ordnung 10 mit 17 Stufen. Die Koeffizienten für ein solches Verfahren müssen 1205 Bestimmungsgleichungen erfüllen.[20] Bei impliziten Runge-Kutta-Verfahren ist die Situation einfacher und übersichtlicher: Für jede Stufenzahl ${\displaystyle s}$ existiert ein Verfahren der Ordnung ${\displaystyle p=2s}$; das ist zugleich die maximal erreichbare Ordnung.[21]

Extrapolationsverfahren

Die Idee der Extrapolation ist nicht auf die Lösung von Anfangswertproblemen mit Einschrittverfahren beschränkt, sondern lässt sich sinngemäß auf alle numerischen Verfahren anwenden, die das zu lösende Problem mit einer Schrittweite ${\displaystyle h}$ diskretisieren. Ein bekanntes Beispiel eines Extrapolationsverfahrens ist etwa die Romberg-Integration zur numerischen Berechnung von Integralen. Sei daher ganz allgemein ${\displaystyle v}$ ein Wert, der numerisch bestimmt werden soll, im Fall dieses Artikels etwa der Wert der Lösungsfunktion eines Anfangswertproblems an einer gegebenen Stelle. Ein numerisches Verfahren, beispielsweise ein Einschrittverfahren, berechnet dafür einen Näherungswert ${\displaystyle {\tilde {v}}(h)}$, der von der Wahl der Schrittweite ${\displaystyle h>0}$ abhängt. Dabei sei angenommen, dass das Verfahren konvergent ist, also, dass ${\displaystyle {\tilde {v}}(h)}$ gegen ${\displaystyle v}$ konvergiert, wenn ${\displaystyle h}$ gegen null konvergiert. Diese Konvergenz ist jedoch nur eine rein theoretische Aussage, da bei der realen Anwendung des Verfahrens zwar Näherungswerte ${\displaystyle {\tilde {v}}(h_{1}),{\tilde {v}}(h_{2}),\dotsc ,{\tilde {v}}(h_{m})}$ für endlich viele verschiedene Schrittweiten ${\displaystyle h_{1}>h_{2}>\ldots >h_{m}}$ berechnet werden können, man aber selbstverständlich nicht die Schrittweite „gegen null konvergieren“ lassen kann. Die berechneten Näherungen für verschiedene Schrittweiten lassen sich jedoch als Information über die (unbekannte) Funktion ${\displaystyle {\tilde {v}}}$ auffassen: Bei den Extrapolationsverfahren wird dabei ${\displaystyle {\tilde {v}}}$ durch ein Interpolationspolynom angenähert, also durch ein Polynom ${\displaystyle P}$ mit

${\displaystyle P(h_{k})={\tilde {v}}(h_{k})}$

für ${\displaystyle k=1,2,\dotsc ,m}$. Der Wert ${\displaystyle P(0)}$ des Polynoms an der Stelle ${\displaystyle h=0}$ wird dann als berechenbare Näherung für den nicht berechenbaren Grenzwert von ${\displaystyle {\tilde {v}}(h)}$ für ${\displaystyle h}$ gegen null verwendet.[22] Einen frühen erfolgreichen Extrapolationsalgorithmus für Anfangswertprobleme veröffentlichten Roland Bulirsch und Josef Stoer im Jahr 1966.[23]

Extrapolation auf ${\displaystyle h=0}$ bei einem Verfahren der Ordnung ${\displaystyle p=2}$

Ein konkretes Beispiel im Falle eines Einschrittverfahrens der Ordnung ${\displaystyle p}$ kann das allgemeine Vorgehen der Extrapolation verständlich machen. Bei einem solchen Verfahren lässt sich die berechnete Näherung für kleine Schrittweiten ${\displaystyle h}$ gut durch ein Polynom der Form

${\displaystyle P(h)=a+bh^{p}}$

mit zunächst unbekannten Parametern ${\displaystyle a}$ und ${\displaystyle b}$ annähern. Berechnet man nun mit dem Verfahren für eine Schrittweite ${\displaystyle h_{1}}$ und für die halbe Schrittweite ${\displaystyle h_{2}={\tfrac {1}{2}}h_{1}}$ zwei Näherungen ${\displaystyle y_{h_{1}}}$ und ${\displaystyle y_{h_{2}}}$, erhält man aus den Interpolationsbedingungen ${\displaystyle P(h_{1})=y_{h_{1}}}$ und ${\displaystyle P(h_{2})=y_{h_{2}}}$ zwei lineare Gleichungen für die Unbekannten ${\displaystyle a}$ und ${\displaystyle b}$. Der auf ${\displaystyle h=0}$ extrapolierte Wert

${\displaystyle P(0)=a=y_{h_{2}}+{\frac {y_{h_{2}}-y_{h_{1}}}{2^{p}-1}}}$

stellt dann im Allgemeinen eine deutlich bessere Näherung dar als die beiden zunächst berechneten Werte. Es lässt sich zeigen, dass die Ordnung des so erhaltenen Einschrittverfahrens mindestens ${\displaystyle p+1}$ ist, also um mindestens 1 größer als beim ursprünglichen Verfahren ist.[24]

Verfahren mit Schrittweitensteuerung

→ Hauptartikel: Schrittweitensteuerung

Ein Vorteil der Einschrittverfahren ist, dass in jedem Schritt ${\displaystyle j}$ unabhängig von den anderen Schritten eine beliebige Schrittweite ${\displaystyle h_{j}}$ verwendet werden kann. In der Praxis stellt sich dabei offensichtlich die Frage, wie ${\displaystyle h_{j}}$ gewählt werden soll. In realen Anwendungen wird stets eine Fehlertoleranz gegeben sein, mit der die Lösung eines Anfangswertproblems berechnet werden soll; zum Beispiel wäre es sinnlos, eine numerische Näherung zu bestimmen, die wesentlich „genauer“ ist als die mit Messabweichungen behafteten Daten für Startwerte und Parameter des gegebenen Problems. Das Ziel wird also sein, die Schrittweiten so zu wählen, dass einerseits die vorgegebenen Fehlertoleranzen eingehalten werden, andererseits aber auch möglichst wenige Schritte zu verwenden, um den Rechenaufwand klein zu halten. Dies wird sich im Allgemeinen nur dann erreichen lassen, wenn die Schrittweiten an den Verlauf der Lösung angepasst werden: kleine Schritte, wo sich die Lösung stark ändert, große Schritte, wo sie nahezu konstant verläuft.[25]

Bei gut konditionierten Anfangswertproblemen lässt sich zeigen, dass der globale Verfahrensfehler ungefähr gleich der Summe der lokalen Abschneidefehler ${\displaystyle \eta _{j}:=\|\eta (t_{j},h_{j})\|}$ in den einzelnen Schritten ist. Daher sollte als Schrittweite ein möglichst großes ${\displaystyle h_{j}}$ gewählt werden, für das ${\displaystyle \eta _{j}}$ unter einer gewählten Toleranzschwelle liegt. Das Problem dabei ist, dass ${\displaystyle \eta _{j}}$ nicht direkt berechnet werden kann, da es ja von der unbekannten exakten Lösung ${\displaystyle y(t_{j})}$ des Anfangswertproblems an der Stelle ${\displaystyle t_{j}}$ abhängt. Die Grundidee der Schrittweitensteuerung ist es daher, ${\displaystyle y(t_{j})}$ mit einem Verfahren anzunähern, das genauer ist als das zugrundeliegende Basisverfahren.[26]

Zwei grundlegende Ideen zur Schrittweitensteuerung sind die Schrittweitenhalbierung und die eingebetteten Verfahren. Bei der Schrittweitenhalbierung wird zusätzlich zum eigentlichen Verfahrensschritt als Vergleichswert das Ergebnis für zwei Schritte mit der halben Schrittweite berechnet. Aus beiden Werten wird dann durch Extrapolation eine genauere Näherung für ${\displaystyle y(t_{j})}$ bestimmt und damit der lokale Fehler ${\displaystyle \eta _{j}}$ geschätzt. Ist dieser zu groß, wird dieser Schritt verworfen und mit einer kleineren Schrittweite wiederholt. Ist er deutlich kleiner als die vorgegebene Toleranz, kann im nächsten Schritt die Schrittweite vergrößert werden.[27] Der zusätzliche Rechenaufwand für dieses Schrittweitenhalbierungsverfahren ist relativ groß; deshalb verwenden moderne Implementierungen meist sogenannte eingebettete Verfahren zur Schrittweitensteuerung. Die Grundidee ist dabei, in jedem Schritt zwei Näherungen für ${\displaystyle y(t_{j})}$ mit zwei Einschrittverfahren zu berechnen, die verschiedene Konvergenzordnungen haben, und damit den lokalen Fehler zu schätzen. Um den Rechenaufwand zu optimieren, sollten die beiden Verfahren möglichst viele Rechenschritte gemeinsam haben: Sie sollten „ineinander eingebettet“ sein. Eingebettete Runge-Kutta-Verfahren verwenden beispielsweise die gleichen Hilfssteigungen und unterscheiden sich nur darin, wie sie diese mitteln. Bekannte eingebettete Verfahren sind unter anderem die Runge-Kutta-Fehlberg-Verfahren (Erwin Fehlberg, 1969) und die Dormand-Prince-Verfahren (J. R. Dormand und P. J. Prince, 1980).[28]

Praxisbeispiel: Lösen von Anfangswertproblemen mit numerischer Software

Für d​ie in diesem Artikel überblicksartig dargestellten mathematischen Konzepte wurden zahlreiche Software-Implementierungen entwickelt, d​ie dem Anwender d​ie Möglichkeit geben, praktische Probleme a​uf einfache Weise numerisch z​u lösen. Als konkretes Beispiel d​azu soll n​un eine Lösung d​er Lotka-Volterra-Gleichungen m​it der verbreiteten numerischen Software Matlab berechnet werden. Die Lotka-Volterra-Gleichungen s​ind ein einfaches Modell a​us der Biologie, d​as die Wechselwirkungen v​on Räuber- u​nd Beutepopulationen beschreibt. Gegeben s​ei dazu d​as Differentialgleichungssystem

${\displaystyle {\begin{aligned}y_{1}'(t)&=ay_{1}(t)-by_{1}(t)y_{2}(t)\\y_{2}'(t)&=cy_{1}(t)y_{2}(t)-dy_{2}(t)\end{aligned}}}$

mit den Parametern ${\displaystyle a=1,b=2,c=1,d=1}$ und der Anfangsbedingung ${\displaystyle y_{1}(0)=3}$, ${\displaystyle y_{2}(0)=1}$. Hierbei entsprechen ${\displaystyle y_{1}}$ und ${\displaystyle y_{2}}$ der zeitlichen Entwicklung der Beute- bzw. der Räuberpopulation. Die Lösung soll auf dem Zeitintervall ${\displaystyle [0,20]}$ berechnet werden.

Für die Berechnung mithilfe von Matlab wird zunächst für die gegebenen Parameterwerte die Funktion ${\displaystyle f}$ auf der rechten Seite der Differentialgleichung ${\displaystyle y'=f(t,y)}$ definiert:

a = 1; b = 2; c = 1; d = 1;
f = @(t,y) [a*y(1) - b*y(1)*y(2); c*y(1)*y(2) - d*y(2)];

Außerdem werden d​as Zeitintervall u​nd die Anfangswerte benötigt:

t_int = [0, 20];
y0 = [3; 1];

Anschließend k​ann die Lösung berechnet werden:

[t, y] = ode45(f, t_int, y0);

Die Matlab-Funktion ode45 implementiert e​in Einschrittverfahren, d​as zwei eingebettete explizite Runge-Kutta-Verfahren m​it den Konvergenzordnungen 4 u​nd 5 z​ur Schrittweitensteuerung verwendet.[29]

Die Lösung kann nun gezeichnet werden, ${\displaystyle y_{1}}$ als blaue Kurve und ${\displaystyle y_{2}}$ als rote; die berechneten Punkte werden durch kleine Kreise markiert:

figure(1)
plot(t, y(:,1), 'b-o', t, y(:,2), 'r-o')

Das Ergebnis i​st unten i​m linken Bild dargestellt. Das rechte Bild z​eigt die v​om Verfahren verwendeten Schrittweiten u​nd wurde erzeugt mit

figure(2)
plot(t(1:end-1), diff(t))

Lösungen der Lotka-Volterra-Gleichungen

Verwendete Schrittweiten

Dieses Beispiel k​ann ohne Änderungen a​uch mit d​er freien numerischen Software GNU Octave ausgeführt werden. Mit d​em dort implementierten Verfahren ergibt s​ich allerdings e​ine etwas andere Schrittweitenfolge.

Literatur

• John C. Butcher: Numerical Methods for Ordinary Differential Equations. John Wiley & Sons, Chichester 2008, ISBN 978-0-470-72335-7.
• Wolfgang Dahmen, Arnold Reusken: Numerik für Ingenieure und Naturwissenschaftler. 2. Auflage. Springer, Berlin/Heidelberg 2008, ISBN 978-3-540-76492-2, Kap. 11: Gewöhnliche Differentialgleichungen.
• Peter Deuflhard, Folkmar Bornemann: Numerische Mathematik 2 – Gewöhnliche Differentialgleichungen. 3. Auflage. Walter de Gruyter, Berlin 2008, ISBN 978-3-11-020356-1.
• David F. Griffiths, Desmond J. Higham: Numerical Methods for Ordinary Differential Equations – Initial Value Problems. Springer, London 2010, ISBN 978-0-85729-147-9.
• Robert Plato: Numerische Mathematik kompakt. 4. Auflage. Vieweg+Teubner, Wiesbaden 2010, ISBN 978-3-8348-1018-2, Kap. 7: Einschrittverfahren für Anfangswertprobleme.
• Hans-Jürgen Reinhardt: Numerik gewöhnlicher Differentialgleichungen. 2. Auflage. Walter de Gruyter, Berlin/Boston 2012, ISBN 978-3-11-028045-6.
• Hans Rudolf Schwarz, Norbert Köckler: Numerische Mathematik. 8. Auflage. Vieweg+Teubner, Wiesbaden 2011, ISBN 978-3-8348-1551-4, Kap. 8: Anfangswertprobleme.
• Karl Strehmel, Rüdiger Weiner, Helmut Podhaisky: Numerik gewöhnlicher Differentialgleichungen. 2. Auflage. Springer Spektrum, Wiesbaden 2012, ISBN 978-3-8348-1847-8.

Weblinks

• Lars Grüne: Numerische Methoden für gewöhnliche Differentialgleichungen (Numerische Mathematik II). (PDF) 2008, abgerufen am 20. August 2018 (Vorlesungsskript der Universität Bayreuth).
• Peter Spellucci: Numerik gewöhnlicher Differentialgleichungen. (PDF) 2007, abgerufen am 20. August 2018 (Vorlesungsskript der TU Darmstadt).
• Hans U. Fuchs: Numerical Methods for Differential Equations. (PDF) 2007, abgerufen am 20. August 2018 (englisch, Vorlesungsskript der Zürcher Hochschule für Angewandte Wissenschaften).
• Mathe Tutorial: Rechner für allgemeine Differentialgleichungen 1. Ordnung. In: Mathe Tutorial. Abgerufen am 20. August 2018.

Einzelnachweise

1. Thomas Sonar: 3000 Jahre Analysis. Springer, Berlin/Heidelberg 2011, ISBN 978-3-642-17203-8, S. 378–388 und 401–426.
2. Jean-Luc Chabert u. a.: A History of Algorithms. Springer, Berlin/Heidelberg 1999, ISBN 978-3-540-63369-3, S. 374–378.
3. Wolfgang Dahmen, Arnold Reusken: Numerik für Ingenieure und Naturwissenschaftler. 2. Auflage. Springer, Berlin/Heidelberg 2008, ISBN 978-3-540-76492-2, S. 386 f.
4. Wolfgang Dahmen, Arnold Reusken: Numerik für Ingenieure und Naturwissenschaftler. 2. Auflage. Springer, Berlin/Heidelberg 2008, ISBN 978-3-540-76492-2, S. 386–392.
5. Hans Rudolf Schwarz, Norbert Köckler: Numerische Mathematik. 8. Auflage. Vieweg+Teubner, Wiesbaden 2011, ISBN 978-3-8348-1551-4, S. 350 f.
6. Robert Plato: Numerische Mathematik kompakt. 4. Auflage. Vieweg+Teubner, Wiesbaden 2010, ISBN 978-3-8348-1018-2, S. 157.
7. Robert Plato: Numerische Mathematik kompakt. 4. Auflage. Vieweg+Teubner, Wiesbaden 2010, ISBN 978-3-8348-1018-2, S. 156.
8. Robert Plato: Numerische Mathematik kompakt. 4. Auflage. Vieweg+Teubner, Wiesbaden 2010, ISBN 978-3-8348-1018-2, S. 157.
9. Hans-Jürgen Reinhardt: Numerik gewöhnlicher Differentialgleichungen. 2. Auflage. Walter de Gruyter, Berlin/Boston 2012, ISBN 978-3-11-028045-6, S. 42 f.
10. John C. Butcher: Numerical Methods for Ordinary Differential Equations. John Wiley & Sons, Chichester 2008, ISBN 978-0-470-72335-7, S. 95–100.
11. J. C. Butcher: Numerical methods for ordinary differential equations in the 20th century. In: Journal of Computational and Applied Mathematics. Band 125, Nr. 1–2, 15. Dezember 2000, S. 21 f. (online).
12. Peter Deuflhard, Folkmar Bornemann: Numerische Mathematik 2 – Gewöhnliche Differentialgleichungen. 3. Auflage. Walter de Gruyter, Berlin 2008, ISBN 978-3-11-020356-1, S. 228 f.
13. Peter Deuflhard, Folkmar Bornemann: Numerische Mathematik 2 – Gewöhnliche Differentialgleichungen. 3. Auflage. Walter de Gruyter, Berlin 2008, ISBN 978-3-11-020356-1, S. 229–231.
14. Wolfgang Dahmen, Arnold Reusken: Numerik für Ingenieure und Naturwissenschaftler. 2. Auflage. Springer, Berlin/Heidelberg 2008, ISBN 978-3-540-76492-2, S. 443 f.
15. Karl Strehmel, Rüdiger Weiner, Helmut Podhaisky: Numerik gewöhnlicher Differentialgleichungen. 2. Auflage. Springer Spektrum, Wiesbaden 2012, ISBN 978-3-8348-1847-8, S. 258 f.
16. Jean-Luc Chabert u. a.: A History of Algorithms. Springer, Berlin/Heidelberg 1999, ISBN 978-3-540-63369-3, S. 378 f.
17. Jean-Luc Chabert u. a.: A History of Algorithms. Springer, Berlin/Heidelberg 1999, ISBN 978-3-540-63369-3, S. 381–388.
18. Wolfgang Dahmen, Arnold Reusken: Numerik für Ingenieure und Naturwissenschaftler. 2. Auflage. Springer, Berlin/Heidelberg 2008, ISBN 978-3-540-76492-2, S. 406 f.
19. J. C. Butcher: Numerical methods for ordinary differential equations in the 20th century. In: Journal of Computational and Applied Mathematics. Band 125, Nr. 1–2, 15. Dezember 2000, S. 4–6 (online).
20. Peter Deuflhard, Folkmar Bornemann: Numerische Mathematik 2 – Gewöhnliche Differentialgleichungen. 3. Auflage. Walter de Gruyter, Berlin 2008, ISBN 978-3-11-020356-1, S. 160–162.
21. Karl Strehmel, Rüdiger Weiner, Helmut Podhaisky: Numerik gewöhnlicher Differentialgleichungen. 2. Auflage. Springer Spektrum, Wiesbaden 2012, ISBN 978-3-8348-1847-8, S. 219–221.
22. Karl Strehmel, Rüdiger Weiner, Helmut Podhaisky: Numerik gewöhnlicher Differentialgleichungen. 2. Auflage. Springer Spektrum, Wiesbaden 2012, ISBN 978-3-8348-1847-8, S. 79 ff.
23. J. C. Butcher: Numerical methods for ordinary differential equations in the 20th century. In: Journal of Computational and Applied Mathematics. Band 125, Nr. 1–2, 15. Dezember 2000, S. 26 (online).
24. Robert Plato: Numerische Mathematik kompakt. 4. Auflage. Vieweg+Teubner, Wiesbaden 2010, ISBN 978-3-8348-1018-2, S. 171–173.
25. Karl Strehmel, Rüdiger Weiner, Helmut Podhaisky: Numerik gewöhnlicher Differentialgleichungen. 2. Auflage. Springer Spektrum, Wiesbaden 2012, ISBN 978-3-8348-1847-8, S. 57–59.
26. Peter Deuflhard, Folkmar Bornemann: Numerische Mathematik 2 – Gewöhnliche Differentialgleichungen. 3. Auflage. Walter de Gruyter, Berlin 2008, ISBN 978-3-11-020356-1, S. 199–204.
27. Robert Plato: Numerische Mathematik kompakt. 4. Auflage. Vieweg+Teubner, Wiesbaden 2010, ISBN 978-3-8348-1018-2, Kap. 7: Einschrittverfahren für Anfangswertprobleme, S. 173–177.
28. Karl Strehmel, Rüdiger Weiner, Helmut Podhaisky: Numerik gewöhnlicher Differentialgleichungen. 2. Auflage. Springer Spektrum, Wiesbaden 2012, ISBN 978-3-8348-1847-8, S. 64–70.
29. ode45: Solve nonstiff differential equations — medium order method. MathWorks, abgerufen am 23. November 2017 (englisch).