Součásti dokumentu 01DIFR
Zdrojový kód
%\wikiskriptum{01DIFR}
\section{Numerické řešení počátečních úloh pro obyčejné
diferenciální rovnice}
Hledáme řešení rovnice $y'=f(x,y)$, $y(x_0)=y_0$. Řekneme, že úloha je
numericky vyřešena, právě když se podaří sestavit řešení ve tvaru
$y(x_0+h)\doteq y(x_0)+\Delta y_0(x_0,y_0,h)$.
\subsection{Analytická metoda}
Provedeme Taylorův rozvoj funkce $y$:
\[y(x_0+h)-y(x_0)=hy'(x_0)+\frac{h^2}{2}y''(x_0)+\cdots\]
Dále platí
\[y'(x_0)=f(x_0,y_0)=f_0\]
Tento vztah můžeme dál derivovat
\[y''(x_0)=\frac{\pd f}{\pd x}(x_0,y_0)+
\frac{\pd f}{\pd y}(x_0,y_0)y'(x_0)=f_x+f_yf_0,\]
\[
\begin{split}
y'''(x_0)&=\frac{\pd^2f}{\pd x^2}(x_0,y_0)+
2\frac{\pd^2f}{\pd x\pd y}(x_0,y_0)y'(x_0)+\\
&\quad +\frac{\pd^2f}{\pd y^2}(x_0,y_0)y'^2(x_0)+
\frac{\pd f}{\pd y}(x_0,y_0)y''(x_0)=\\
&=f_{x^2}+2f_{xy}f_0+f_{y^2}f_0^2+f_y(f_x+f_yf_0)
\end{split}
\]
a tak lze pokračovat libovolně dlouho (za předpokladu, že $f$ má
derivace dostatečně vysokého řádu).
\subsection{Runge-Kuttovy metody}
Předchozí metoda je výpočetně náročná, proto se v~praxi používají
následující metody, kde přírůstek hledáme ve tvaru
\[\Delta y_0=p_1k_1(h)+p_2k_2(h)+\dots+p_rk_r(h),\]
kde $k_i(h)=hf(\xi_i(h),\eta_i(h))$, $\xi_i(h)=x_0+\alpha_i h$,
$\eta_i=y_0+\beta_{i1}k_1(h)+\beta_{i2}k_2(h)+\dots+\beta_{i,i-1}k_{i-1}(h)$, $\alpha_1=0$.
\begin{align*}
k_1(h)&=hf(x_0,y_0)\\
k_2(h)&=hf(x_0+\alpha_2h,y_0+\beta_{21}k_1(h))\\
k_3(h)&=hf(x_0+\alpha_3h,y_0+\beta_{31}k_1(h)+\beta_{32}k_2(h))\\
&\vdots\\
k_r(h)&=hf(x_0+\alpha_rh,y_0+\beta_{r1}k_1(h)+\dots+\beta_{r,r-1}k_{r-1}(h))
\end{align*}
\[y(x_0+h)\doteq y(x_0)+\underbrace{
p_1k_1(h)+p_2k_2(h)+\dots+p_rk_r(h)
}_{\text{Runge-Kuttovský přírůstek}}.\]
Pod pojmem \uv{skutečný přírůstek} rozumíme $y(x_0+h)-y(x_0)$. Zbývá
vyřešit volbu $\alpha,\beta,p$.
Pokud rozvineme Runge-Kuttovský a skutečný přírůstek v~mocninách $h$,
chceme, aby se rozvoje shodovaly do co možná nejvyšší mocniny
$h$. Budou-li se shodovat až do $h^p$, je chyba řádu $h^{p+1}$. Toho
chceme dosáhnout pro libovolnou volbu pravé strany. Z~toho budeme
vycházet při volbě koeficientů $\alpha,\beta,p$.
Označmerozdíl mezi správnou a spočtenou hodnotou
\[\phi_r(h)=[y(x_0+h)-y(x_0)]-[p_1k_1(h)+p_2k_2(h)+\dots+p_rk_r(h)].\]
Výše uvedená podmínka pak odpovídá podmínce
\[\phi_r(0)=\phi_r'(0)=\dots=\phi_r^{(s)}(0)=0.\]
Pro $r=1$:
\[\phi_1(h)=[y(x_0+h)-y(x_0)]-[p_1hf(x_0,y_0)]\]
\[\phi_1'(0)=y'(x_0)-p_1f_0=f_0-p_1f_0=(1-p_1)f_0\]
z~toho vychází podmínka $p_1=1$.
Výsledná metoda
\[y(x_0+h)\doteq y(x_0)+hf(x_0,y_0)\]
se označuje jako Eulerova.
\index{numerická metoda, Eulerova}
Pro $r=2$:
\[\phi_2(h)=[y(x_0+h)-y(x_0)]-[p_1k_1(h)+p_2k_2(h)]\]
\begin{align*}
k_1(h)&=hf(x_0,y_0)\\
k_2(h)&=hf(x_0+\alpha_2h,y_0+\beta_{21}k_1(h))\\
k_1'(h)&=f(x_0,y_0)\\
k_2'(h)&=f(x_0+\alpha_2h,y_0+\beta_{21}k_1(h))+h\left[
\frac{\pd f}{\pd x}\alpha_2+\frac{\pd f}{\pd y}\beta_{21}k_1'(h)
\right]\\
k_1''(h)&=0\\
k_2''(h)&=2\left[
\frac{\pd f}{\pd x}(x_0+\alpha_2h,y_0+\beta_{21}k_1(h))+
\frac{\pd f}{\pd y}\beta_{21}k_1'(h)
\right]+h[\cdots]'
\end{align*}
\begin{align*}
k_1'(0)&=f_0 & k_1''(0)&=0\\
k_2'(0)&=f_0 & k_2''(0)&=2(\alpha_2f_x+\beta_{21}f_yf_0)
\end{align*}
\begin{align*}
y'(x)&=f(x,y) & y'(x_0)&=f_0\\
y''(x)&=\frac{\pd f}{\pd x}(x,y)+\frac{\pd f}{\pd y}(x,y)y'(x)
& y''(x_0)&=f_x+f_yf_0
\end{align*}
\begin{align*}
0 = \phi_2'(0)&=f_0-[p_1f_0+p_2f_0]=[1-p_1-p_2]f_0\\
0 = \phi_2''(0)&=f_x+f_yf_0-2p_2(\alpha_2f_x+\beta_{21}f_yf_0)=
[1-2\alpha_2p_2]f_x+[1-2\beta_{21}p_2]f_yf_0
\end{align*}
Z~toho vyplývá podmínka
\begin{align*}
p_1+p_2&=1\\
2\alpha_2p_2&=1\\
2\beta_{21}p_2&=1.
\end{align*}
V~praxi se užívají následující volby:
\begin{align*}
\alpha_2&=\beta_{21}=1 & \alpha_2&=\beta_{21}=\frac12\\
p_1&=p_2=\frac12 & p_1&=0 \quad p_2=1\\
k_1&=hf(x_0,y_0) & k_1&=hf(x_0,y_0)\\
k_2&=hf(x_0+h,y_0+h) & k_2&=hf\left(x_0+\frac h2,y_0+\frac{k_1}{2}\right)\\
y(x_0+h)&\doteq y(x_0)+\frac12(k_1+k_2) & y(x_0+h)&\doteq y(x_0)+k_2
\end{align*}
Pro $r=3$:
\[\phi_3(h)=[y(x_0+h)-y(x_0)]-[p_1k_1(h)+p_2k_2(h)+p_3k_3(h)]\]
\begin{align*}
k_1(h)&=hf(x_0,y_0)\\
k_2(h)&=hf(x_0+\alpha_2h,y_0+\beta_{21}k_1)\\
k_3(h)&=hf(x_0+\alpha_2h,y_0+\beta_{31}k_1+\beta_{32}k_2)\\
k_1'(h)&=f(x_0,y_0)\\
k_2'(h)&=f(x_0+\alpha_2h,y_0+\beta_{21}k_1)+h\left[
\frac{\pd f}{\pd x}\alpha_2+\frac{\pd f}{\pd y}\beta_{21}k_1'
\right]\\
k_3'(h)&=f(x_0+\alpha_2h,y_0+\beta_{31}k_1+\beta_{32}k_2)+
h\left[\frac{\pd f}{\pd x}+\frac{\pd f}{\pd y}
(\beta_{31}k_1'+\beta_{32}k_2')\right]\\
k_1''(h)&=0\\
k_2''(h)&=2\left[
\frac{\pd f}{\pd x}(x_0+\alpha_2h,y_0+\beta_{21}k_1)\alpha_2+
\frac{\pd f}{\pd y}\beta_{21}k'
\right]+h[\cdots]'\\
k_3''(h)&=2\left[
\frac{\pd f}{\pd x}
(x_0+\alpha_2h,y_0+\beta_{31}k_1+\beta_{32}k_2)
\alpha_3+\frac{\pd f}{\pd y}(\beta_{31}k_1'+\beta_{32}k_2')
\right]+h[\cdots]'\\
k_1'''(h)&=0\\
k_2'''(h)&=3\left[
\frac{\pd^2f}{\pd x^2}(x_0+\alpha_2h,y_0+\beta_{21}k_1)
\alpha_2^2+2\frac{\pd^2f}{\pd x\pd y}\alpha_2\beta_{21}k_1'+
\frac{\pd^2f}{\pd y^2}(\beta_{21}k_1')^2
\right]+\\
&\quad +h[\cdots]''\\
k_3'''(h)&=3\left[
\frac{\pd^2f}{\pd x^2}
(x_0+\alpha_2h,y_0+\beta_{31}k_1+\beta_{32}k_2)\alpha_3^2+
2\frac{\pd^2f}{\pd x\pd y}\alpha_3(\beta_{31}k_1'+\beta_{32}k_2')+\right.\\
&\quad\left.+\frac{\pd^2f}{\pd y^2}(\beta_{31}k_1'+\beta_{32}k_2')+
\frac{\pd f}{\pd y}(\beta_{31}k_1''+\beta_{32}k_2'')
\right]+h[\cdots]''
\end{align*}
\begin{align*}
k_1'(0)&=f_0\\
k_1''(0)&=0\\
k_1'''(0)&=0\\
k_2'(0)&=f_0\\
k_2''(0)&=2(\alpha_2f_x+\beta_{21}f_yf_0)\\
k_2'''(0)&=3(\alpha_2^2f_{x^2}+2\alpha_2\beta_{21}f_{xy}f_0+
\beta_{21}^2f_{y^2}f_0^2)\\
k_3'(0)&=0\\
k_3''(0)&=2(\alpha_3f_x+(\beta_{31}+\beta_{32})f_yf_0\\
k_3'''(0)&=3(\alpha_3^2f_{x^2}+2\alpha_3(\beta_{31}+\beta_{32})
f_{xy}f_0+(\beta_{31}+\beta_{32})^2f_{y^2}f_0^2+\\
&\quad+2\beta_{32}f_y(\alpha_2f_x+\beta_{21}f_yf_0))
\end{align*}
\begin{align*}
y'(x)&=f(x,y)\\
y''(x)&=\frac{\pd f}{\pd x}(x,y)+\frac{\pd f}{\pd y}(x,y)y'(x)\\
y'''(x)&=\frac{\pd^2 f}{\pd x^2}(x,y)+
2\frac{\pd^2 f}{\pd x\pd y}(x,y)y'(x)+
\frac{\pd^2 f}{\pd y^2}(x,y)y'^2(x)+\frac{\pd f}{\pd y}(x,y)y''(x)\\
y'(x_0)&=f_0\\
y''(x_0)&=f_x+f_yf_0\\
y'''(x_0)&=f_{x^2}+2f_{xy}f_0+f_{y^2}f_0^2+f_y(f_x+f_yf_0)
\end{align*}
\begin{align*}
\phi_3'(0)&=f_0-[p_1f_0+p_2f_0+p_3f_0]=[1-p_1-p_2-p_3]\\
\phi_3''(0)&=f_x+f_yf_0-[2p_2(\alpha_2f_x+\beta_{21}f_yf_0)+
2p_3(\alpha_3f_x+(\beta_{31}+\beta_{32})f_yf_0)]=\\
&=f_x[1-2\alpha_2p_2-2\alpha_3 p_3]+
f_yf_0[1-2p_2\beta_{21}-2p_3(\beta_{31}+\beta_{32})]\\
\phi_3'''(0)&=f_{x^2}[1-3p_2\alpha_2^2-3p_3\alpha_3^2]+
f_{xy}f_0[2-6p_2\alpha_2\beta_{21}-6p_3\alpha_3(\beta_{31}+\beta_{32})]+\\
&\quad+f_0^2f_{y^2}[1-3p_2\beta_{21}^2-3p_3(\beta_{31}+\beta_{32})^2]+
f_yf_x[1-6p_3\alpha_2\beta_{32}]+\\
&\quad+f_y^2f_0[1-6p_3\beta_{21}\beta_{32}],
\end{align*}
z~toho dostáváme
\begin{align*}
1-p_1-p_2-p_3&=0\\
1-2\alpha_2p_2-2\alpha_3 p_3&=0\\
1-2p_2\beta_{21}-2p_3(\beta_{31}+\beta_{32})&=0\\
1-3p_2\alpha_2^2-3p_3\alpha_3^2&=0\\
2-6p_2\alpha_2\beta_{21}-6p_3\alpha_3(\beta_{31}+\beta_{32})&=0\\
1-3p_2\beta_{21}^2-3p_3(\beta_{31}+\beta_{32})^2&=0\\
1-6p_3\alpha_2\beta_{32}&=0\\
1-6p_3\beta_{21}\beta_{32}&=0
\end{align*}
Tyto rovnice jsou ovšem závislé a jsou ekvivalentní s~následující
soustavou:
\begin{align*}
1&=p_1+p_2+p_3\\
\frac12&=p_2\alpha_2+p_3\alpha_3\\
1&=3p_2\alpha_2^2+3p_3\alpha_3^2\\
\alpha_2&=\beta_{21}\\
\alpha_3&=\beta_{31}+\beta_{32}\\
\frac16&=p_3\beta_{32}\alpha_2.
\end{align*}
Každé řešení této soustavy dává Runge-Kuttovu metodu. Čtvrtá derivace
nulovat nejde. V~praxi se používají následující metody:
\begin{align*}
k_1&=hf(x_0,y_0)\\
k_2&=hf\left(x_0+\frac h2,y_0+\frac{k_1}{2}\right)\\
k_3&=hf(x_0+h,y_0-k_1+2k_2)\\
y(x_0+h)&\doteq y(x_0)+\frac16(k_1+4k_2+k_3)
\end{align*}
nebo
\begin{align*}
k_1&=hf(x_0,y_0)\\
k_2&=hf\left(x_0+\frac h3,y_0+\frac{k_1}{3}\right)\\
k_3&=hf\left(x_0+\frac{2h}{3},y_0+\frac{2k_2}{3}\right)\\
y(x_0+h)&\doteq y(x_0)+\frac14(k_1+3k_3).
\end{align*}
Pro $r=4$: Pokud sestavíme polynomy pro $r=4$ a požadovali
$\phi^{(4)}=0$, dojdeme k~následujícím 11 podmínkám:
\begin{align*}
\alpha_2&=\beta_{21}\\
\alpha_3&=\beta_{31}+\beta_{32}\\
\alpha_4&=\beta_{41}+\beta_{42}+\beta_{43}\\
p_1+p_2+p_3+p_4&=1\\
p_2\alpha_2+p_3\alpha_3+p_4\alpha_4&=\frac12\\
p_2\alpha_2^2+p_3\alpha_3^2+p_4\alpha_4^2&=\frac14\\
p_2\alpha_2^3+p_3\alpha_3^3+p_4\alpha_4^3&=\frac14\\
p_3\beta_{32}\alpha_2+p_4\beta_{42}\alpha_2+p_4\beta_{43}\alpha_3&=\frac16\\
p_3\beta_{32}\alpha_2\alpha_3+p_4\beta_{42}\alpha_2\alpha_4
+p_4\beta_{43}\alpha_3\alpha_4&=\frac18\\
p_3\beta_{32}\alpha_2^2+p_4\beta_{42}\alpha_2^2+p_4\beta_{43}\alpha_3^2&=\frac1{12}\\
p_4\beta_{43}\beta_{32}\alpha_2&=\frac1{24}
\end{align*}
Máme 13 neznámých, ale $\phi_4^{(5)}$ už nejde nulovat. Uvedeme
následující tři metody.
\begin{enumerate}
\item Standardní Runge-Kuttova metoda:
\index{numerická metoda, standardní Runge-Kuttova}
\begin{align*}
k_1&=hf(x_0,y_0)\\
k_2&=hf\left(x_0+\frac h2,y_0+\frac{k_1}2\right)\\
k_3&=hf\left(x_0+\frac h2,y_0+\frac{k_2}2\right)\\
k_4&=hf(x_0+h,y_0+k_3)\\
y(x_0+h)&\doteq y(x_0)+\frac16(k_1+2k_2+2k_3+k_4)
\end{align*}
\item Tříosminové pravidlo:
\begin{align*}
k_1&=hf(x_0,y_0)\\
k_2&=hf\left(x_0+\frac h3,y_0+\frac{k_1}3\right)\\
k_3&=hf\left(x_0+\frac{2h}3,y_0-\frac{k_1}3+k_2\right)\\
k_4&=hf(x_0+h,y_0+k_1-k_2+k_3)\\
y(x_0+h)&\doteq y(x_0)+\frac18(k_1+3k_2+3k_3+k_4)
\end{align*}
\item
\begin{align*}
k_1&=hf(x_0,y_0)\\
k_2&=hf\left(x_0+\frac h4,y_0+\frac{k_1}4\right)\\
k_3&=hf\left(x_0+\frac h2,y_0-\frac{k_2}2\right)\\
k_4&=hf(x_0+h,y_0+k_1-2k_2+2k_3)\\
y(x_0+h)&\doteq y(x_0)+\frac16(k_1+4k_3+4k_4)
\end{align*}
\end{enumerate}
Abychom dosáhli maximální chyby $O(h^6)$, nevystačíme s~$r=5$, ale až
s~$r=6$. Proto se používají zejména metody s~$r=4$.
V~praxi se nejčastěji používají metody s~automatickým výběrem
kroku. Nejjednodušší způsob je spočítat následující hodnotu s~krokem
$h$ a $h/2$, je-li relativní rozdíl větší než nějaké $\epsilon$,
interval se dále půlí. Případně se krok může zase prodlužovat ---
např. mám počitadlo, kolikrát odchylka vyhovovala a když dosáhne
určité hodnoty, zase zkusím krok prodloužit.
\begin{example}
Řešme rovnici $y'=y$, $y(0)=1$. Víme, že řešením je $y(x)=e^x$. Pokud
použijeme standardní Runge-Kuttovu metodu, máme
\begin{align*}
k_1&=0.1*1=0.1\\
k_2&=0.1*1.05=0.105\\
k_3&=0.1*1.0525=0.10525\\
k_4&=0.1*1.10525=0.110525\\
y(0.1)&=1+\frac16[0.1+0.21+0.2105+0.110525]=1+\frac16*0.631025\doteq
1.1051708333,
\end{align*}
přičemž $e^{0.1}=1.105170918$.
\end{example}
Kteroukoli z~Runge-Kuttových metod pro rovnici lze použít i pro systém
rovnic.
\begin{example}
Mějme systém rovnic $y'=f(x,y,z)$, $z'=g(x,y,z)$ s~počátečními
podmínkami $y(x_0)=y_0$, $z(x_0)=z_0$. Standardní Runge Kuttova metoda
pro tento systém je
\begin{align*}
k_1&=hf(x_0,y_0,z_0) & l_1&=hg(x_0,y_0,z_0)\\
k_2&=hf\left(x_0+\frac h2,y_0+\frac{k_1}2,z_0+\frac{l_1}2\right)
& l_2&=hg\left(x_0+\frac h2,y_0+\frac{k_1}2,z_0+\frac{l_1}2\right)\\
k_3&=hf\left(x_0+\frac h2,y_0+\frac{k_2}2,z_0+\frac{l_2}2\right)
& l_3&=hg\left(x_0+\frac h2,y_0+\frac{k_2}2,z_0+\frac{l_2}2\right)\\
k_4&=hf(x_0+h,y_0+k_3,z_0+l_3) & l_4&=hg(x_0+h,y_0+k_3,z_0+l_3)
\end{align*}
\begin{align*}
y(x_0+h)&\doteq y(x_0)+\frac16(k_1+2k_2+2k_3+k_4)\\
z(x_0+h)&\doteq z(x_0)+\frac16(l_1+2l_2+2l_3+l_4)
\end{align*}
\end{example}
Rovnici $n$-tého řádu lze převést na systém a pak řešit Runge-Kuttovou
metodou. Pro rovnici druhého řádu
$y''=f(x,y,y')$
existuje přímo Runge-Kuttův vzorec
\begin{align*}
k_1&=\frac{h^2}2 f(x_0,y_0,y_0')\\
k_2&=\frac{h^2}2 f\left(x_0+\frac h2,
y_0+\frac h2y_0'+\frac{k_1}4,
y_0'+\frac{k_1}h\right)\\
k_3&=\frac{h^2}4 f\left(x_0+\frac h2,
y_0+\frac h2y_0'+\frac{k_2}4,
y_0'+\frac{k_2}h\right)\\
k_4&=\frac{h^2}2 f\left(x_0+h,
y_0+hy_0'+k_3,y_0'+\frac{2k_3}h
\right)\\
y(x_0+h)&\doteq y(x_0)+hy'(x_0)+\frac13(k_1+k_2+k_3)\\
y'(x_0+h)&\doteq y'(x_0)+\frac1{3h}(k_1+2k_2+2k_3+k_4)
\end{align*}
\subsection{Řešení lineárních diferenčních rovnic}
Lineární diferenční rovnicí $k$-tého řádu nazýváme rovnici
\[a_k(n)y_{n+k}+a_{k-1}(n)y_{n+k-1}+\dots+a_0(n)y_n=b_n,\]
kde $a_i(n)$, $b(n)$ jsou posloupnosti, $a_k(n)\not=0$,
$a_0(n)\not=0$.
Řešením lineární diferenční rovnice je posloupnost $y_n$, která je
jednoznačně dáno hodnotami $y_0,y_1,\dots,y_{k-1}$, neboť pro $y_i$
pak dostáváme jednoznačný rekurentní vztah. Pokud je $b(n)=0$,
nazýváme rovnici rovnicí bez pravé strany.
\begin{enumerate}
\item Jsou-li $y_n^{(1)},y_n^{(2)},\dots,y_n^{(k)}$ řešení rovnice bez
pravé strany, pak jejich lineární kombinace
\[\sum_{i=1}^k c_iy_n^{(i)}\]
je také řešení.
\item Jsou-li $y_n^{(1)},y_n^{(2)},\dots,y_n^{(k)}$ řešení rovnice bez
pravé strany a platí-li
\[
\begin{vmatrix}
y_0^{(1)} & y_0^{(2)} & \hdots & y_0^{(k)}\\
y_1^{(1)} & y_1^{(2)} & \hdots & y_1^{(k)}\\
\vdots & \vdots & & \vdots\\
y_{k-1}^{(1)} & y_{k-1}^{(2)} & \hdots & y_{k-1}^{(k)}
\end{vmatrix}\not=0
\]
a $Y_n$ je řešení rovnice bez pravé strany, pak existují jednoznačná
$c_1,c_2,\dots,c_k$ taková, že
\[Y_n=\sum_{i=1}^k c_iy_n^{(i)}.\]
\begin{proof}
Matice soustavy rovnic
\begin{align*}
c_1y_0^{(1)}+c_2y_0^{(2)}+\dots+c_ky_0^{(k)}&=Y_0\\
c_1y_1^{(1)}+c_2y_1^{(2)}+\dots+c_ky_1^{(k)}&=Y_1\\
&\vdots\\
c_1y_{k-1}^{(1)}+c_2y_{k-1}^{(2)}+\dots+c_ky_{k-1}^{(k)}&=Y_{k-1}\\
\end{align*}
je regulární, tudíž má jednoznačné řešení.
\end{proof}
\item Obecné řešení rovnice s~pravou stranou má tvar
\[Y_n=\sum_{i=1}^n c_iy_n^{(i)}+p_n,\]
kde $p_n$ je partikulární řešení.
\end{enumerate}
\subsection{Diferenční rovnice s~konstantními koeficienty}
Řešení rovnice tvaru
\[a_ky_{n+k}+a_{k-1}y_{n+k-1}+\dots+a_0y_n=0\]
hledáme ve tvaru $y_n=z^n$. Po dosazení
\[a_kz^{n+k}+a_{k-1}z^{n+k-1}+\dots+a_0z^n=0\]
a po vykrácení $z^n$
\[a_kz^k+a_{k-1}z^{k-1}+\dots+a_0=0\]
dostaneme {\bf charakteristickou rovnici}.
\index{rovnice, charakteristická}
Je-li $z_i$ $r$-násobný kořen, řeší rovnici $z_i^n$ a dále také
$nz_i^n,n^2z_i^n,\dots,n^{r-1}z_i^n$.
\begin{proof}[Nástin důkazu]
Dosazením do diferenční rovnice dostaneme
\[a_k(n+k)z^{n+k}+a_{k-1}(n+k-1)z^{n+k-1}+\dots+a_0nz^n,\]
po úpravě
\[n[\underbrace{a_kz^k+a_{k-1}z^{k-1}+\dots+a_0}_{\text{charakteristický
polynom}}]
+z[\underbrace{ka_kz^{k-1}+(k-1)a_{k-1}z^{k-2}+\dots+a_1}
_{\text{derivace charakteristického polynomu}}]=0.\]
Protože $z$ je $r$-násobný kořen, je kořenem i $1.$, $2.$ až $r-1$-té
derivace.
\end{proof}
\begin{theorem}
Nechť charakteristická rovnice má kořeny $z_1,z_2,\dots,z_s$ vzájemně
různé, násobností $r_1,r_2,\dots,r_s$. Potom řešení
\[z_1^n,nz_1^n,\dots,n^{r-1}_1z_1^n,\dots,
z_s^n,nz_s^n,\dots,n^{r_s-1}z_s^n\]
tvoří fundamentální systém.
\begin{proof}
Předpokládejme, že řešení jsou lineárně závislá. Potom pro každé $n$
musí platit
\[p_1(n)z_1^n+p_2(n)z_2^n+\dots+p_s(n)z_s^n=0,\]
kde $p_1,p_2,\dots,p_s$ jsou polynomy v~$n$ a alespoň jeden z~nich je
nenulový. Dokážeme, že to nemůže platit. To provedeme indukcí podle
počtu nenulových polynomů $s$.
Pro $s=1$: Buď $p_1(n)z_1^n=0$. Po vykrácení $z_1^n$ je $p_1(n)=0$,
což je spor.
\[p_1(n)z_1^n=-p_2(n)z_2^n-\dots-p_s(n)z_s^n\]
Po vydělení $z_1^n$
\[p_1(n)=-p_2(n)\zeta_2^n-\dots-p_s(n)\zeta_s^n,\]
kde $\zeta_i\not=1$ a současně jsou vzájemně různá, neboť $z_i$ jsou
vzájemně různá. Tento vztah musí platit i pro $n+1$.
\[p_1(n+1)=-p_2(n+1)\zeta_2^{n+1}-\dots-p_s(n)\zeta_s^{n+1}.\]
Odečtením dostaneme
\[p_1(n+1)-p_1(n)=-[p_2(n+1)\zeta_2-p_2(n)]\zeta_2^n-\dots
-[p_s(n+1)\zeta_s-p_s(n+1)]\zeta_s^n.\]
Na levé straně máme teď polynom (ostře) nižšího stupně než $p_1$,
protože nejvyšší mocniny se odečetly, naopak stupně polynomů na pravé
straně se díky $\zeta_i\not=1$ nezmění. Tímto způsobem lze postupně
snížit stupeň levé strany až k~nulovému polynomu. Podle indukčního
předpokladu pak jsou všechny polynomy na pravé straně nulové.
\end{proof}
\end{theorem}
\subsection{Jednokrokové metody}
\index{numerické metody, jednokrokové}
\begin{define}
Obecnou jednokrokovou metodou nazveme metodu danou formulí tvaru
\begin{equation}
\label{objednokr}
y_{n+1}=y_n+h\Phi_f(x_n,y_n,h).
\end{equation}
\end{define}
\begin{define}
\index{formule, regulární}
Formule \eqref{objednokr} se nazývá {\bf regulární}, jestliže funkce
$\Phi_f(x,y,h)$ je definována a spojitá na množině, kde $x_0\le x\le
a$, $-\infty<y<+\infty$, $h\in\la 0,h_0\ra$ $(h_0>0)$ a existuje-li
konstanta $M$ (nezávislá na $x$ a $h$) tak, že
\begin{equation}
\label{oj1}
\abs{\Phi_f(x,y,h)-\Phi_f(x,z,h)}\le M\abs{y-z}
\end{equation}
pro každé $x\in\la x_0,a\ra$, $y,z\in(-\infty,+\infty)$ a $h\in\la
0,h_0\ra$.
\end{define}
\begin{define}
Formule \eqref{objednokr} se nazývá {\bf stupně $p$}, existuje-li
konstanta $L$ a $p\in\N$ tak, že
\begin{equation}
\label{oj2}
\abs{\Phi_f(t,y,h)-\frac{r(t+h)-y}{h}}\le Lh^p
\end{equation}
pro každé $t\in\la x_0,a\ra$, $y\in(-\infty,+\infty)$, $h\in\la
0,h_0\ra$, kde $r(x)$ je řešení rovnice $y'=f(x,y)$ na intervalu $\la
t,t+h\ra$ s~počáteční podmínkou $r(t)=y$.
\end{define}
\begin{remark}
$\abs{h\Phi_f(t,y,h)-(r(t+h)-r(t))}\le Lh^{p+1}$.
\end{remark}
\begin{theorem}
Nechť je dána regulární obecná jednokroková formule \eqref{objednokr},
která je stupně $p\in\N$. Nechť $y(x)$ je řešení rovnice $y'=f(x,y)$ a
$y_0,y_1,\dots,y_N$ jsou hodnoty splňující vztah
$y_{n+1}=y_n+h\Phi_f(x_n,y_n,h)+\delta_n$ pro $n=0,1,2,\dots,N-1$. Pak
pro $n=0,1,2,\dots,N$ platí
\[\abs{y_n-y(x_n)}\le\abs{y_0-y(x_0)}e^{M(x_n-x_0)}+
\left(Lh^p+\frac{\delta}{h}\right)
\frac{e^{M(x_n-x_0)}-1}{M},\]
kde
\[\delta=\max_{n=0,\dots,N-1}\abs{\delta_n}\]
a $M$ a $L$ jsou konstanty definované \eqref{oj1} a \eqref{oj2}.
\begin{proof}
Buď $r_n=y_n-y(x_n)$ pro $n=0,1,\dots,N$. Platí, že
\[
\begin{split}
r_{n+1}&=y_{n+1}-y(x_{n+1})=y_n+h\Phi_f(x_n,y_n,h)+\delta_n-y(x_{n+1})=\\
&=\underbrace{(y_n-y(x_n))}_{r_n}+h[\Phi_f(x_n,y_n,h)-\Phi_f(x_n,y(x_n),h)]+\\
&\quad +h\left[\Phi_f(x_n,y(x_n),h)-\frac{y(x_n+h)-y(x_n)}{h}\right]+\delta_n,
\end{split}
\]
s~použitím $\delta\ge\delta_n$ a konstant $M$, $L$ můžeme psát
$\abs{r_{n+1}}\le\abs{r_n}+Mh\abs{r_n}+Lh^{p+1}+\delta$.
Zavedeme $R_0=\abs{r_0}$, $R_{n+1}=(1+hM)R_n+Lh^{p+1}+\delta$,
Platí, že $\abs{r_n}\le R_n$. To dokážeme indukcí:
$\abs{r_{n+1}}\le R_n+MhR_n+Lh^{p+1}+\delta=R_{n+1}$.
Dostali jsme tak diferenční rovnici pro $R_n$. Rovnice bez pravé
strany má tvar $R_{n+1}=(1+hM)R_n$, charakteristická rovnice je
$z-(1+hM)=0$, obecné řešení rovnice bez pravé strany má tvar
$R_n=C(1+hM)^n$.
Hledáme partikulární řešení ve tvaru $P=(1+hM)P+Lh^{p+1}+\delta$, z
čehož dostáváme
\[P=\frac{Lh^{p+1}+\delta}{hM}.\]
Nakonec doladíme konstantu $C$:
\[R_n=C(1+hM)^n-\left(Lh^p+\frac{\delta}{h}\right)\frac1M,\]
\[\abs{r_0}=C-\left(Lh^p+\frac{\delta}h\right)\frac1M,\]
\[C=\abs{r_0}+\left(Lh^p+\frac{\delta}h\right)\frac1M.\]
Obecné řešení rovnice s~pravou stranou má tvar
\[R_n=\abs{r_0}(1+hM)^n+\left(Lh^p+\frac{\delta}h\right)
\frac{(1+hM)^n-1}M.\]
Protože $1+x<e^x$, je
\[\begin{split}
R_n&\le\abs{r_0}e^{nhM}+\left(Lh^p+\frac{\delta}h\right)
\frac{e^{nhM}-1}{M}=\\
&=\abs{y_0-y(x_0)}e^{(x_n-x_0)M}+
\left(Lh^p+\frac{\delta}h\right)\frac{e^{M(x_n-x_0)}-1}{M}.\qed
\end{split}\]
\noqed
\end{proof}
\end{theorem}
\begin{remark}
První člen je způsoben počáteční podmínkou, druhý je chyba metody,
$\delta$ je zaokrouhlovací chyba počítače. Nemá cenu jít s~$h$ pod
jistou mez, protože jinak se bude zvětšovat chyba $\delta/h$. Jediným
způsobem jak zvýšit přesnost je pak metoda vyššího řádu.
\end{remark}
\begin{remark}
Závislost chyby je dost špatná, dá se zkonstruovat rovnice, kdy to
bude nejhorší --- chyba bude exponenciálně závislá.
\end{remark}
\begin{define}
nechť rovnice $y'=f(x,y)$ na intervalu $\la x_0,+\infty)$ má řešení
$y(x)=0$, tj. $0=f(x,0)$. Toto nulové řešení se nazývá {\bf stabilní
vzhledem k~soustavným poruchám}, jestliže pro každé $\epsilon>0$
existuje $\delta>0$ tak, že pro každou spojitou funkci $\eta(x)$,
která pro $x\in\la x_0,+\infty)$ s~výjimkou spočetného množství
izolovaných bodů splňuje rovnici $\eta'=f(x,\eta)+g(x)$, kde
$\abs{\eta(0)}<\delta$, $\abs{g(x)}\le\delta$ pro $x\in\la
x_0,+\infty)$ a $g(x)$ je libovolná měřitelná funkce, platí nerovnost
$\abs{\eta(x)}<\epsilon$ pro $x\in\la x_0,+\infty)$, tj. poškodím-li
počáteční podmínku a pravou stranu o~$\delta$, změní se řešení o
$\epsilon$.
\end{define}
\begin{define}
\index{formule, úplně regulární}
Obecná jednokroková formule se nazývá {\bf úplně regulární}, je-li
funkce $\Phi_f(x,y,h)$ omezená, spojitá ve všech svých proměnných,
stejnoměrně spojitá v~$x$, lipschitzovská v~$y$ s~konstantou $M$
(nezávislou na $x$ a $h$) a stejnoměrně spojitá v~$h$.
\end{define}
\begin{theorem}
Nechť je dána rovnice $y'=f(x,y)$ a na $\la x_0,+\infty)$ platí
$f(x,0)=0$. Nechť nulové řešení je stabilní vzhledem k~soustavným
poruchám. Nechť je dána úplně regulární jednokroková metoda
\eqref{objednokr} stupně $p\in\N$. Pak pro každé $\epsilon>0$ existují
$h_1,\delta>0$ tak, že pro každé řešení diferenční rovnice
$y_{n+1}=y_n+h\Phi_f(x_n,y_n,h)+\delta_n$ pro $n=0,1,2,\dots$, pro
které $\abs{y_0}<\delta$, $\delta_n\le h\delta$, $h<h_1$ splňuje
nerovnost $\abs{y_n}<\epsilon$.
\begin{proof}
Buď
\[\eta(x)=y_n+\frac{y_{n+1}-y_n}{h}(x-x_n)=
y_n+\left(\Phi_f(x_n,y_n,h)+\frac{\delta_n}h\right)(x-x_n)\]
pro $x\in\la x_n,x_{n+1}\ra$, $n=0,1,\dots$. Buď $g(x)=\eta'-f(x,\eta)$,
$\abs{g(x)}\le\delta_1$.
\[
\begin{split}
\abs{\eta'-f(x,\eta)}&\le
\abs{\Phi_f(x_n,y_n,h)+\frac{\delta_n}h-f(x,\eta)}\le
\abs{\Phi_f(x_n,y_n,h)-\Phi_f(x_n,\eta,h)}+\\
&\quad+\abs{\Phi_f(x_n,\eta,h)-\Phi_f(x,\eta,h)}+
\abs{\Phi_f(x,\eta,h)-f(x,\eta)}+\frac{\abs{\delta_n}}h
\end{split}
\]
Protože $\Phi_f$ je lipschitzovské v~$y$, je
\[\abs{\Phi_f(x_n,y_n,h)-\Phi_f(x_n,\eta,h)}\le
M\abs{\eta-y_n}\le M\abs{\Phi_f(x_n,y_n,h)+\frac{\delta_n}h}
\abs{x-x_0}.\] To lze libovolně zmenšit volbou $h_1$ a $\delta$.
Člen $\abs{\Phi_f(x_n,\eta,h)-\Phi_f(x,\eta,h)}$ lze libovolně zmenšit
volbou $h_1$ díky stejnoměrné spojitosti v~$x$.
Buď $r'(x)=f(x,r(x))$, $r(x)=\eta(x)$. Protože
\[\lim_{h\to 0}\abs{\Phi_f(x,\eta,h)-\frac{r(x+h)-\eta}h}=0,\]
a $\Phi_f$ je spojitá v~$x$, je $\Phi_f(x,\eta,0)=r'(x)=f(x,\eta)$. Díky
stejnoměrné spojitosti v~$h$ jde pro dostatečně malé $h$ třetí člen
k~nule.
\end{proof}
\end{theorem}
\subsection{Mnohokrokové (diferenční) metody}
\index{numerické metody, mnohokroké (diferenční)}
Zabýváme se úlohou nalézt přibližné hodnoty řešení rovnice
$y'=f(x,y)$, $y(x_0)=y_0$ v~bodech $x_0,x_1,\dots$, přičemž
$x_i=x_0+ih$. K~výpočtu $y_{n+1}$ potřebujeme hodnoty
$y_{m-k},y_{m-k+1},\dots,y_m$, ale do pravé strany dosazujeme pouze
jednou. Problém je pouze na počátku --- tam se musí použít Runge-Kutta
nebo Taylor.
Vycházíme z~toho, že pro přesné $y_i$ platí
\begin{equation}
\label{mnohokr}
y_{m+1}=y_{m-j}+\int_{x_{m-j}}^{x_{m+1}}f(x,y(x))\d x.
\end{equation}
Při numerickém výpočtu funkci $f(x,y(x))$ nahradíme interpolačním
polynomem k~uzlům $x_{m-k},x_{m-k+1},\dots,x_m$ nebo
$x_{m-k},x_{m-k+1},\dots,x_m,x_{m+1}$. V~prvním případě se jedná o
{\bf explicitní metody} --- dostaneme explicitní vyjádření hodnoty
$y_{m+1}$. V~druhém případě jde o~{\bf implicitní metody} --- ve
vzorci vystoupí i $f(x_{m+1},y_{m+1})$. Poté se buď $y_{m+1}$ vypočte
přímo z~rovnice (méně obvyklé) nebo se vypočítá přibližná hodnota
iteračně.
\index{numerická metoda, explicitní} \index{numerická metoda, implicitní}
Volbou $j=0$ dostaneme tzv. {\bf Adamsovy formule}.
\index{formule, Adamsonova}
\subsubsection{Explicitní Adamsovy formule}
Označme $f_i=f(x_i,y_i)=y'(x_i)=y_i'$. Použijeme Newtonův vzorec pro
interpolaci vpřed, $f(x,y(x))=L_{m,k}(x)+R_{m,k}(x)$, kde
\[L_{m,k}(x)=f_m+tf_{m-\frac12}^1+\frac{t(t+1)}{2!}f_{m-1}^2+
\dots+\frac{t(t+1)\cdots(t+k-1)}{k!}f_{m-\frac{k}2}^k,
\]
\[
\begin{split}
R_{m,k}&=\frac{(x-x_m)(x-x_{m-1})\cdots(x-x_{m-k})}{(k+1)!}
\frac{\d^{k+1}}{\d x^{k+1}}f(\xi,y(\xi))=\\
&=h^{k+1}\frac{t(t+1)\cdots(t+k)}{(k+1)!}
\frac{\d^{k+1}}{\d x^{k+1}}f(\xi,y(\xi)).
\end{split}
\]
Po dosazení do \eqref{mnohokr} máme
\[\underbrace{y_{m+1}=y_m+\int_{x_m}^{x_{m+1}}L_{m,k}(x)\d x}_
{\text{přibližný vzorec pro výpočet $y_{m+1}$}}
+\underbrace{\int_{x_m}^{x_{m+1}}R_{m,k}(x)\d x}_
{\text{chyba v~jednom kroku}},
\]
po substituci $x=x_m+th$
\[y_{m+1}=y_m+h\int_0^1\left[
f_m+tf_{m-\frac12}^1+\dots+\frac{t(t+1)\cdots(t+k-1)}{k!}
f_{m-\frac k2}^k\right]\d t,\]
\[y_{m+1}=y_m+h[f_m+a_1f_{m-\frac12}^2+\dots+a_kf_{m-\frac
k2}^k]+l_{m,k},\]
kde
\[a_i=\int_0^1\frac{t(t+1)\cdots(t+i-1)}{i!}\d t.\]
Některá $a$:
\[
\begin{split}
a_1&=\int_0^1 t\d t=\frac12, a_2=\int_0^1\frac{t(t+1)}{2}\d t=
\left[\frac{t^3}{6}+\frac{t^2}{4}\right]_0^1=\frac{5}{12},\
a_3=\frac38,\ a_4=\frac{251}{720},\\a_5&=\frac{95}{288},\
a_6=\frac{19087}{60480},\ a_7=\frac{5275}{17280},\
a_8=\frac{1070017}{3628800}.
\end{split}\]
Užití Adamsových vzorců:
\[f_i^k=f_{i+\frac k2}-\binom k1f_{i+\frac k2-1}+
\binom k2f_{i+\frac k2-2} -\dots+(-1)^k\binom kkf_{i-\frac k2}\]
\[
\begin{split}
y_{m+1}&\doteq y_m+h\left[y_m'+\frac12(y_m'-y_{m-1}')
+\frac 5{12}(y_m'-2y_{m-1}'+y_{m-2}')+\right.\\
&\left.\quad+\frac 38(y_m'-3y_{m-1}'+3y_{m-2}'-y_{m-3}')\cdots\right]
\end{split}
\]
\[
\begin{split}
k&=0:\ y_{m+1}\doteq y_m+hy_m'\\
k&=1:\ y_{m+1}\doteq y_m+\frac h2[3y_m'-y_{m-1}]\\
k&=2:\ y_{m+1}\doteq y_m+\frac h{12}[23y_m'-16y_{m-1}'+5y_{m-2}']\\
k&=3:\ y_{m+1}\doteq y_m+\frac h{24}[55y_m'-59y_{m-1}'+
37y_{m-2}'-9y_{m-3}'].
\end{split}
\]
Chyba $l_{m,k}$ bude
\[
\begin{split}
l_{m,k}&=\int_{x_m}^{x_{m+1}}R_{m,k}(x)\d x=h^{k+2}\int_0^1
\frac{t(t+1)\cdots(t+k)}{(k+1)!}\frac{\d^{k+1}}{\d x^{k+1}}
f(\xi,y(\xi))\d t=\\&=h^{k+2}y^{(k+2)}(\eta)a_{k+2}.
\end{split}
\]
Zkráceně $\abs{l_{m+k}}\le h^{k+2}a_{k+1}M_{k+2}$.
\subsubsection{Implicitní Adamsovy formule}
Opět použijeme vzorec pro interpolaci vpřed,
$f(x,y(x))=L_{m,k}(x)+R_{m,k}(x)$,
\[L_{m,k}(x)=f_{m+1}+tf_{m+\frac12}^1+\frac{t(t+1)}{2!}f_m^2+
\cdots+
\frac{t(t+1)\cdots(t+k)}{(k+1)!}f_{m-\frac{k-1}2}^{k+1},\]
\[R_{m,k}(x)h^{k+2}\frac{t(t+1)\cdots(t+k+1)}{(k+2)!}
\frac{\d^{k+1}}{\d x^{k+2}}f(\xi,y(\xi)).\]
Dosadíme do \eqref{mnohokr},
\[y_{m+1}=y_m+\int_{x_m}^{x_{m+1}}L_{m,k}(x)\d x+
\int_{x_m}^{x_{m+1}}R_{m,k}(x)\d x.\]
Zavedeme substituci,
\[
\begin{split}
y_{m+1}&=y_m+\int_{-1}^0 L_{m,k}(x)\d t+l_{m,k}=
\\&=y_m+h\left[
f_{m+1}+b_1f_{m+\frac12}^1+b_2f_m^2+\dots+b_{k+1}
f_{m-\frac{k-1}2^{k+1}}
\right]+l_{m,k},
\end{split}
\]
kde
\[b_i=\int_{-1}^0\frac{t(t+1)\cdots(t+i-1)}{i!}\d t.\]
Uvedeme některá $b$:
\[
\begin{split}
b_1&=\int_{-1}^0 t\d t=-\frac12,\ b_2=-\frac{1}{12},\
b_3=-\frac{1}{24}, b_4=-\frac{13}{720},\\b_5&=-\frac{9}{160},\
b_6=-\frac{863}{60480}, b_7=-\frac{275}{24135},\
b_8=-\frac{33953}{3628800}.
\end{split}
\]
Konkrétní metody:
\[
\begin{split}
k&=-1:\ y_{m+1}\doteq y_m+hy'_{m+1}\\
k&=0:\ y_{m+1}\doteq y_m+\frac h2[y_{m+1}'+y_m']\\
k&=1:\ y_{m+1}\doteq y_m+\frac h{12}[5y_{m+1}'+8y_{m}'-1y_{m-1}']\\
k&=2:\ y_{m+1}\doteq y_m+\frac h{24}[9y_{m+1}'+19y_{m}'-5y_{m-1}'+y_{m-2}']\\
k&=3:\ y_{m+1}\doteq y_m+\frac h{720}
[251y_{m+1}'+646y_{m}'-264y_{m-1}'+106y_{m-2}'-19y_{m-3}']
\end{split}
\]
Chyba aproximace je
\[
\begin{split}
l_{m,k}&=\int_{x_m}^{x_{m+1}}R_{m,k}(x)\d x=
h^{k+3}\int_{-1}^0
\frac{t(t+1)\cdots(t+k+1)}{(k+2)!}\frac{\d^{k+2}}{\d x^{k+2}}
f(\xi,y(\xi))\d t=\\
&=h^{k+3}y^{(k+3)}(\eta)b_{k+2},
\end{split}
\]
celkem $\abs{l_{n,k}}\le h^{k+3}\abs{b_{k+2}}M_{k+3}$.
Implicitní metody se používají v~tzv.\index{numerická metoda, prediktor-korektor}
{\bf metodách prediktor-korektor}. Nejprve se vypočte $y_{m+1}$ pomocí explicitního
vzorce, pak se získaná hodnota dosadí do pravé strany implicitního
vzorce. Lze dokázat, že nová hodnota $y_{m+1}$ je přesnější.
\begin{define}
{\bf Obecná diferenční formule $k$-tého řádu} pro řešení rovnice
$y'=f(x,y)$ je formule tvaru
\begin{equation}
\label{odiff}
\sum_{i=0}^k\alpha_i y_{n+i}=
h\sum_{i=0}^k\beta_i f(x_{n+i},y_{n+i}),
\end{equation}
kde $n\in\{0,1,\dots,N-k\}$ a $a_k\not=0$. Pro $\beta_k=0$ je to {\bf
explicitní formule}, pro $\beta_k\not=0$ je to {\bf implicitní
formule}.
\end{define}
\begin{define}
Diferenční formule \eqref{odiff} se nazývá {\bf stupně $p\in\No$},
jestliže platí následujících $p+1$ podmínek:
\[
\begin{split}
\sum_{i=0}^k\alpha_i&=0\\
\sum_{i=0}^k\frac{i^s\alpha_i}{s!}&=
\sum_{i=0}^k\frac{i^{s-1}\beta_i}{(s-1)!}\quad
\text{pro $s=1,\dots,p$}.
\end{split}
\]
\end{define}
Známe-li hodnoty $y_n,y_{n+1},\dots,y_{n+k-1}$ přesně a vypočítáme-li
$y_{n+k}$, je to s přesností $O(h^{p+1})$.
\begin{proof}
\[
\begin{split}
&\sum_{i=0}^k[\alpha_i y(x_n+ih)-h\beta_i f(x_n+ih,y(x_n+ih))]=\\
&\quad=\sum_{i=0}^k[\alpha_iy(x_n+ih)-h\beta_iy'(x_n+ih)]
=c_{p+1}h^{p+1}y^{(p+1)}(x_n)+O(h^{p+1})
\end{split}
\]
\begin{align*}
y(x_n+ih)&=y(x_n)+ihy'(x_n)+\frac{(ih)^2}{2!}y''(x_n)+\dots+
\frac{(ih)^{p+1}}{(p+1)!}y^{(p+1)}(x_n)+\\
&\quad+\frac{(ih)^{p+2}}{(p+2)!}y^{(p+2)}(\xi_1)\\
y'(x_n+ih)&=y'(x_n)+ihy''(x_n)+\dots+
\frac{(ih)^p}{p!}y^{(p+1)}(x_n)+
\frac{(ih)^{p+1}}{(p+1)!}y^{(p+2)}(\xi_2)
\end{align*}
\end{proof}
Uvedené požadavky ještě nestačí k tomu, aby daly rozumnou metodu.
\begin{example}
Zkonstruujeme explicitní formuli 2.~řádu, co nejpřesnější:
\[y_{n+2}+\alpha_1y_{n+1}+\alpha_0y_n=
h(\beta_1y_{n+1}'+\beta_0y_n').\]
Z podmínek pro stupeň formule máme
\begin{align*}
1+\alpha_1+\alpha_0&=0\\
2+\alpha_1&=\beta_1+\beta_2\\
\frac12(4+\alpha_1)&=\beta_1\\
\frac16(8+\alpha_1)&=\frac12\beta_1.
\end{align*}
Víc podmínek si klást nelze, protože je u formule 2.~řádu nelze
splnit. Metoda má v jednom kroku přesnost $O(n^4)$. Řešení soustavy je
$\alpha_1=4$, $\alpha_0=-5$, $\beta_1=4$, $\beta_0=2$. Po dosazení
\[y_{n+2}=-4y_{n+1}+5y_n+4(4y'_{n+1}+2y_n').\]
Zkusíme řešit rovnici $y'=-y$, $y(0)=1$. Víme, že řešení je
$y(x)=e^{-x}$. Výše uvedená formule bude mít pro tuto rovnici tvar
\[y_{n+2}=-4(1+h)y_{n+1}+(5-2h)y_n.\] Zvolme krok $h=0.1$, počáteční hodnotu
$y_1=e^{-0.1}\doteq 0.904837$. Dostaneme následující
výsledky (IEEE Double, zaokrouhleno):
\begin{tabular}{d|dr|dr}
\multicolumn{1}{c}{$x$} & \multicolumn{1}{c}{$y_i$} &
\multicolumn{1}{c}{$(y_i-y(x_i))\cdot10^6$} &
\multicolumn{1}{c}{$y_i'$} &
\multicolumn{1}{c}{$(y_i'-y(x_i))\cdot10^6$} \\
\hline
0.0 & 1.000000 & 0 & 1.000000 & 0 \\
0.1 & 0.904837 & 0 & 0.904835 & -2 \\
0.2 & 0.818715 & -15 & 0.818726 & -4 \\
0.3 & 0.740872 & 54 & 0.740812 & -6 \\
0.4 & 0.669997 & -323 & 0.670313 & -7 \\
0.5 & 0.608200 & 1669 & 0.606522 & -8 \\
0.6 & 0.539907 & -8905 & 0.548803 & -9 \\
0.7 & 0.543769 & 47183 & 0.496576 & -10 \\
0.8 & 0.198971 & -250358 & 0.449319 & -10 \\
0.9 & 1.734618 & 1328048 & 0.406560 & -10 \\
1.0 & -6.677259 & -7045138 & 0.367869 & -10 \\
2.5 & & & 30.828821 & 30746736 \\
%0 & 1.0 & 0 & 1.0 & 0\\
%0.1 & 0.904837 & 0 & 0.904835 & -2\\
%0.2 & 0.818717 & -14 & & -5\\
%0.3 & 0.740863 & -276 & & -4\\
%0.4 & 0.670044 & 1418 & & -17\\
%0.5 & 0.607949 & & & 43\\
%1.0 & -5.624569 & -5992448 & & -216788
\end{tabular}
\end{example}
\begin{define}
Formule \eqref{odiff} se nazývá {\bf stabilní podle Dahlquista},
jestliže všechny kořeny polynomu
$\alpha_k\lambda^k+\alpha_{k-1}\lambda^{k-1}+\dots+\alpha_0$ jsou v absolutní
hodnotě menší nebo rovny $1$ a ty, které jsou v absolutní hodnotě
rovny 1, jsou jednoduché.
\end{define}
\begin{lemma}
\label{lemma1}
Nechť $y(x)$ je řešení rovnice $y'=f(x,y)$ v intervalu $\la
x_0,a\ra$, nechť $f$ je definována, spojitá a lipschitzovská vzhledem
k $y$ s konstantou $M$ na $\la x_0,a\ra\times\la
-\infty,+\infty\ra$. Pak pro každé celé $n$: $0\le n\le N$ označme
$u_n=y(x_n)-\lambda f(x_n,y(x_n))$, kde $\lambda$ je libovolně zvolená
hodnota taková, že $\abs{\lambda}M<1$.
Rovnice $z-\lambda f(x_n,z)=\tilde u_n$ má jediné řešení $z=\tilde y_n$
a platí
\[\abs{\tilde y_n-y_n}\le\frac{1}{1-\abs{\lambda}M}
\abs{\tilde u_n-u_n}.\]
\begin{proof}
Zavedeme zobrazení $F(z)=\lambda f(x_n,z)+\tilde u_n$ $\R\mapsto\R$,
dokážeme, že $F(z)$ je kontrahující, tedy $F(z)=z$ má právě jedno
řešení. Protože platí
\[\abs{F(z)-F(y)}=\abs{\lambda}\abs{f(x,z)-f(x,y)}\le
M\abs{\lambda}\abs{z-y}\] a $M\abs{\lambda}<1$, je $F(x)$
kontrahující a výše uvedená rovnice má jediné řešení. Dále platí
\[\abs{\tilde y_n-y_n}=\abs{\lambda f(x_n,\tilde y_n)+\tilde u_n-
\lambda f(x_n,y_n)-u_n}\le\abs{\tilde u_n-u_n}+
\abs{\lambda}M\abs{\tilde y_n-y_n},\]
z toho po úpravě
\[\abs{\tilde y_n-y_n}\le\frac{1}{1-\abs{\lambda}M}\abs{\tilde
u_n-u_n}.\]
\end{proof}
\end{lemma}
\begin{lemma}
\label{lemma2}
Nechť je dán polynom
\[\rho(t)=\sum_{i=0}^k\alpha_it^i\]
a matice
\[
\mathbf A=
\begin{pmatrix}
0 & 1 & \\
0 & 0 & 1 \\
& & & \ddots \\
& & & & 1 \\
-\frac{\alpha_0}{\alpha_k} & -\frac{\alpha_0}{\alpha_k}
& & & -\frac{\alpha_{k-1}}{\alpha_k} \\
\end{pmatrix}.
\]
Zřejmě $\rho(t)$ je násobkem charakteristického polynomu matice
$\mathbf A$. Nechť pro každý kořen $t_i$ polynomu $\rho(t)$ platí
$\abs{t_i}\le 1$ a nechť ty kořeny, pro které $\abs{t_i}=1$, jsou
jednoduché.
Zvolme nějakou maticovou normu. Pak existuje konstanta $G$, která
závisí pouze na koeficientech $\rho(t)$ tak, že platí $\norm{\mathbf
A^n}\le G$ pro $n\in\N$.
Jestliže pro všechny kořeny $\rho(t)$ platí $\abs{t_i}<1$, pak
existují konstanty $G_1$ a $0<\gamma<1$ tak, že $\norm{\mathbf A^n}\le
G_1\gamma^n$ pro $n\in\N$.
\begin{proof}
$\mathbf A$ je regulární, tudíž ji lze zapsat jako
\[\mathbf A=\mathbf T^{-1}
\begin{pmatrix}
J_1\\
& J_2\\
& & \ddots
\end{pmatrix}
\mathbf T,\]
\[
\mathbf A^n=\mathbf T^{-1}
\begin{pmatrix}
J_1^n\\
& J_2^n\\
& & \ddots
\end{pmatrix}
\mathbf T.
\]
Pro maximovou normu platí
$\mnorm{\mathbf A^n}\le
\mnorm{\mathbf T^{-1}}\tilde G
\mnorm{\mathbf T}\le
\Tilde{\Tilde G}$, neboť vlastní čísla jsou buď menší než $1$, v tom
případě jsou prvky $J_i^n$ k nule, vlastní číslo $1$ může být pouze
jednonásobné, v tom případě ale $J=(1)$, takže $J_i^n$ jsou omezeny
konstantou. Pro jiné normy to platí díky topologické ekvivalenci.
Je-li $\rho(\mathbf A)<1$, pak existuje $\epsilon>0$ takové, že
i $\rho(\mathbf A)+\epsilon<1$. Pak
\[
\mathbf A=(\rho(\mathbf A)+\epsilon)\mathbf T^{-1}
\begin{pmatrix}
\frac{1}{\rho(\mathbf A)+\epsilon} J_1\\
& \frac{1}{\rho(\mathbf A)+\epsilon} J_2\\
& & \ddots
\end{pmatrix}
\frac{1}{\rho(\mathbf A)+\epsilon} T
\]
a
\[\norm{\mathbf A^n}\le\norm{\mathbf T^{-1}}
\norm{
\begin{matrix}
\left(\frac{1}{\rho(\mathbf A)+\epsilon} J_1\right)^n\\
& \left(\frac{1}{\rho(\mathbf A)+\epsilon} J_2\right)^n\\
& & \ddots
\end{matrix}}
\norm{\mathbf T}(\rho(\mathbf A)+\epsilon)^n.
\]
\end{proof}
\end{lemma}
\begin{lemma}
\label{lemma3}
Nechť $\phi_k,\psi_k,\chi_k$ jsou konečné posloupnosti čísel pro
$k=0,1,\dots,n$, $\chi_k\ge 0$. Nechť
\[\phi_k\le\psi_k+\sum_{i=0}^{k-1}\chi_i\phi_i\quad
\text{pro $k=0,1,\dots,n$}.\]
Potom platí
\[\phi_n\le\psi_n+\sum_{i=0}^{n-1}\chi_i\psi_i
\prod_{j=i+1}^{n-1}(1+\chi_j).\]
\begin{proof}
Zavedeme další posloupnost
\[\Phi_k=\psi_k+\sum_{i=0}^{k-1}\chi_i\Phi_i,\]
$\Phi_0=\psi_0$. Indukcí dokážeme, že $\phi_k\le\Phi_k$. Pro $k=0$ to
z definice platí, pokud to platí až do $k-1$, je
\[\sum_{i=0}^{k-1}\chi_i\phi_i\le\sum_{i=0}^{k-1}\chi_i\Phi_i,\]
neboť $\chi_i\ge 0$, což bylo dokázat.
Dále dokážeme, že
\[\Phi_n=\psi_n+\sum_{i=0}^{n-1}\chi_i\psi_i
\prod_{j=i+1}^{n-1}(1+\chi_j).\]
Pro $k=0$ to opět z definice platí, předpokládejme platnost až do
$k-1$. Z definice je
\[
\begin{split}
\Phi_k&=\psi_k+\chi_0\Phi_0+\chi_1\Phi_1+\dots+\chi_{k-1}\Phi_{k-1}=\\
&=\psi_k+\chi_0\psi_0+\\
&\qquad+\chi_1[\psi_1+\chi_0\psi_0]+\\
&\qquad+\chi_2[\psi_2+\chi_0\psi_0(1+\chi_1)+\chi_1\psi_1]+\\
&\qquad+\chi_3[\psi_3+\chi_0\psi_0(1+\chi_1)(1+\chi_2)+
\chi_1\psi_1(1+\chi_2)+\chi_2\psi_2]+\\
&\qquad+\chi_4[\psi_4+\chi_0\psi_0(1+\chi_1)(1+\chi_2)(1+\chi_3)+
\chi_1\psi_1(1+\chi_2)(1+\chi_3)+\\
&\qquad\quad+\chi_2\psi_2(1+\chi_3)+\chi_3\phi_3]+\\
&\qquad+\chi_{k-1}[\psi_{k-1}+
\chi_0\psi_0(1+\chi_1)(1+\chi_2)\cdots(1+\chi_{k-1})+\\
&\qquad\quad+\chi_1\psi_1(1+\chi_2)\cdots(1+\chi_{k-1})+\cdots]=\\
&=\psi_k+\chi_0\psi_0[1+\chi_1+\chi_2(1+\chi_1)+
\chi_3(1+\chi_1)(1+\chi_2)+\\
&\qquad\quad+\chi_4(1+\chi_1)(1+\chi_2)(1+\chi_3)+\dots
+\chi_0\psi_0(1+\chi_1)\cdots(1+\chi_{k-1})]=\cdots
\end{split}
\]
Budeme vytýkat $(1+\chi_1)$, $(1+\chi_2)$ až $(1+\chi_{k-1})$.
\end{proof}
\end{lemma}
\begin{theorem}
Nechť pravá strana diferenciální rovnice $y'=f(x,y)$ je definována a
spojitá v intervalu $x_0\le x\le a$, $-\infty<y<+\infty$ a splňuje
vzhledem k $y$ Lipschitzovu podmínku s konstantou $M$. Nechť $y(x)$ je
řešení rovnice na intervalu $\la x_0,a\ra$ a nechť má spojité derivace
do řádu $p+1$, kde $p\ge 1$. Nechť uvažovaná diferenční formule řádu
$k$ je stupně $p$ a stabilní podle Dahlquista. Nechť
$y_0,y_1,\dots,y_n$ jsou hodnoty vypočítané podle vzorců
\[\sum_{i=0}^k\alpha_i y_{n+i}=
h\sum_{i=0}^k\beta_i f(x_{n+i},y_{n+i})+\delta_n
\text{ pro $n=0,1,\dots,N-k$}\]
a $y_i$ pro $i=0,1,\dots,k-1$ je dáno. Pak (pro dostatečně malá $h$)
platí
\[
\abs{y_n-y(x_n)}\le\frac{G}{1-\abs{\lambda}M}\left[
(1+\abs{\lambda}M)\vartheta+\frac{\abs{x_n-x_0}}{\abs{\alpha_n}}
\left(\frac{\delta}{h}+Kh^p\right)
\right]e^{G\tilde M(x_n-x_0)},
\]
kde $n=0,1,\dots,N$,
\[\vartheta=\max_{i=0,1,\dots,k-1}\abs{y_i-y(x_i)},\quad
\delta=\max_{i=0,1,\dots,N-k}\abs{\delta_i},\quad
\lambda=h\frac{\beta_k}{\alpha_k}\]
a $\tilde G$, $\tilde M$, $K$ jsou konstanty, které závisí jen na
koeficientech formule a na diferenciální rovnici a jsou nezávislé na
$h$ pro dostatečně malá $h$.
\begin{proof}
Označme $y_n-y(x_n)=r_n$. Odečtením vztahů
\[
\sum_{i=0}^k\alpha_iy_{n+i}=h\sum_{i=0}^k\beta_i
f(x_{n+i},y_{n+i})+\delta_n
\]
\[
\sum_{i=0}^k\alpha_iy(x_{n+i})=h\sum_{i=0}^k\beta_i
f(x_{n+i},y(x_{n+i}))-l_n,
\]
přičemž $l_n\le Kh^{p+1}$, dostaneme
\begin{equation}
\label{chd1}
\sum_{i=0}^k\alpha_ir_{n+i}=h\sum_{i=0}^k\beta_i
[f(x_{n+i},y_{n+i})-f(x_{n+i},y(x_{n+i}))]+\delta_n+l_n.
\end{equation}
Označme $u_n=y(x_n)-\lambda f(x_n,y(x_n))$, $\tilde u_n=y_n-\lambda
f(x_n,y_n)$. Podle lemmatu \ref{lemma1} se dá odhadnout
\[\abs{r_n}\le\frac{1}{1-\abs{\lambda}M}\abs{\tilde u_n-u_n}.\]
Místo $r_n$ budeme dále odhadovat $\tilde u_n-u_n$. Platí, že
\[
\begin{split}
&\sum_{i=0}^k\alpha_i(\tilde u_{n+i}-u_{n+i})=
\sum_{i=0}^k\alpha_ir_{n+i}-\lambda\sum_{i=0}^k\alpha_i
[f(x_{n+i},y_{n+i})-f(x_{n+i},y(x_{n+i}))]=\\
&\quad=h\sum_{i=0}^k(\beta_i-\frac{\beta_k}{\alpha_k}\alpha_i)
[f(x_{n+i},y_{n+i})-f(x_{n+i},y(x_{n+i}))]+\delta_n+l_n=q_n.
\end{split}
\]
Pro $i=k$ dostaneme $\beta_k-\frac{\beta_k}{\alpha_k}\alpha_k=0$,
takže stačí sčítat od 0 do $k-1$. Zavedeme vektory
\[
\vec u^{(n)}=
\begin{pmatrix}
\tilde u_n-u_n\\
\tilde u_{n+1}-u_{n+1}\\
\vdots\\
\tilde u_{n+k-1}-u_{n+k-1}
\end{pmatrix}
,\quad
\vec q^{(n)}=
\begin{pmatrix}
0\\0\\0\\\vdots\\\frac{q_n}{\alpha_k}
\end{pmatrix}
\]
a matici
\[
\mathbf A=
\begin{pmatrix}
0 & 1 \\
0 & 0 & 1 \\
& & & \ddots\\
0 & 0 & 0 & \cdots & 1 \\
-\frac{\alpha_0}{\alpha_k} & -\frac{\alpha_1}{\alpha_k} & & \cdots &
-\frac{\alpha_{k-1}}{\alpha_k}
\end{pmatrix}.
\]
Potom platí, že $\vec u^{(n+1)}=A\vec u^{(n)}+\vec q^{(n)}$. Až do
$k-1$-té složky je to jasné, pro poslední složku to vyplývá z
předchozí rovnosti.
Dále je
\[\vec u^{(n)}=\mathbf A\vec u^{(n-1)}+\vec q^{(n-1)}=
\mathbf A^2\vec u^{(n-2)}+\mathbf A\vec q^{(n-2)}+\vec q^{(n-1)},\]
\[\vec u^{(n)}=\mathbf A^n\vec u^{(0)}+\sum_{i=0}^{n-1}
\mathbf A^{n-1-i}\vec q^{(i)}.\]
Pomocí vektoru $\vec u^{(n)}$ provedeme odhad $r_{n+i}$ nezávislý na
$i$:
\[\abs{r_{n+i}}\le\frac{1}{1-\abs{\lambda}M}\abs{\tilde u_{n+i}-u_{n+i}}\le
\frac{1}{1-\abs{\lambda}M}\mnorm{\vec u^{(n)}},\]
neboť pro $i=0,\dots,k-1$ jsou to složky $\vec u^{(n)}$.
\[
\begin{split}
\mnorm{\vec q^{(n)}}&=\frac1{\abs{\alpha_k}}\abs{q_n}\le
\frac{Mh}{\abs{\alpha_k}}\sum_{i=0}^{k-1}\abs{\beta_i-\frac{\beta_k}{\alpha_k}\alpha_i}
\abs{r_{n+i}}+\frac1{\abs{\alpha_k}}\left(\delta+Kh^{p+1}\right)\le\\
&\le h\underbrace{\left[
\frac{M}{\abs{\alpha_k}}\sum_{i=0}^{k-1}\abs{\beta_i-\frac{\beta_k}{\alpha_k}\alpha_i}
\right]
\frac{1}{1-\abs{\lambda}M}}_{\tilde M}\mnorm{\vec u^{(n)}}
+\frac1{\abs{\alpha_k}}
(\delta+Kh^{p+1})
\end{split}
\]
\[
\mnorm{\vec q^{(n)}}\le\tilde M h\mnorm{\vec u^{(n)}}+
\frac{1}{\alpha_k}(\delta+Kh^{p+1})
\]
Protože podle lemmatu \ref{lemma2} mají všechny mocniny $\mathbf A$
společný odhad $G$, platí
\[
\begin{split}
\mnorm{\vec u^{(n)}}&\le\mnorm{\mathbf A^n}\mnorm{\vec u^{(0)}}+
\tilde M h\sum_{i=0}^{n-1}\mnorm{\mathbf A^{n-i-1}}
\mnorm{\vec u^{(i)}}+\\
&\quad+\frac{\delta+Kh^{p+1}}{\abs{\alpha_k}}
\sum_{i=0}^{n-1}\mnorm{\mathbf A^{n-i-1}}\le\\
&\le\underbrace{G\tilde M h}_{\chi_i}\sum_{i=0}^{n-1}
\underbrace{\mnorm{\vec u^{(i)}}}_{\phi_i}+
\underbrace{Gn\frac{\delta+Kh^{p+1}}{\abs{\alpha_k}}+
G\mnorm{\vec u^{(0)}}}_{\psi_n}.
\end{split}
\]
Podle lemmatu \ref{lemma3} z toho, že platí
\[\phi_k\le\psi_k+\sum_{i=0}^{k-1}\chi_i\phi_i,\]
vyplývá, že
\[\phi_n\le\psi_n+\sum_{i=0}^{n-1}\chi_i\psi_i\prod_{j=i+1}^{n-1}(1+\chi_j)\]
\[
\begin{split}
\mnorm{\vec u^{(n)}}&\le G\tilde M h\sum_{i=0}^{n-1}\left[
Gi\frac{\delta+Kh^{p+1}}{\abs{\alpha_k}}+G\mnorm{\vec u^{(0)}}
\right](1+G\tilde Mh)^{n-i-1}+\\
&\quad+Gn\frac{\delta+Kh^{p+1}}{\abs{\alpha_k}}+G\mnorm{\vec u^{(0)}}
\end{split}
\]
\[
\begin{split}
\mnorm{\vec u^{(n)}}&\le G\mnorm{\vec u^{(0)}}\left[
1+G\tilde M h\sum_{i=0}^{n-1}(1+G\tilde Mh)^{n-i-1}
\right]+\\
&\quad+\frac{\delta+Kh^{p+1}}{\abs{\alpha_k}}\left[
Gn+G\tilde M h\sum_{i=0}^{n-1}Gi(1+G\tilde M h)^{n-i-1}
\right]
\end{split}
\]
Protože
\[
\sum_{i=0}^{n-1}iq^{n-1-i}=\sum_{j=1}^{n-1}\sum_{i=0}^{j-1}q^i,
\]
je druhá hranatá závorka rovna
\[
\begin{split}
&Gn+G\tilde Mh\sum_{j=1}^{n-1}G\frac{(1+G\tilde Mh)^j -1}
{G\tilde Mh}=Gn+G\sum_{j=1}^{n-1}(1+G\tilde Mh)^j-G(n-1)=\\
&\quad=G+G(1+G\tilde Mh)\frac{(1+G\tilde Mh)^{n-1}-1}{G\tilde Mh}=
\frac{1}{\tilde Mh}[(1+G\tilde Mh)^n-1].
\end{split}
\]
\[
\mnorm{\vec u^{(n)}}\le G\mnorm{\vec u^{(0)}}[1+(1+G\tilde Mh)^n-1]+
\frac{\delta+Kh^{p+1}}{\abs{\alpha_k}}
\frac{1}{\tilde Mh}[(1+G\tilde Mh)^n-1],
\]
to lze dále upravit:
\[
\mnorm{\vec u^{(n)}}\le G\mnorm{\vec u^{(0)}}(1+G\tilde Mh)^n+
\frac{1}{\abs{\alpha_k}\tilde M}\left(
\frac\delta h+Kh^p
\right)[(1+G\tilde Mh)^n-1]
\]
Využijeme nyní toho, že $1+G\tilde Mh<e^{G\tilde Mh}$, tedy
$(1+G\tilde Mh)^n<e^{nG\tilde Mh}=e^{G\tilde M(x_n-x_0)}$ a
dále $e^{G\tilde M(x_n-x_0)}-1<G\tilde M(x_n-x_0)e^{G\tilde
M(x_n-x_0)}$. To první je jasné, to druhé vyplývá třeba z Taylorova
rozvoje $e^x$.
\[\mnorm{\vec u^{(n)}}\le G\mnorm{\vec u^{(0)}}e^{G\tilde M(x_n-x_0)}+
\frac{G(x_n-x_0)}{\abs{\alpha_k}}\left(
\frac{\delta}{h}+Kh^p
\right)e^{G\tilde M(x_n-x_0)}\]
Po vynásobení $1/(1-\abs{\lambda}M)$, dostáváme hledanou
nerovnost. Stačí si uvědomit, že
\[
\begin{split}
\abs{\tilde u_i-u_i}&=
\abs{y_i-y(x_i)-\lambda[f(x_i,y_i)-f(x_i,y(x_i))]}\le\\
&=\abs{y_i-y(x_i)}+\abs{\lambda}M\abs{y_i-y(x_i)}=
(1+\abs{\lambda}M)\abs{y_i-y(x_i)}=\\
&=(1+\abs{\lambda}M)\vartheta.
\end{split}
\]
\end{proof}
\end{theorem}
\begin{remark}
$\delta_n$ je zaokrouhlovací chyba počítače a ani počáteční podmínky
nejsou dány přesně. Opět nelze jít s $h$ libovolně nízko, pro další
zvýšení přesnosti lze použít jedině metodu vyššího řádu.
\end{remark}
