1. Úvod do numerické matematiky Tento učební text byl podpořen z Operačního programu Praha - Adaptabilita
Irena Sýkorová
Zdroje chyb Reálné problémy, např. fyzikální, ekonomické, biologické, vyjadřujeme matematickým modelem. K řešení matematické úlohy pak užijeme nějakou numerickou metodu. Každá numerická úloha má konečný počet vstupních i výstupních dat a při jejím řešení můžeme provést jen konečný počet operací. Obsahem numerické matematiky jsou především metody, kterými se matematická úloha převádí na numerickou. Při procesu řešení nějakého problému se vždy dopouštíme určitého zjednodušení, z čehož vznikají různé chyby. a) Matematický model je zjednodušením reality, nepřesnosti se říká chyba matematického modelu. b) Místo přesného řešení matematického modelu získáme často jen jeho aproximaci, to je tzv. chyba numerické metody (nekonečný proces nahradíme konečným). c) Výsledek ovlivňují i chyby ve vstupních datech, což jsou například chyby měření. d) Je třeba vzít v úvahu i zaokrouhlovací chyby; ty jsou způsobené uložením čísel v počítači a další zaokrouhlovací chyby vznikají při provádění početních operací. Pokud vybíráme vhodnou numerickou úlohu k řešení nějakého problému, musíme vzít v úvahu otázku konvergence, složitost algoritmu (počet operací), nároky na paměť počítače a stabilitu. Algoritmus je stabilní, jestliže je dobře podmíněný, tj. málo citlivý na chyby ve vstupních datech a numericky stabilní, tj. málo citlivý na zaokrouhlovací chyby vznikající v průběhu výpočtu.
Zobrazení čísel v počítači Jsme zvyklí počítat s čísly vyjádřenými v desítkové poziční soustavě. V počítači jsou však čísla uložena ve dvojkové soustavě (pouze 0 a 1) a pro každé číslo je vyhrazen jen konečný počet míst. Z toho mohou pramenit různé nepřesnosti výpočtu: (1) Převedením čísla z desítkové do dvojkové soustavy se může stát, že číslo s konečným rozvojem je nahrazeno číslem s nekonečným rozvojem (např. číslo 0, 4 je ve dvojkové soustavě periodické, (0, 4)10 = (0, 0110)2 ). (2) Číslo s nekonečným rozvojem je zobrazeno jen s konečným počtem číslic. (3) Při provádění operací mohou výsledek ovlivnit zaokrouhlovací chyby. (4) Pro počítačové operace neplatí asociativní a distributivní zákon. Nemusí platit (x + y) + z = x + (y + z),
(x · y) · z = x · (y · z),
x·
1 = 1. x
Neplatnost asociativního zákona ukážeme na jednoduchém příkladu. Příklad 1: Předpokládejme, že můžeme zobrazit maximálně dvojciferná čísla, a máme za úkol vypočítat 82 + 46 − 31. Budeme-li počítat zleva, (82 + 46) − 31 = 128 − 31 nastane chyba, protože číslo 128 už nemůžeme zobrazit. Zatímco, při druhém způsobu uzávorkování, tj. 82 + (46 − 31) = 82 + 15 = 97 výpočet proběhne bez problémů. Každé nenulové reálné číslo může být jednoznačně vyjádřeno v tzv. normalizovaném tvaru ! ∞ X −i E ai 2 x = sgn(x) · 2 · 1 + , kde E ∈ Z , ai ∈ {0, 1}. i=1
2
Výraz M =1+
∞ X
ai 2−i
i=1
se nazývá mantissa a celému číslu E se říká exponent, sgn(x) označuje znaménko čísla x. Reálná čísla jsou v počítači zobrazena v pohyblivé řádové čárce, tj. každé nenulové číslo je zobrazeno ve tvaru E
x = sgn(x) · 2 ·
1+
k X
−i
ai 2
i=1
!
kde E ∈ Z, Emin ≤ E ≤ Emax , k ∈ N, ai ∈ {0, 1}.
,
Podle standardu IEEE 754 (Institute of Electrical and Electronic Engineers) mohou být reálná čísla jsou zobrazována: • v jednoduché přesnosti (32 bitů), • v dvojnásobné přesnosti (64 bitů), • v rozšířeném formátu (80 bitů).
Pro jednoduchou přesnost se používá slovo o délce 32 bitů, kde levý bit je určen pro znaménko (0 pro kladná čísla, 1 pro záporná), následujících 8 bitů je pro exponent a zbývajících 23 bitů slouží k uložení mantisy. ±
E
a1 , a2 , · · · , a23
Protože pro mantisu platí 1 ≤ M < 2, tedy nejvyšší bit je vždy roven 1, není potřeba tento bit zobrazovat. V jednoduché přesnosti je možné zobrazit exponenty od Emin = −126 do Emax = 127 (jejich počítačová reprezentace obsahuje vždy jedničky i nuly), tedy lze zobrazit normalizovaná čísla x, pro která platí Nmin ≤ |x| ≤ Nmax
-Nmax
-Nmin 0 Nmin
Nmax
Nejmenší kladné zobrazitelné normalizované číslo je 1, 0 · 2−126 ≈ 1, 2 · 10−38 . Množina zobrazitelných normalizovaných čísel je rozložena nerovnoměrně, kolem nuly je mezera, proto byla rozšířena i o tzv. subnormalní čísla (nelze je vyjádřit v normalizovaném tvaru, mantisa začíná nulou), to jsou čísla s nulovým exponentem. Nejmenší kladné subnormální číslo je 1, 4 · 10−45 . Ve dvojnásobné přesnosti jsou čísla zobrazena podobným způsobem, ale v 64 bitovém slově, kde je jeden bit určen pro znaménko, 11 bitů pro exponent a zbývajících 52 bitů pro mantisu. Exponenty 00000000 a 11111111 jsou určeny pro speciální hodnoty, kromě subnormálních čísel slouží také k označení nuly (tu je možné zobrazit dvojím způsobem jako +0 nebo −0) a čísla přesahující rozsah zobrazitelých čísel, tj. |x| > Nmax (tzv. přetečení). K označení nedefinovaných operací slouží tzv. N aN (not a number), nezobrazitelné, nedefinované číslo. Znaménkový bit
exponent E
mantisa
význam
0 nebo 1
000 0000 0
000 0000 0000 0000 0000 0000
0 resp. 1
111 1111 1
000 0000 0000 0000 0000 0000
±0 (dva způsoby zobrazení)
0 resp. 1
111 1111 1
nenulový řetězec
0 nebo 1
000 0000 0
nenulový řetězec
±∞ (dvě různé hodnoty) N aN – chyba
subnormální čísla |x| ≤ Nmin
Rozsah kladných zobrazitelných normalizovaných čísel a počet platných dekadických číslic pro jednoduchou a dvojnásobnou přesnost je uveden v následující tabulce.
3
Přesnost jednoduchá dvojnásobná
od
do −38
1, 2 · 10
−308
2, 2 · 10
počet platných dekadických číslic 38
3, 4 · 10
7
308
16
1, 8 · 10
V jednoduché přesnosti můžeme tedy počítat zhruba se sedmi platnými dekadickými číslicemi. Abychom omezili vliv zaokrouhovacích chyb, během výpočtů nezaokrouhlujeme mezivýsledky, počítáme s maximální možnou přesností.
Absolutní a relativní chyba Označíme-li x∗ hledanou přesnou veličinu a x její aproximaci, pak absolutní chyba ea je dána absolutní hodnotou jejich rozdílu ea = |x∗ − x|. Protože však přesnou hodnotu x∗ většinou neznáme, hledáme odhad absolutní chyby ε(x) přibližné hodnoty x |x∗ − x| ≤ ε(x). Přesnost lépe charakterizuje relativní chyba er , poměr absolutní chyby k absolutní hodnotě přesné veličiny |x∗ − x| . er = |x∗ | 0dhad relativní chyby δ(x) přibližné hodnoty x je |x∗ − x| ≤ δ(x). |x∗ |
Norma vektoru Aritmetický vektor x je uspořádaná n-tice reálných čísel, tj. x = (x1 , . . . , xn ). Norma aritmetického vektoru x ∈ Rn je taková nezáporná funkce (značíme k · k), která splňuje tyto podmínky: (i) kxk > 0 pro každý nenulový vektor a kxk = 0 pouze pro x = o,
(ii) kc · xk = |c| · kxk pro každé reálné číslo c,
(iii) kx + yk ≤ kxk + kyk pro každé dva vektory z Rn . Nejčastěji užívané normy vektoru jsou: a) kxk∞ = max{|x1 |, |x2 |, . . . , |xn |}, tzv. maximová norma, n P |xi |, tzv. oktaedrická norma, b) kxk1 = i=1 s n P c) kxk2 = x2i , tzv. euklidovská norma. i=1
Příklad 2: Vypočítáme všechny normy vektoru x = (1, −2, 0, 3).
kxk∞ = max{| 1 |, | − 2|, | 0 |, | 3 |} = max{1, 2, 0, 3} = 3
kxk1 = | 1 | + | − 2| + | 0 | + | 3 | = 1 + 2 + 0 + 3 = 6 p √ √ kxk2 = 12 + (−2)2 + 02 + 32 = 1 + 4 + 0 + 9 = 14
4
MATLAB – norma vektoru v: maximová: >> norm(v,inf) oktaedrická: >> norm(v,1) euklidovská: >> norm(v,2) >> norm(v) Pojem norma můžeme zavést i v jiných prostorech než jen v Rn : i) V prostoru reálných čísel R je norma totéž co absolutní hodnota kxk = |x|,
x ∈ R,
ii) v prostoru spojitých funkcí na uzavřeném intervalu Cha,bi může být norma definovaná takto: kf k = max |f (x)|, a≤x≤b
f ∈ Cha,bi .
Prostory, ve kterých je definovaná norma, se nazývají normované. P o z n á m k a : Pro u, v z normovaného prostoru M se číslo ku − vk nazývá vzdálenost prvků u a v.
Konvergence Uvažujeme nyní nějaký normovaný prostor M a posloupnost prvků vn ∈ M, n = 1, 2, . . . . Řekneme, že posloupnost prvků (vn ) konverguje k v ∈ M, tj. vn −→ v pro n −→ ∞ (resp. lim vn = v), jestliže
n→∞
kvn − vk −→ 0
pro n −→ ∞.
Řekneme, že posloupnost prvků (vn ) je Cauchyovská jestliže kvn − vm k −→ 0
pro n, m −→ ∞.
P o z n á m k a : Jestliže je posloupnost konvergentní v prostoru M, pak je i Cauchyovská v M. Normovaný prostor M se nazývá úplný, jestliže každá Cauchyovská posloupnost je konvergentní v M. Následující normované prostory jsou úplné: 1) M = R s normou kxk = |x|, 2) M = Rn s normami kxk∞ , kxk1 , kxk2 , 3) M = Cha,bi s normou kf k = max |f (x)| (prostor funkcí spojitých na uzavřeném intervalu a≤x≤b
ha, bi ).
Příklad 3: Uvažujeme M = (0, 1) s normou kxk = |x|. Pak M není úplný. Posloupnost xn = n1 je Cauchyovská v M = (0, 1), protože pro každé n ∈ N je n1 ∈ (0, 1) a 1 1 |n − m | → 0 pro n, m → ∞. Tato posloupnost však není konvergentní v M = (0, 1), protože její / M). limita 0 není v M obsažená, tj. lim n1 = 0 (∈ n→∞
5
Banachova věta Některé metody, které budeme používat, jsou založené na následující větě. Věta (Banachova). Nechť M je úplný normovaný prostor, číslo α ∈ (0, 1), f je takové zobrazení M do M, pro které platí kf (x) − f (y)k ≤ α kx − yk
pro každé
x, y ∈ M.
Pak a) Existuje jediný bod xp , pro který f (xp ) = xp . b) Pro libovolný bod x0 ∈ M posloupnost bodů xn = f (xn−1 )
pro
n = 1, 2, . . .
konverguje k xp , tj. xp = lim xn . n→∞
c) Platí odhady: αk kx0 − x1 k 1−α α kxk−1 − xk k kxk − xp k ≤ 1−α kxk − xp k ≤
P o z n á m k a : Bod xp se nazývá pevný bod zobrazení f .
6
pro pro
k > 1, k > 1.
