Kurs:Gewöhnliche Differentialgleichungen/3.2 Schrittweitensteuerung bei Einschrittverfahren

Schrittweitensteuerung bei Einschrittverfahren

Die Anpassung der Schrittweite $h$ während der numerischen Berechnung der Lösung spielt eine wichtige Rolle bei den Anfangswertaufgaben, deren rechte Seite eine steile Funktion ist. Wird die Schrittweite in solchen Fällen fest gehalten, wird der numerische Fehler in den Bereichen des steilen Anstiegs von $f$ stark wachsen und das numerische Verfahren versagen. Schrittweitensteuerung ist essenziel auch für die sogenannten steifen Aufgaben mit sehr unterschiedlichen Wachstumsraten der Lösungskomponenten. Die steifen Probleme werden noch später näher studiert. In diesem Abschnitt wird beschrieben, wie man den genauen Ausdruck für den lokalen Fehler eines Verfahrens für die Anpassung der Schrittweite $h$ einsetzen kann. Wir beschränken uns hier auf das explizite Euerverfahren.

Der lokale Fehler beschreibt, welchen Fehler das ESV in jedem einzelnen Berechnungsschritt zulässt, siehe Definition 3.8, oder wie gut die Verfahrensfunktion die rechte Seite der AWA, $f(t,y)$, annähert. Im Fall des expliziten Eulerverfahrens ist der lokaler Fehler, ausgewertet an der Stelle $t_i$ $\begin{array}{ccc}\tau (t_i,h)& =& {\frac{y}{\text{'}}}^{\prime }(\xi ),\xi \in (t_i,t_i+h),\text{bzw.,}\\ \tau (t_i,h)& =& {\frac{y}{\text{'}}}^{\prime }(t_i)+𝒪(h^2).\end{array}(3.16)$

Je steiler die Funktion $f$ ist, (je größer $f^{\prime }(t_i,y(t_i))=y^{″}(t_i)$), desto größer wird der lokale Fehler $\tau (t_i,h)$. Das Ziel der Schrittweitensteuerung ist die Schrittweite an den Verlauf der Funktion $f$ anzupassen, um den lokalen Fehler in jedem Berechnungschritt $i=1,\dots N_h-1$ unter einer gewissen Toleranz zu halten.

Um die passende Schrittweite in jedem Schritt zu finden, müssen wir die numerische Lösung und den zukünftigen Fehler vorab ”abtasten”. Bei explizitem Eulerverfahren bietet sich für diese Abschätzung der Abgleich der numerischen Lösung mit der Lösung mithilfe zweier ,,halber” Euler-Schritte (mit Schritweite $\frac{h}{2}$). Um dies zu erläutern, untersuchen wir zuerst den lokalen Fehler des Eulerverfahrens mit halber Schrittweite: $\begin{array}{ccc}\text{1. Schritt:}& & {\bar{\eta }}_{i+\frac{1}{2}}={\bar{\eta }}_i+\frac{h}{2}f(t_i,{\bar{\eta }}_i)),\\ \text{2. Schritt:}& & {\bar{\eta }}_i={\bar{\eta }}_{i+\frac{1}{2}}+\frac{h}{2}f(t_i+\frac{h}{2},{\bar{\eta }}_{i+\frac{1}{2}})={\bar{\eta }}_i+\frac{h}{2}f(t_i,{\bar{\eta }}_i)+\frac{h}{2}f(t_i+\frac{h}{2},{\bar{\eta }}_{i+\frac{1}{2}}).\end{array}(3.17)$

Der lokale Fehler dieses Verfahrens ist zweimal kleiner als der des Eulerfervahrens, denn $\begin{array}{ccc}h\bar{\tau }(t_i,h)& =& \bar{\epsilon }(t_i,h)=y(t_i+h)-[y(t_i)+\frac{h}{2}f(t_i,y(t_i))+\frac{h}{2}f(t_i+\frac{h}{2},y(t_i)+\frac{h}{2}f(t_i,y(t_i))\left)\right]\\ \stackrel{Taylor}{=}& & y(t_i)+h{y}^{\prime }(t_i)+\frac{h^2}{2}{y}^{″}(t_i)\\ & & -[y(t_i)+\frac{h}{2}f(t_i,y(t_i))+\frac{h}{2}(f+\frac{h}{2}(f_t+f{f}_y))(t_i,y(t_i))]+𝒪(h^3)\\ & =& \frac{h^2}{2}{y}^{″}(t_i)-\frac{h^2}{4}\underset{=y^{″}(t_i)}{\underbrace{(f_t+f{f}_y)(t_i,y(t_i))}}+𝒪(h^3)=\frac{h^2}{4}{y}^{″}(t_i)+𝒪(h^3),\end{array}(3.18)$

vergleiche (3.16). Nach der Definition des lokalen Fehlers ist der lokale Fehler gleich dem Unterschied zwischen der exakten und der numerichen Lösung (berechnet aus exaktem Startwert), $\epsilon (t_i,h)=y(t_i+h)-{\eta }_{i+1}$, siehe Bemerkung 3.1 Punkt ii). Demnach erhalten wir für das Eulerverahren $${\eta }_{i+1}\stackrel{(3.16)}{=}y(t_i+h)-\frac{h^2}{2}{y}^{″}(t_i)+𝒪(h^3),$$

und für das Eulerverfahren mit halber Schrittweite $${\bar{\eta }}_{i+1}\stackrel{(3.18)}{=}y(t_i+h)-\frac{h^2}{4}{y}^{″}(t_i)+𝒪(h^3).$$

Diese Erkenntnis lässt zwei Folgerungen schließen:
a) Der lokale Fehler vom Eulerverfahren lässt sich mithilfe des Unterschiedes der numerischen Lösung mit Schrittweite $h$ und halber Schrittweite vorhersagen, $$2({\bar{\eta }}_{i+1}-{\eta }_{i+1})=\frac{h^2}{2}{y}^{″}(t_i)+\widetilde{𝒪}(h^3).(3.19)$$

b) Mit einer geeigneten Kombination der zwei Lösungen ${\bar{\eta }}_{i+1},{\eta }_{i+1}$ ist es möglich das führende Fehlerglied $\frac{h^2}{2}{y}^{″}(t_i)$ zu eliminieren und eine numerische Lösung zweiter Ordnung zu rekonstruieren, denn es gilt $$2{\bar{\eta }}_{i+1}-{\eta }_{i+1}=y(t_i+h)+𝒪(h^3).$$

Die Gleichung (3.19) ist die Grundlage der Schrittweitenanpassung des Eulerverfahrens. Der lokaler Fehler $\frac{h^2}{2}{y}^{″}(t_i)+𝒪(h^3)$ lässt sich vorab bestimmen, indem man den Unterschied der numerischen Lösungen mit ganzer und halber Schrittweite $2({\bar{\eta }}_{i+1}-{\eta }_{i+1})\equiv \epsilon (t_i,h)$ bildet. Die Schrittweite $h$ wird so gewählt, dass der lokaler Fehler die vorher gewählte Toleranz $tol$ nicht überschreitet:

ALGORITMUS: Schrittweitensteuerung

1. $t=t_i$, wähle eine Probeschrittweite $H$ und die Toleranz $tol$ für den lokalen Fehler.
   

Berechne die numerische Lösung im nächsten Schritt $t_i+H:$

1. • mithilfe eines Schrittes des Eulerverfahrens mit Schrittweite $H$: ${\eta }_{i+1}^H$,
   • mithilfe von zwei Schritten des Eulerverfahrens mit halber Schrittweite $\frac{H}{2}:$, ${\bar{\eta }}_{i+1}^H$, siehe (3.17).
2. Schätze die... führende Fehlerkonstante mittels (3.19), $$|\epsilon (t_i,H)|=2|{\bar{\eta }}_{i+1}^H-{\eta }_{i+1}^H|\approx \gamma H^2,H\in {ℝ}_{+}⟹\gamma \approx \frac{|\epsilon (t_i,H)|}{H^2}\equiv \frac{2|{\bar{\eta }}_{i+1}^H-{\eta }_{i+1}^H|}{H^2}.$$
3. Wähle die optimale Schrittweite $h_{opt}$ so, dass der lokale Fehler $\epsilon (t_i,h_{opt})$ die vorgegebene Toleranz $tol$ nicht überschreitet: $$tol\ge |\epsilon (t_i,h_{opt})|\approx \gamma h_{opt}^2\equiv \frac{|\epsilon (t_i,H)|}{H^2}{h}_{opt}^2⟹$$ $$h_{opt}\le \alpha \sqrt{\frac{tol}{|\epsilon (t_i,H)|}}H=\alpha \sqrt{\frac{tol}{2|{\bar{\eta }}_{i+1}^H-{\eta }_{i+1}^H|}}H.$$ Gleichung (3.20) Hier kann der Parameter $\alpha$ gewählt werden, um die optimale Schrittweite zusätzlich noch zu skalieren. In der Praxis wird of $\alpha =0.9$ angewendet.
4. Bestimme die numerische Lösung ${\eta }_{i+1}$ im nächsten Schritt $t_i+h_{opt}$ mithilfe der Schrittweite $h_{opt}$ und des Verfahrens zweiter Ordnung, ${\eta }_{i+1}:={\bar{\eta }}_{i+1}^{h_{opt}}$, siehe (3.17).
   

Dieser Schritt ist optional, man kann hier auch das klassische explizite Eulerverfahren anwenden, allerdings mit Verlust der Genauigkeit (Konsistenzordnung).

Anstelle der Toleranz für den absoluten lokalen Fehler, wird in der Schrittweitensteuerung oft die Toleranz für den relativen Fehler (%) angenommen. In dem Fall wird statt $\epsilon (t_i,H)$ er relativer Fehler ${\epsilon }_{rel}(t_i,H)=\frac{\epsilon (t_i,H)}{{\eta }_{i+1}}$ für die Schrittberechnung benutzt, $$h_{opt}\le \alpha \sqrt{\frac{tol}{|{\epsilon }_{rel}(t_i,H)|}}H=\alpha \sqrt{\frac{tol.|{\eta }_{i+1}|}{2|{\bar{\eta }}_{i+1}^H-{\eta }_{i+1}^H|}}H.$$