Obyčajné differenciálne rovnice/Kontrola kroku pre jednokrokové metódy

Kontrola kroku pre jednokrokové metódy

Nastavenie veľkosti kroku $h$ počas numerického výpočtu riešenia zohráva dôležitú úlohu v úlohách s počiatočnou hodnotou, ktorých pravá strana je strmá funkcia. Ak sa v takýchto prípadoch zachová pevná veľkosť kroku, numerická chyba v oblastiach strmého nárastu z$f$ stark wachsen und das numerische Verfahren versagen. Kontrola veľkosti kroku je nevyhnutná aj pre takzvané tuhé problémy s veľmi rozdielnymi rýchlosťami rastu zložiek riešenia. Problémy s tuhosťou sa budú podrobnejšie skúmať neskôr. V tejto časti opisuje, ako použiť presný výraz pre lokálnu chybu metódy na nastavenie veľkosti kroku $h$. Obmedzíme sa tu na explicitnú Euerovu metódu.

Lokálna chyba opisuje chybu, ktorú ESV pripúšťa v každom jednotlivom kroku výpočtu, pozri definíciu 3.8, alebo ako dobre sa funkcia metódy aproximuje pravej strane AWA, $f(t,y)$. V prípade explicitnej Eulerovej metódy je lokálna chyba, vyhodnotená v bode $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)$

Čím je funkcia $f$ strmšia (čím väčšia $f^{\prime }(t_i,y(t_i))=y^{″}(t_i)$), tým väčšia je lokálna chyba $\tau (t_i,h)$. Cieľom riadenia veľkosti kroku je prispôsobiť veľkosť kroku priebehu funkcie $f$ tak, aby sa lokálna chyba v každom kroku výpočtu $i=1,\dots N_h-1$ udržala pod určitou toleranciou.

Aby sme našli vhodnú veľkosť kroku v každom kroku, musíme vopred "skenovať" numerické riešenie a budúcu chybu. Pri explicitnej Eulerovej metóde sa tento odhad dá urobiť porovnaním numerického riešenia s riešením pomocou dvoch "polovičných" Eulerových krokov (s veľkosťou kroku $\frac{h}{2}$). Aby sme to vysvetlili, najprv preskúmame lokálnu chybu Eulerovej metódy s polovičnou veľkosťou kroku: $\begin{array}{ccc}\text{1. krok:}& & {\bar{\eta }}_{i+\frac{1}{2}}={\bar{\eta }}_i+\frac{h}{2}f(t_i,{\bar{\eta }}_i)),\\ \text{2. Krok: }& & {\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)$

Lokálna chyba tejto metódy je dvakrát menšia ako chyba Eulerovho vzorca, pretože $\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))+\frac{h}{2}f(t_i+\frac{h}{2},y(t_i)+\frac{h}{2}f(t_i,y(t_i)))]\\ \\ stackrelTaylor=& & 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)$

porovnaj (3.16). Podľa definície lokálnej chyby je lokálna chyba rovná rozdielu medzi presným a numerickým riešením (vypočítaným z presnej počiatočnej hodnoty), $\epsilon (t_i,h)=y(t_i+h)-{\eta }_{i+1}$, pozri poznámku 3.1 bod ii). V súlade s tým pre Eulerovu metódu dostaneme $${\eta }_{i+1}\stackrel{(3.16)}{=}y(t_i+h)-\frac{h^2}{2}{y}^{″}(t_i)+𝒪(h^3),$$

a pre Eulerovu metódu s polovičnou veľkosťou kroku $${\bar{\eta }}_{i+1}\stackrel{(3.18)}{=}y(t_i+h)-\frac{h^2}{4}{y}^{″}(t_i)+𝒪(h^3).$$

Toto zistenie vedie k dvom záverom:
a) Lokálnu chybu Eulerovej metódy možno predpovedať pomocou rozdielu numerického riešenia s veľkosťou kroku $h$ a polovičnou veľkosťou kroku, $$2({\bar{\eta }}_{i+1}-{\eta }_{i+1})=\frac{h^2}{2}{y}^{″}(t_i)+\widetilde{𝒪}(h^3).(3.19)$$

b) S vhodnou kombináciou dvoch riešení ${\bar{\eta }}_{i+1},{\eta }_{i+1}$ je možné odstrániť vedúci chybový prvok $\frac{h^2}{2}{y}^{″}(t_i)$ a rekonštruovať numerické riešenie druhého rádu, pretože $$2{\bar{\eta }}_{i+1}-{\eta }_{i+1}=y(t_i+h)+𝒪(h^3).$$

Rovnica (3.19) je základom pre úpravu veľkosti kroku Eulerovej metódy. Lokálnu chybu $\frac{h^2}{2}{y}^{″}(t_i)+𝒪(h^3)$ možno vopred určiť výpočtom rozdielu medzi numerickými riešeniami s plnou a polovičnou veľkosťou kroku $2({\bar{\eta }}_{i+1}-{\eta }_{i+1})\equiv \epsilon (t_i,h)$. Veľkosť kroku $h$ sa volí tak, aby lokálna chyba neprekročila predtým zvolenú toleranciu $tol$:

ALGORITM: Regulácia šírky kroku

1. $t=t_i$, vyberte veľkosť skúšobného kroku $H$ a toleranciu $tol$ pre lokálnu chybu.
   

V ďalšom kroku vypočítajte numerické riešenie $t_i+H:$.

1. • pomocou jedného kroku Eulerovej metódy s veľkosťou kroku $H$: ${\eta }_{i+1}^H$,
   • pomocou dvoch krokov Eulerovej metódy s polovičnou veľkosťou kroku $\frac{H}{2}:$, ${\bar{\eta }}_{i+1}^H$, pozri (3.17).
2. Odhadnite... vedúca konštanta chyby pomocou (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. Zvoľte optimálnu veľkosť kroku $h_{opt}$ tak, aby lokálna chyba $\epsilon (t_i,h_{opt})$ neprekročila zadanú toleranciu $tol$: $$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.$$ Rovnica (3.20) Parameter $\alpha$ sa tu môže zvoliť na dodatočné škálovanie optimálnej veľkosti kroku. V praxi sa často používa $\alpha =0,9$.
4. Určite numerické riešenie ${\eta }_{i+1}$ v ďalšom kroku $t_i+h_{opt}$ pomocou veľkosti kroku $h_{opt}$ a metódy druhého rádu, ${\eta }_{i+1}:={\bar{\eta }}_{i+1}^{h_{opt}}$, siehe (3. 17).
   

Tento krok nie je povinný; možno tu použiť aj klasickú explicitnú Eulerovu metódu, ale so stratou presnosti (poradie konzistencie).

Namiesto tolerancie absolútnej lokálnej chyby sa pri riadení veľkosti kroku často predpokladá tolerancia relatívnej chyby (%). V tomto prípade sa namiesto $\epsilon (t_i,H)$ použije relatívna chyba ${\epsilon }_{rel}(t_i,H)=\frac{\epsilon (t_i,H)}{{\eta }_{i+1}}$ sa používa na výpočet kroku, $$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.$$