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

Kontrola kroku pre jednokrokové metódy

Nastavenie veľkosti kroku ${\textstyle 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${\textstyle 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 ${\textstyle 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, ${\textstyle f(t,y)}$. V prípade explicitnej Eulerovej metódy je lokálna chyba, vyhodnotená v bode ${\textstyle t_{i}}$ ${\begin{aligned}\tau (t_{i},h)&=&{\frac {y}{'}}'(\xi ),\quad \xi \in (t_{i},t_{i}+h),\quad {\mbox{bzw.,}}\\\tau (t_{i},h)&=&{\frac {y}{'}}'(t_{i})+{\cal {O}}(h^{2}).\quad \end{aligned}}\qquad (3.16)$

Čím je funkcia ${\textstyle f}$ strmšia (čím väčšia ${\textstyle f'(t_{i},y(t_{i}))=y''(t_{i})}$), tým väčšia je lokálna chyba ${\textstyle \tau (t_{i},h)}$. Cieľom riadenia veľkosti kroku je prispôsobiť veľkosť kroku priebehu funkcie ${\textstyle f}$ tak, aby sa lokálna chyba v každom kroku výpočtu ${\textstyle i=1,\ldots 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{aligned}{\mbox{1. krok:}}&&{\overline {\eta }}_{i+{\frac {1}{2}}}={\overline {\eta }}_{i}+{\frac {h}{2}}f(t_{i},{\overline {\eta }}_{i})),\\{\mbox{2. Krok: }}&&{\overline {\eta }}_{i}={\overline {\eta }}_{i+{\frac {1}{2}}}+{\frac {h}{2}}f(t_{i}+{\frac {h}{2}},{\overline {\eta }}_{i+{\frac {1}{2}}})={\overline {\eta }}_{i}+{\frac {h}{2}}f(t_{i},{\overline {\eta }}_{i})+{\frac {h}{2}}f{\big (}t_{i}+{\frac {h}{2}},{\overline {\eta }}_{i+{\frac {1}{2}}}{\big )}.\quad \end{aligned}}\qquad (3.17)$

Lokálna chyba tejto metódy je dvakrát menšia ako chyba Eulerovho vzorca, pretože ${\begin{aligned}h{\overline {\tau }}(t_{i},h)&=&{\overline {\varepsilon }}(t_{i},h)=y(t_{i}+h)-{\big [}y(t_{i})+{\frac {h}{2}}f(t_{i}))+{\frac {h}{2}}f{\Big (}t_{i}+{\frac {h}{2}},y(t_{i})+{\frac {h}{2}}f(t_{i},y(t_{i})){\Big )}{\big ]}\\\\stackrel{Taylor}{=}&&y(t_{i})+hy'(t_{i})+{\frac {h^{2}}{2}}y''(t_{i})\\&&-{\big [}y(t_{i})+{\frac {h}{2}}f(t_{i},y(t_{i}))+{\frac {h}{2}}{\big (}f+{\frac {h}{2}}(f_{t}+ff_{y}){\big )}(t_{i},y(t_{i}))){\big ]}+{\cal {O}}(h^{3})\\&=&{\frac {h^{2}}{2}}y''(t_{i})-{\frac {h^{2}}{4}}\underbrace {(f_{t}+ff_{y})(t_{i},y(t_{i}))} _{=y''(t_{i})}+{\cal {O}}(h^{3})={\frac {h^{2}}{4}}y''(t_{i})+{\cal {O}}(h^{3}),\end{aligned}}\qquad (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), ${\textstyle \varepsilon (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})+{\cal {O}}(h^{3}),$$

a pre Eulerovu metódu s polovičnou veľkosťou kroku $${\overline {\eta }}_{i+1}{\stackrel {({3.18})}{=}}y(t_{i}+h)-{\frac {h^{2}}{4}}y''(t_{i})+{\cal {O}}(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 ${\textstyle h}$ a polovičnou veľkosťou kroku, $$2({\overline {\eta }}_{i+1}-\eta _{i+1})={\frac {h^{2}}{2}}y''(t_{i})+{\tilde {\cal {O}}}(h^{3}).\qquad (3.19)$$

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

Rovnica (3.19) je základom pre úpravu veľkosti kroku Eulerovej metódy. Lokálnu chybu ${\textstyle {\frac {h^{2}}{2}}y''(t_{i})+{\cal {O}}(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 ${\textstyle 2({\overline {\eta }}_{i+1}-\eta _{i+1})\equiv \varepsilon (t_{i},h)}$. Veľkosť kroku ${\textstyle h}$ sa volí tak, aby lokálna chyba neprekročila predtým zvolenú toleranciu ${\textstyle tol}$:

ALGORITM: Regulácia šírky kroku

1. ${\textstyle t=t_{i}}$, vyberte veľkosť skúšobného kroku ${\textstyle H}$ a toleranciu ${\textstyle tol}$ pre lokálnu chybu.
   

V ďalšom kroku vypočítajte numerické riešenie ${\textstyle t_{i}+H:}$.

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