mirror of
https://github.com/vale981/TUD_MATH_BA
synced 2025-03-06 01:51:38 -05:00
361 lines
No EOL
19 KiB
TeX
Executable file
361 lines
No EOL
19 KiB
TeX
Executable file
\section{Einschrittverfahren}
|
|
\subsection{Grundlagen}
|
|
|
|
Anstelle der gesuchten Lösungsfunktion $y$: $[a,b]\to\real^m$ einer AWA ist man an möglichst guten Näherungen $y^k\in\real^m$ ($k=0,1,2,...,N$) für die Funktionswerte $y(x_k)\in\real^m$ der Funktion $y$ an \begriff{Gitterpunkten} $x_k\in [a,b]$ interessiert. Auf Grundlage der Paare $(x_k,y^k)$ ($k=0,1,...,N-1$) ist es auch möglich, eine Näherungsfunktion $y$ zu erzeugen (etwa durch Interpolation).
|
|
|
|
Einschrittverfahren bilden eine Klasse von Verfahren, die Näherungen $y^k$ zu erzeugen. Das \begriff{Gitter} $\{x_0,...,x_N\}$ is so gewählt, dass
|
|
\begin{align}
|
|
x_0=a<x_1<x_2<\dots < x_{N-1} < x_N = b\notag
|
|
\end{align}
|
|
Außerdem setzen wir
|
|
\begin{align}
|
|
h_k = x_{k+1}-x_k\quad\text{für} k=0,...,N-1\notag
|
|
\end{align}
|
|
und bezeichnen $h_k$ als \begriff{Schrittweite}. Falls $h=h_0=\dots=h_{N-1}$, so heißen die Gitterpunkte bzw. das Gitter \begriff[Gitter!]{gleichabständig} oder \begriff[Gitter!]{äquidistant}.
|
|
|
|
Ein Verfahren zur Erzeugung einer Folge $y^0,...,y^N$ heißt \begriff{Einschrittverfahren} für das AWA \cref{3_1_1}, wenn
|
|
\begin{align}
|
|
\label{3_1_3}
|
|
y^{k+1} = y_k + h_k\Phi(x_k,y_k,y^{k+1},h_k) \quad\text{für} k=0,...,N-1
|
|
\end{align}
|
|
Dabei bezeichnet $\Phi(x,y,z,h)$ den Funktionswert einer \begriff{Verfahrensfunktion}
|
|
\begin{align}
|
|
\Phi: [a,b]\times\real^m\times \real^m\times (0,b-a]\to \real^m\notag
|
|
\end{align}
|
|
die das jeweilige Einschrittverfahren definiert. Man beachte, dass $y^0$ bereits durch die Anfangsbedingung in \cref{3_1_1} gegeben ist. Ein Einschrittverfahren heißt \begriff[Einschrittverfahren!]{implizit}, falls $\Phi$ tatsächlich von $z$ abhängt. Dann ist zur Bestimmung von $y^{k+1}$ aus \cref{3_1_3} die Lösung eines im Allgemeinen nichtlinearen Gleichungssystems erforderlich. Falls $\Phi$ nicht von $z$ abhängt, heißt das Einschrittverfahren \begriff[Einschrittverfahren!]{explizit}. Das explizite \begriff{\person{Euler}-Verfahren} (auch \begriff{Polygonzugverfahren} genannt) ist gegeben durch
|
|
\begin{align}
|
|
\label{3_1_4}
|
|
\Phi(x,y,z,h) = f(x,y)
|
|
\end{align}
|
|
das heißt
|
|
\begin{align}
|
|
y^{k+1} = y^k + h_kf(x_k,y^k) \notag
|
|
\end{align}
|
|
Für das implizite \person{Euler}-Verfahren gilt die Vorschrift
|
|
\begin{align}
|
|
y^{k+1} = y^k + h_kf(x_k + h_k,y^{k+1}) \notag
|
|
\end{align}
|
|
Um die Güte der Näherungen $y^k$ zu beurteilen, untersuchen wir zunächst den lokalen Diskretisierungsfehler eines Einschrittverfahrens.
|
|
|
|
\subsection{Lokaler Diskretisierungsfehler und Konsistenz}
|
|
|
|
\begin{definition}[lokaler Diskretisierungsfehler]
|
|
Seien $y$: $[a,b]\to \real^m$ Lösung des Differentialgleichung $y'=f(x,y)$ und $\Phi$ die Verfahrensfunktion eines Einschrittverfahrens. Für $x\in[a,b)$ und $h>0$ mit $x+h\le b$ heißt
|
|
\begin{align}
|
|
\label{3_1_5}
|
|
\Delta(x,h) = y(x+h) - \bigg( y(x) + h\Phi\big(x,y(x),y(x+h),h\big)\bigg)
|
|
\end{align}
|
|
\begriff{lokaler Diskretisierungsfehler} und
|
|
\begin{align}
|
|
\frac{\Delta(x,h)}{h} = \frac{y(x+h)-y(x)}{h} - \Phi(x,y(x),y(x+h),h)
|
|
\end{align}
|
|
relativer lokaler Diskretisierungsfehler des Einschrittverfahrens.
|
|
\end{definition}
|
|
|
|
Der lokale Diskretisierungsfehler gibt also die Abweichung zwischen exakter Lösung $y(x+h)$ an der Stelle $x+h$ und der Näherung an dieser Stelle an, wobei angenommen wird, dass die Näherung unter Verwendung der exakten Lösung $y(x)$ (und ggf. $y(x+h)$) berechnet wird. Die Bezeichnung relativer Diskretisierungsfehler ist bezüglich der Schrittweite $h$ zu verstehen.
|
|
|
|
\begin{definition}[konsistent, Konsistenzordnung]
|
|
Ein Einschrittverfahren heißt \begriff{konsistent} zur Differentialgleichung $y'=f(x,y)$, wenn
|
|
\begin{align}
|
|
\lim\limits_{h\downarrow 0} \left\Vert\frac{\Delta(x,h)}{h}\right\Vert =0\quad\forall x\in [a,b) \notag
|
|
\end{align}
|
|
für jede Lösung $y$: $[a,b]\to\real^m$ der Differentialgleichung gilt. Gibt es außerdem $p\ge 1$, $M>0$, $\tilde{h}>0$, so dass
|
|
\begin{align}
|
|
\left\Vert\frac{\Delta(x,h)}{h}\right\Vert \le Mh^p\quad\forall (x,h)\in [a,b)\times (0,\tilde{h})\text{ mit } x+h\le b \notag
|
|
\end{align}
|
|
für jede Lösung $y$: $[a,b]\to\real^m$ der Differentialgleichung gilt, so hat das Einschrittverfahren (für diese Differentialgleichung) die \begriff{Konsistenzordnung} $p$.
|
|
\end{definition}
|
|
|
|
\begin{proposition}
|
|
\proplbl{3_2_3}
|
|
Sei $f$: $[a,b]\times \real^m\to\real^m$ stetig differenzierbar. Dann hat das explizite \person{Euler}-Verfahren die Konsistenzordnung 1.
|
|
\end{proposition}
|
|
\begin{proof}
|
|
Mit \cref{3_1_4} folgt
|
|
\begin{align}
|
|
\Delta(x,h) = y(x+h) - y(x) - hf(x,y(x)) \notag
|
|
\end{align}
|
|
Da $y$ die Differentialgleichung $y'=f(x,y)$ löst und $f$ stetig differenzierbar ist, muss $y$ zweimal stetig differenzierbar sein. Aus der \person{Taylor}-Formel erhält man für $i\in\{1,...,m\}$
|
|
\begin{align}
|
|
\Delta(x,h)_i &= y'_i(x)h + \frac{1}{2}y''_i(\xi_i(x,h))h^2 - hf_i(x,y(x)) \notag \\
|
|
&= \frac{1}{2} y''_i(\xi_i(x,h))h^2 \notag
|
|
\end{align}
|
|
für ein $\xi_i(x,h)\in (x,x+h)$. Die Stetigkeit von $y''$ auf $[a,b]$ und Division durch $h$ liefert die Behauptung mit $M=\frac{1}{2}\max_{1\le i\le m}\max_{\xi\in[a,b]}\Vert y''_i(\xi)\Vert$ und $\tilde{h}=b-a$.
|
|
\end{proof}
|
|
|
|
\subsection{Konvergenz von Einschrittverfahren}
|
|
|
|
Zum Gitter $G=\{x_0,...,x_N\}\subset [a,b]$ mit $x_0=a$ und $x_N=b$ seien $y^0,...,y^N\in\R^m$ durch ein Einschrittverfahren erzeugt. Weiter bezeichne $y$: $[a,b]\to\R^m$ die eindeutige Lösung der AWA \cref{3_1_1}. Dann seien
|
|
\begin{align}
|
|
e(x_k) &= y(x_k)-y^k \notag \\
|
|
e(G) &= \max_{x\in G}\norm{e(x)} \notag
|
|
\end{align}
|
|
sowie
|
|
\begin{align}
|
|
h_{max}(G) = \max_{k=0,...,N-1} h_k\notag
|
|
\end{align}
|
|
definiert.
|
|
|
|
\begin{definition}[konvergent]
|
|
Die AWA \cref{3_1_1} besitzt die eindeutige Lösung $y$: $[a,b]\to\R^m$. Ein Einschrittverfahren für diese AWA heißt dann \begriff[Einschrittverfahren!]{konvergent}, falls
|
|
\begin{align}
|
|
\lim_{l\to\infty} e(G_l)=0\notag
|
|
\end{align}
|
|
für alle Gitterfolgen $\{G_t\}$ gilt, für die $\lim_{l\to\infty} h_{max}(G_l)=0$. Gibt es außerdem $p\ge 1$, $C>0$, $\tilde{h}>0$, so dass
|
|
\begin{align}
|
|
e(G) \le C\cdot h_{max}(G)^p\notag
|
|
\end{align}
|
|
für jedes Gitter mit $h_{max}(G)\le\tilde{h}$, so hat das Einschrittverfahren für die gegebene AWA die \begriff[Einschrittverfahren!]{Konvergenzordnung} $p$.
|
|
\end{definition}
|
|
|
|
\begin{lemma}[diskretes \person{Grönwall}sches Lemma]
|
|
\proplbl{3_2_5}
|
|
Falls die Zahlenfolgen $\{\alpha_k\}$, $\{\beta_k\}$, $\{v_k\}\subset [0,\infty)$ den Bedingungen
|
|
\begin{align}
|
|
v_0=0\quad\text{und}\quad v_{k+1}=(1+\alpha_k)v_k+\beta_k\quad\forall k=0,...,N-1\notag
|
|
\end{align}
|
|
genügen, dann folgt
|
|
\begin{align}
|
|
v_{k+1} \le \sum_{i=0}^{k} \beta_i\cdot \exp\left(\sum_{j=i+1}^{k} \alpha_j\right)\quad\text{für } k=0,...,N-1\notag
|
|
\end{align}
|
|
gilt zusätzlich $\alpha_k=\alpha>0$ und $\beta_k=\beta>0$ für jedes $k=0,...,N-1$, dann folgt
|
|
\begin{align}
|
|
v_k \le \frac{\beta}{\alpha}(\exp(k\alpha)-1)\quad\text{für } k=0,...,N-1\notag
|
|
\end{align}
|
|
\end{lemma}
|
|
\begin{proof}
|
|
Zum Beispiel durch vollständige Induktion (vgl. Übungsaufgabe).
|
|
\end{proof}
|
|
|
|
In der Literatur findet man für vorstehende und ähnliche Aussagen die Bezeichnung \textit{diskretes \person{Grönwall}sches Lemma}.
|
|
|
|
\begin{proposition}
|
|
Die AWA \cref{3_1_1} besitze die eindeutige Lösung $y$: $[a,b]\to\R^m$. Ein Einschrittverfahren mit der Verfahrensfunktion $\Phi$ habe für die Differentialgleichung $y'=f(x,y)$ die Konsistenzordnung $p$. Es gebe ferner $L_\Phi>0$ und $H>0$, so dass die Lipschitz-Bedingung
|
|
\begin{align}
|
|
\label{3_1_6}
|
|
\norm{\Phi(x,y,z,h) - \Phi(x,\tilde{y},\tilde{z},h)} \le L_\Phi(\norm{y-\tilde{y}} + \norm{z-\tilde{z}})
|
|
\end{align}
|
|
für alle $(x,y,z,h),(x,\tilde{y},\tilde{z},h)\in [a,b]\times \R^m\times \R^m\times (0,H]$ gilt. Dann besitzt das Einschrittverfahren die Konvergenzordnung $p$.
|
|
\end{proposition}
|
|
\begin{proof}
|
|
Entsprechend \cref{3_1_3} und \cref{3_1_5} gilt
|
|
\begin{align}
|
|
y^{k+1} = y^k + h_k\Phi(x_k,y^k,y^{k+1},h_k) \notag
|
|
\end{align}
|
|
und
|
|
\begin{align}
|
|
y(x_{k+1}) = y(x_k) + h_k\Phi(x_k,y(x_k),y(x_k+h_k),h_k) + \Delta(x_k,h_k) \notag
|
|
\end{align}
|
|
also folgt
|
|
\begin{align}
|
|
e(x_{k+1}) &= y(x_{k+1}) - y^{k+1} \notag \\
|
|
&= y(x_k) - y^k + h_k\big(\Phi((x_k,y(x_k),y(x_k+h_k),h_k) - \Phi(x_k,y^k,y^{k+1},h_k)\big) + \Delta(x_k,h_k) \notag
|
|
\end{align}
|
|
und weiter mit \cref{3_1_6} für $0<h_k\le \tilde{h}=\min\{H,\tilde{h},\frac{1}{2L_\Phi}\}$
|
|
\begin{align}
|
|
\norm{e(x_{k+1})} \le \norm{e(x_k)} + \norm{\Delta(x_k,h_k)} + h_kL_\Phi(\norm{e(x_k)} + \norm{e(x_{k+1})}) \notag
|
|
\end{align}
|
|
Durch Umstellen und Beachtung der Konsistenzordnung ergibt sich
|
|
\begin{align}
|
|
\label{3_1_7}
|
|
\norm{e(x_{k+1})} \le \frac{1+h_kL_\Phi}{1-h_kL_\Phi}\norm{e(x_k)} + \frac{M}{1-h_kL_\Phi}h_k^{p+1}
|
|
\end{align}
|
|
Mit $\alpha_k=4h_kL_\Phi$ hat man (wegen $2h_kL\Phi\le 1$)
|
|
\begin{align}
|
|
\frac{1+h_kL_\Phi}{1-h_kL_\Phi} = 1 + \frac{2h_kL_\Phi}{1-h_kL_\Phi} \le 1+\alpha_k \notag
|
|
\end{align}
|
|
Setzt man weiter $v_k=\norm{e(x_k)}$ und $\beta_k=2Mh_k^{p+1}$, so erhält man aus \cref{3_1_7}
|
|
\begin{align}
|
|
v_{k+1} \le (1+\alpha_k)v_k + \beta_k\quad k=0,...,N-1\notag
|
|
\end{align}
|
|
Nach \propref{3_2_5} folgt daraus (beachte $v_0=\norm{e(x_0)} = \norm{y(x_0) - y^0} = 0$)
|
|
\begin{align}
|
|
v_{k+1} \le \left(\sum_{i=0}^k \beta_i\right)\exp\left(\sum_{i=0}^k \alpha_i\right) \quad\text{für } k=0,...,N-1\notag
|
|
\end{align}
|
|
und damit
|
|
\begin{align}
|
|
\norm{e(x_{k+1})} = v_{k+1} &\le 2M\left(\sum_{i=0}^k h_i^{p+1}\right)\exp\left(4L_\Phi\sum_{i=0}^k h_i\right) \notag \\
|
|
&\le 2Mh_{max}(G)^p(x_{k+1}-x_0)\exp(4L_\Phi(x_{k+1} - x_0)) \notag
|
|
\end{align}
|
|
für $k=0,...,N-1$ sowie
|
|
\begin{align}
|
|
e(G) \le 2M(b-a)\exp(4L_\Phi(b-a))h_{max}(G)^p \notag
|
|
\end{align}
|
|
Also besitzt das Einschrittverfahren die Konvergenzordnung $p$.
|
|
\end{proof}
|
|
|
|
\subsection{Stabilität gegenüber Rundungsfehlern}
|
|
|
|
Wir betrachten das Einschrittverfahren \cref{3_1_3} für ein gleichabständiges Gitter ($h_k=h$) bei exakter Rechnung, das heißt
|
|
\begin{align}
|
|
\label{3_1_8}
|
|
y^{k+1} = y^k + h\Phi(x_k,y^k,y^{k+1},h)\quad\text{für } k=0,...,N-1
|
|
\end{align}
|
|
Weiter beschreibe
|
|
\begin{align}
|
|
\label{3_1_9}
|
|
\tilde{y}^0 = y^0\quad\text{und}\quad \tilde{y}^{k+1}=\tilde{y}^k + h\Phi(x_k,\tilde{y}^k,\tilde{y}^{k+1},h)+\epsilon_k \quad\text{für } k=0,...,N-1
|
|
\end{align}
|
|
ein gestörtes Verfahren, das heißt $\tilde{y}^1,...,\tilde{y}^N$ sind die tatsächlich im Computer berechneten Größen.
|
|
|
|
\begin{proposition}
|
|
\proplbl{3_2_7}
|
|
Sei $y^o\in\R^m$ gegeben und $y^1,...,y^N$ sowie $\tilde{y}^1,...,\tilde{y}^N$ entsprechend \cref{3_1_8} und \cref{3_1_9} berechnet, wobei $\norm{\epsilon_k}<\epsilon$ für alle $k=0,...,N-1$ mit einem $\epsilon>0$ gelte. Außerdem sei für gewisse $L_\Phi,H>0$ die Lipschitz-Bedingung \cref{3_1_6} für alle $(x,y,z,h),(x,\tilde{y},\tilde{z},h)\in [a,b]\times \R^m\times \R^m\times (0,H]$ erfüllt. Dann gibt es $\tilde{h}>0$, so dass
|
|
\begin{align}
|
|
\norm{y^k-\tilde{y}^k} \le \frac{\epsilon}{2hL_\Phi}(\exp(4L_\Phi(x_k-a)) - 1)\quad\text{für } k=0,...,N\notag
|
|
\end{align}
|
|
falls $0>h<\tilde{h}$.
|
|
\end{proposition}
|
|
\begin{proof}
|
|
Für $z^k=y^k-\tilde{y}^k$ folgt
|
|
\begin{align}
|
|
z^{k+1} &= y^{k+1} - \tilde{y}^{k+1} \notag \\
|
|
&= y^k - \tilde{y}^k + h\big(\Phi(x_k,y^k,y^{k+1},h) - \Phi(x_k,\tilde{y}^k,\tilde{y}^{k+1},h)\big) - \epsilon_k \notag
|
|
\end{align}
|
|
und weiter
|
|
\begin{align}
|
|
\norm{z^{k+1}} \le \norm{z^k} + hL_\Phi(\norm{z^k} + \norm{z^{k+1}}) + \epsilon \notag
|
|
\end{align}
|
|
Mit $v_k=\norm{z^k}$, $\alpha=4hL_\Phi$ und $\beta=2\epsilon$ hat man für $0<h\le\tilde{h}=\min{H,\frac{1}{2L_\Phi}}$ die Differenzenungleichung
|
|
\begin{align}
|
|
v_{k+1} \le (1+\alpha)v_k + \beta \notag
|
|
\end{align}
|
|
für $k=0,...,N-1$. \propref{3_2_5} liefert
|
|
\begin{align}
|
|
\norm{y^k-\tilde{y}^k} = \norm{z^k} &= v_k\notag \\
|
|
&\le \frac{\epsilon}{2hL_\Phi}(\exp(k4hL_\Phi) - 1) \notag \\
|
|
&= \frac{\epsilon}{2hL_\Phi}(\exp(4L_\Phi(x_k-a)) - 1) \notag
|
|
\end{align}
|
|
\end{proof}
|
|
|
|
Selbst wenn die Abschätzung in \propref{3_2_7} nicht scharf ist, muss man damit rechnen, dass der Rundungsfehler wie $\sfrac{1}{h}$ wächst. Der Gesamtfehler eines Einschrittverfahrens ans einer Stelle $x_k$ setzt sich aus dem Verfahrensfehler $\norm{y(x_k)-y^k}$ und dem Rundungsfehler $\norm{y^k-\tilde{y}^k}$ zusammen. Für ein Verfahren der Konvergenzordnung $p$ ergibt sich also (bei äquidistantem Gitter) für den Gesamtfehler
|
|
\begin{align}
|
|
\norm{y(x_k)-\tilde{y}^k} &\le \norm{y(x_k)-y^k} + \norm{y^k-\tilde{y}^k} \notag \\
|
|
&\le Ch^p + \tilde{C}\frac{\epsilon}{h} \notag
|
|
\end{align}
|
|
Minimiert man die rechte Seite der Abschätzung in Abhängigkeit von $h$, so folgt, dass man $h$ nicht kleiner als
|
|
$\sim \sqrt[p+1]{\epsilon}$ wählen sollte. Setzt man speziell $h=\sqrt[p+1]{\epsilon}$, dann folgt
|
|
\begin{align}
|
|
Ch^p + \tilde{C}\frac{\epsilon}{h} = C\exp\left(\frac{p}{p+1}\right) + \tilde{C}\exp\left(\frac{p}{p+1}\right) \notag
|
|
\end{align}
|
|
Durch Erhöhung der Konvergenzordnung $p$ kann man also versuchen, mit einer größeren Schrittweite einen kleineren Gesamtfehler zu erreichen. Ein weiterer Grund für das Interesse an Verfahren mit höherer Konvergenzordnung liegt in der Möglichkeit, die Gesamtzahl der erforderlichen Funktionswertbestimmungen der Funktion $f$ zu verringern.
|
|
|
|
\subsection{\person{Runge-Kutta}-Verfahren}
|
|
|
|
Die Klasse der \person{Runge-Kutta}-Verfahren (RKV) ist eine Möglichkeit, Einschrittverfahren mit höheren Konsistenz- bzw. Konvergenzordnungen zu konstruieren. Betrachten wir folgende Idee, eine Näherung $y^{k+1}$ für $y(x_{k+1})$ aus einer Näherung $y^k$ für $y(x_k)$ zu erzeugen.
|
|
|
|
Wegen $y'=f(x,y)$ liefert der Hauptsatz der Differential- und Integralrechnung
|
|
\begin{align}
|
|
\label{3_1_10}
|
|
y(x_{k+1}) = y(x_k) + \int_{x_k}^{x_{k+1}} f(x,y(x))\diff x
|
|
\end{align}
|
|
Approximiert man das Integral durch eine gewichtet Summe von Funktionswerten (vgl. \person{Newton-Cotes} Formeln), so ergibt sich die folgende Verfahrensidee
|
|
\begin{align}
|
|
\label{3_1_11}
|
|
y^{k+1} = y^k + h_k\sum_{i=1}^s c_if(s_i,y(s_i))
|
|
\end{align}
|
|
wobei $c_1,...,c_s$ die Gewichte und $s_1,...,s_s$ Stützstellen bezeichnen. Zur Darstellung der Stützstellen sei
|
|
\begin{align}
|
|
s_i = x_k + \alpha_ih_k\quad\text{für } i=1,...,s\notag
|
|
\end{align}
|
|
mit $\alpha_1=0$ und den Parametern $\alpha_2,...,\alpha_s$. Da $y(s_i)$ unbekannt ist, ersetzt man $f(s_i,y(s_i))$ zunächst durch einen Parameter $k^i$, wobei $k^1=f(x_k,y^k)\approx f(x_k,y(x_k))$ gesetzt wird. Um $y(s_i)$ und damit $f(s_i,y(s_i))$ zu approximieren, verwendet man (bei expliziten RKV) den Ansatz
|
|
\begin{align}
|
|
y(s_i) \approx y^k + h_k\sum_{j=1}^{i-1} \beta_{ij}k^j\notag
|
|
\end{align}
|
|
mit Parametern $\beta_{ij}$. Bei sogenannten impliziten RKV läuft die Summation von $j=1$ bis $j=s$ (und mindestens ein $\beta_{ij}$ mit $j\ge 1$ ist ungleich 0). Für die Parameter $\alpha_i$, $k^i$, $\beta_{ij}$ ergibt sich (im expliziten Fall) somit das folgende Gleichungssystem
|
|
\begin{align}
|
|
\label{3_1_12}
|
|
\begin{split}
|
|
k^1 &= f(x_k,y^k) \\
|
|
k^2 &= f(x_k+\alpha_2h_k,y^k+h_k\beta_{21}k^1) \\
|
|
k^3 &= f(x_k+\alpha_3h_k,y_k+h_k(\beta_{31}k^1 + \beta_{32}k^2)) \\
|
|
&\vdots \\
|
|
k^s &= f(x_k+\alpha_sh_k,y^k+h_k(\beta_{s1}k^1+\dots+\beta_{s,s-1}k^{s-1}))
|
|
\end{split}
|
|
\end{align}
|
|
Ersetzt man in \cref{3_1_11} die unbekannten Vektoren $f(s_i,y(s_i))$ durch die Näherungen $k^i$, so hat man das $s$-stufige \begriff{\person{Runge-Kutta}-Verfahren}
|
|
\begin{align}
|
|
\label{3_1_13}
|
|
y^{k+1} = y^k + h_k\sum_{i=1}^s c_ik^i
|
|
\end{align}
|
|
mit den Parametern $c_1,...,c_s$. Die Verfahrensfunktion eines expliziten RKV ist damit gegeben durch
|
|
\begin{align}
|
|
\Phi(x,y,h) = \sum_{i=1}^s c_if\left(x+\alpha_ih_iy+h\sum_{j=1}^{i-1}\beta_{ij}k^j(x,y,h)\right) \notag
|
|
\end{align}
|
|
wobei $k^i=k^i(x,y,h)$ entsprechend \cref{3_1_12} verwendet wird (die bei expliziten Verfahren nicht vorhandene Abhängigkeit der Funktion $\Phi$ von $z$ wurde weggelassen). Zum Beispiel ist das explizite \person{Euler}-Verfahren $y^{k+1}=y^k+h_kf(x_k,y^k)$ ein einstufiges RKV mit $c_1=1$.
|
|
|
|
\begin{proposition}
|
|
Sei $f$: $[a,b]\times \R^m\to\R^m$ stetig differenzierbar. Ein explizites RKV \cref{3_1_13} mit
|
|
\begin{align}
|
|
\label{3_1_14}
|
|
\sum_{i=1}^s c_i=1
|
|
\end{align}
|
|
hat dann (mindestens) die Konsistenzordnung 1.
|
|
\end{proposition}
|
|
\begin{proof}
|
|
Sei $x\in\R$ fest gegeben. Weiter sei $\eta=y(x)$. Da $f$ stetig differenzierbar ist, gibt es $L_f>0$, so dass
|
|
\begin{align}
|
|
\norm{f(x,\eta) - f(x+\delta x,\eta + \delta\eta)} \le L_f(\abs{\delta x} + \norm{\delta\eta}) \notag
|
|
\end{align}
|
|
für alle $(\delta x,\delta\eta)\in\R\times\R^m$ mit $\abs{\delta x} + \norm{\delta\eta}\le 1$. Induktiv folgt damit, dass $\tilde{h}>0$ existiert, so dass
|
|
\begin{align}
|
|
\norm{k^i(x,\eta,h) - f(x,\eta)} = \mathcal{O}(h) \quad\forall h\in [0,\tilde{h}] \notag
|
|
\end{align}
|
|
für alle $i=1,...,s$. Also gilt wegen \cref{3_1_14}
|
|
\begin{align}
|
|
\norm{\Phi(x,\eta,h) - f(x,\eta)} = \mathcal{O}(h) \quad\forall h\in [0,\tilde{h}] \notag
|
|
\end{align}
|
|
Daraus erhält man (da $f$ in $[a,b]$ stetig differenzierbar und somit $y$ zweimal stetig differenzierbar ist, vgl. Beweis zu \propref{3_2_3})
|
|
\begin{align}
|
|
\norm{\frac{\Delta(x,h)}{h}} &= \norm{\frac{y(x+h) - y(x)}{h} - f(x,y(x)) + f(x,y(x)) - \Phi(x,y(x),h)} \notag \\
|
|
&\le \mathcal{O}(h) \notag
|
|
\end{align}
|
|
\end{proof}
|
|
|
|
\begin{proposition}
|
|
Sei $f$: $[a,b]\times\R^m\to\R^m$ zweimal stetig differenzierbar. Ein explizites RKV \cref{3_1_13} mit \cref{3_1_14},
|
|
\begin{align}
|
|
\label{3_1_15}
|
|
\sum_{j=1}^{i-1}\beta_{ij}=\alpha_i\quad\text{für } i=2,...,s
|
|
\end{align}
|
|
und
|
|
\begin{align}
|
|
\label{3_1_16}
|
|
\sum_{i=2}^s c_i\alpha_i = ^\frac{1}{2}
|
|
\end{align}
|
|
hat dann (mindestens) die Konvergenzordnung 2.
|
|
\end{proposition}
|
|
\begin{proof}
|
|
Übungsaufgabe
|
|
\end{proof}
|
|
|
|
Verwendet man zur Approximation des bestimmten Integrals in \cref{3_1_10} die Trapezregel, das heißt
|
|
\begin{align}
|
|
\int_{x_k}^{x_{k+1}} f(x,y(x))\diff x \approx \frac{1}{2}\bigg(f(x_k,y(x_k)) + f(x_k+h_k,y(x_k+h_k))\bigg) \notag
|
|
\end{align}
|
|
und ersetzt man $f(x_k,y(x_k))$ und $f(x_k+h_k,y(x_k+h_k))$ im RKV durch $k^1=f(x_k,y^k)$ bzw. $k^2=f(x_k+h_k,y^k+h_kk^1)$, dann ergibt sich ein 2-stufiges RKV mit
|
|
\begin{align}
|
|
c_1=c_2=\frac{1}{2},\quad \alpha_2=1,\quad \beta_{21}=1\notag
|
|
\end{align}
|
|
das heißt die Bedingungen \cref{3_1_14}, \cref{3_1_15} und \cref{3_1_16} sind erfüllt. Also besitzt dieses explizite RKV die Konsistenzordnung 2. es ist als \begriff{Verfahren von \person{Heun}} bekannt.
|
|
|
|
Verwendet man zur Quadratur des Integrals in \cref{3_1_10} die \person{Simpson}-Regel, das erhält man ein 4-stufiges RKV mit folgenden Parametern (im sogenannten \begriff{\person{Butcher}-Schema})
|
|
\begin{align}
|
|
\begin{array}{c|cccc}
|
|
0 & & & & \\
|
|
\alpha_2 & \beta_{21} &&& \\
|
|
\alpha_3 & \beta_{31} & \beta_{32} && \\
|
|
\alpha_4 & \beta_{41} & \beta_{42} & \beta_{43} & \\
|
|
\hline
|
|
& c_1 & c_2 & c_3 & c_4
|
|
\end{array} =
|
|
\begin{array}{c|cccc}
|
|
0 & & & & \\
|
|
\sfrac{1}{2} & \sfrac{1}{2} &&& \\
|
|
\sfrac{1}{2} & 0 & \sfrac{1}{2} && \\
|
|
1 & 0 & 0 & 1 & \\
|
|
\hline
|
|
& \sfrac{1}{6} & \sfrac{1}{3} & \sfrac{1}{3} & \sfrac{1}{6}
|
|
\end{array}\notag
|
|
\end{align}
|
|
Dieses Verfahren hat die Konsistenzordnung 4 (sofern $f$ hinreichend glatt). |