Príklady z hodiny

Na prvých dvoch prednáškach ste preberali základy numerických metód a vyjasnili ste si pojmy ako reprezentácia čísel v počítači, typy chýb pri riešení numerických úloh, korektnosť, podmienenosť a numerická stabilita. Na tomto cvičení si na príkladoch demonštrujeme niektoré z týchto poznatkov.

Na začiatok si pripomeňme, že zdroje chýb pri numerických metódach sa všeobecne delia na

  • chyby vstupných dát (napr. chyby nameraných dát, ktoré potom zadávame do kódu, alebo chyba modelu reality),
  • chyby metódy (truncation errors),
  • zaokrúhľovacie chyby (roundoff errors).

Rád metódy
Odvodenie rádu metódy pre metódy aproximácie derivácie ste si ukázali na prednáške. V rámci opakovania si ukážeme, ako odvodiť rád metódy doprednej aproximácie prvej derivácie funkcie $f$:

$f’(x)\sim\frac{f(x+h)-f(x)}{h}.$

Rád metódy je rovný najmenšej mocnine konečneho krátkeho kroku $h$ (ktorý nahrádza nekonečne krátky krok), ktorá ostala po aproximácii analytickej derivácie. Urobíme Taylorov rozvoj členu $f(x+h)$:

$f’(x)\sim\frac{f(x+h)-f(x)}{h}\sim\frac{f(x)+h\frac{f’(x)}{1!}+h^2\frac{f\prime\prime(x)}{2!}+h^3\frac{f\prime\prime\prime(x)}{3!}-f(x)}{h}$,

$f’(x)\sim\frac{f(x+h)-f(x)}{h}\sim f’(x)+h\frac{f\prime\prime(x)}{2!}+h^2\frac{f\prime\prime\prime(x)}{3!}$.

V tomto prípade vidíme, že metóda je prvého rádu, pretože člen s najnižšou mocninou, ktorý nám po Taylorovom rozvoji ostal, je úmerný $h^1$. Je potrebné si uvedomiť, že krok $h$ je veľmi malý, preto členy s vyššími mocninami $h$ prispievajú pri aproximácii menšou chybou.
TIP: Pri rozvoji si dajte pozor, do koľkých členov daný výraz rozviniete. Ak sa Vám všetky výrazy s $h$ pri aritmetických úpravach vyrušia, je potrebné sa vrátiť na začiatok, a rozvinúť Taylora do viacerých členov.

Zaokrúhľovacie chyby
Na nasledujúcom príklade si ukážeme, ako sa zvýši relatívna chyba výpočtu, ak od seba odčítame podobne veľké čísla. Majme čísla

$x = 1.32483726$,

$y = 1.32483357$.

Predstavme si, že to mohli byť pôvodne čísla s vyšším počtom desatinných miest, my sme však mali miesto len na to, aby sme ich ukladali presnosťou na 9 platných číslic. Vezmime si napr. číslo $x$. To sme mohli získať zaokrúhlením čísel $1.324837255$, $1.324837256$, $1.324837257$, $1.324837258$, $1.324837259$, $1.324837260$,$1.324837261$,$1.324837262$, $1.324837263$, $1.324837264$ (a ďalšími číslami, ktoré by mali väčší počet desatinných miest než 9). Najvǎčší rozdiel by vznikol pri zaokrúhľovaní čísla $1.324837255$. (Keďže robíme len odhad, zanedbáme zaokrúhlenie čísel s väčším počtom desatinných miest, ktoré by bolo pri našom odhade chyby zanedbateľné.) Absolútnu chybu odhadneme tým najväčśim rozdielom, teda

$a(x) = \lvert1.32483726-1.324837255\rvert=5\times10^{-9}$.

Odhad relatívnej chyby potom môžeme vypočítať ako

$r(x) = \frac{a(x)}{\lvert x \rvert}=\frac{\lvert 1.324837255-1.32483726 \rvert}{\lvert 1.32483726 \rvert}=4\times10^{-9}$.

To isté platí pre $y$, t.j. $a(y) = 5\times10^{-9}$, $r(y) = 4 \times10^{-9}$. Vykonajme teraz rozdiel

$\lvert x - y \rvert = \lvert 1.32483726 - 1.32483357 \rvert = 4 \times 10^{-6}$.

Odhad relatívnej chyby tohto rozdielu sa dá vypočítať ako

$r(\lvert x - y \rvert)=\frac{a(x)+a(y)}{\lvert x - y\rvert}=\frac{5\times10^{-9}+5\times10^{-9}}{\lvert1.32483726 - 1.32483357\rvert}=3\times 10^{-3} $.

Všimnime si, že v rámci jedinej operácie sme prešli z relatívnej presnosti $\sim 10^{-9}$ na oveľa nižšiu relatívnu presnosť $\sim 10^{-3}$, posunuli sme sa až o 6 rádov!

Numerická (ne)stabilita

Príklad stabilnej a nestabilnej metódy si ukážeme na príklade padajúceho meteoritu v Matlabe. Úlohou je riešiť diferenciálnu rovnicu

$\frac{\mathrm{d}v}{\mathrm{d}t}=-v$,

s počiatočnou podmienkou $v(0)=1$. Táto rovnica (zjednodušene) opisuje brzdnú dráhu meteoritu v atmosfére, určuje teda zmenu rýchlosti $v$ v čase $t$. Analytické riešenie tejto rovnice poznáme:

$v(t) =\mathrm{exp}(-t)$.

My ju budeme riešiť pomocou dvoch rozdielnych numerických metód, a potom budeme mať možnosť si porovnať, ktorá z nich je numericky stabilnejšia.

Eulerova metóda
Nekonečne krátky časový krok $\mathrm{d}t$ nahraďme krátkym konečným krokom $h$:

$\frac{v(t+h)-v(t)}{h} = -v(t)$

Rýchlosť v čase $t+h$ budeme teda počítať ako

$v(t+h) = -v(t)h+v(t)$

Dvojkroková metóda

Pri dvojkrokovej metóde sa derivácia aproximuje pomocou dvojnásobného kroku $2h$:

$\frac{v(t+h)-v(t-h)}{2h} = -v(t)$

Rýchlosť v čase $t+h$ budeme teda počítať ako

$v(t+h)= -2h v(t)+v(t-h)$

Všimnite si, že pri výpočte kroku v bode $v(t+h)$ sme využili hodnoty $v(t)$ aj $v(t-h)$, narozdiel od Eulerovej metódy, kde sa využívalo len $v(t)$. Prečítajte si ->analytické vysvetlenie<- stability týchto algoritmov. Vidíte, že pri dvojkrokovej metóde chyba rastie s počtom krokov - to je numerická definícia nestability. Vo videu si teraz názorne ukážeme rast tejto chyby:


Na záver si prečítajte aj ->stručný prehľad o numerickej stabilite a výbere numerickej metódy<-.