PA081: Programování numerických výpočtů 3. Systémy nelineárních rovnic Aleš Křenek jaro 2019 Řešený problém ► formulace problému Fi(xi,...,xN) = 0 pro i = l,...iV PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen 2/25 Řešený problém ► formulace problému Fi(xi,.. .,Xn) = 0 pro i = 1,...N ► není to jednoduché, neexistuje univerzální řešení PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen Řešený problém ► formulace problému Fi(xi,.m.,xN) = 0 pro i = l,...iV ► není to jednoduché, neexistuje univerzální řešení ► pro dvě dimenze, funkce f(x,y),g(x,y) PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen g pos 2/25 Neřešitelný problém ► funkce nemají obecně nic společného ► řešení celého systému se objevují v podstatě náhodně ► v obecném případě hledáme společné body hyperpovrchů dimenze N - 1 ► oblasti, kde je konkrétní Fi nulová ► konkrétní Fi jich může vygenerovat několik ► ani nevíme kolik jich je ► bez podrobnější znalosti konkrétního problému se neobejdeme PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen Metoda prosté iterace ► zobecnění metody pro jednu proměnnou ► problém přeformulujeme do podoby xi = Gi(x) x2 = G2(x) PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen xn = GN(x) ► potřebujeme nějakou metriku na RN ► je-li G(x) kontrakce na nějaké uzavřené omezené podmnožině RN, pak v ní existuje pevný bod G(§) = ► posloupnost Xi+i = G(xí) konverguje k § Metoda prosté iterace ► praktické kritérium konvergence dGi(x) dx N ► resp. slabší n 1 J = l dGjjx) dx pro všechna í,j = l,...,iVaO(F,Xi,n,u) vypočet xí + ôx Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen ei < udlxíll + 116x11) 10/25 Newtonova metoda pro systémy rovnic Odhady přesnosti - Tisseurova věta ► J(x*) je dobře podmíněná uk(J(X*)) < 1/8 ► chyba výpočtu Jakobiánu a lineárních rovnic je malá ^IIJ(Xi)"1|l0(F,Xi,n,u)|| < 1/8 ► J je lipschitzovsky spojitá funkce PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen IJ(v) - J(w)|| < j8||v- w ► počáteční odhad je dostatečně blízko j6||J(x*) -i x0 -x*|| < 1/8 Newtonova metoda pro systémy rovnic Odhady přesnosti - Tisseurova věta PA081: Programování numerických výpočtů A. Křenek ► pak relativní chyba řešení Newtonovou metodou klesá až dosáhne IIJ(x#) || , . -—c//(F,x*,u) + u X Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen Newtonova metoda pro systémy rovnic Odhady přesnosti - Tisseurova věta PA081: Programování numerických výpočtů A. Křenek ► pak relativní chyba řešení Newtonovou metodou klesá až dosáhne llJ(X*) ll . /r s. x Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen přesnost závisí na podmíněnosti J 12/25 Newtonova metoda pro systémy rovnic Odhady přesnosti - Tisseurova věta ► pak relativní chyba řešení Newtonovou metodou klesá až dosáhne IIJ(x#) || , . -—c//(F,x*,u) + u x ► přesnost závisí na podmíněnosti J ► nezávisí na 4>, tj. nepřesnosti výpočtu J a řešení lineárních rovnic ► pouze jsou na ně kladeny jisté podmínky ► lze použít nepřesné řešení lineárních rovnic ► Jakobián není třeba počítat v každém kroku PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen Newtonova metoda pro systémy rovnic Odhady přesnosti - Tisseurova věta ► pak relativní chyba řešení Newtonovou metodou klesá až dosáhne IIJ(x#) || , . -—c//(F,x*,u) + u x ► přesnost závisí na podmíněnosti J ► nezávisí na 4>, tj. nepřesnosti výpočtu J a řešení lineárních rovnic ► pouze jsou na ně kladeny jisté podmínky ► lze použít nepřesné řešení lineárních rovnic ► Jakobián není třeba počítat v každém kroku ► nepřesnosti mají vliv na rychlost konvergence PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen Newtonova metoda pro systémy rovnic Shrnutí vlastností ► vlastnosti analogické jednodimenzionální verzi ► rychlá (kvadratická) konvergence ► začne stagnovat při dosažení strojové přesnosti v x nebo F ► kritéria kovergence zpravidla Yj\xí\<£x nebo ^ |Fí | < eF (které nastane dříve) ► špatné globální vlastnosti ► při nepřesném odhadu může snadno zabloudit ► nutná kontrola řešení ► potřebné výpočty derivací ► není tolik jiných metod na výběr ► v nouzi má smysl i aproximace konečnou diferencí ► dopředný (recyklace hodnoty funkce) nebo centrální výpočet PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen Newtonova metoda pro systémy rovnic Shrnutí vlastností ► problém /' — 0 „obohacen" na singularity J ► mnoho řešení - určitě se netrefíme zadne resem - ?? Newtonova metoda pro systémy rovnic Shrnutí vlastností ► problém /' — 0 „obohacen" na singularity J ► mnoho řešení - určitě se netrefíme zadne reseni - ?? ► závažnější numerické potíže ► např. opakované přičtení malé hodnoty by nasměrovalo metodu jinam ► přetečení a podtečení při násobení ► maskovány v maticových operacích - špatně odhalitelné PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen Newtonova metoda pro systémy rovnic Shrnutí vlastností ► problém /' — 0 „obohacen" na singularity J ► mnoho řešení - určitě se netrefíme zadne reseni - ?? ► závažnější numerické potíže ► např. opakované přičtení malé hodnoty by nasměrovalo metodu jinam ► přetečení a podtečení při násobení ► maskovány v maticových operacích - špatně odhalitelné ► správné škálování ► řádově srovnatelné hodnoty nezávislých proměnných ► dtto. funkční hodnoty, i vůči proměnným ► nejlépe v řádu jednotek ► odpovídající nastavení intervalu pro diference místo derivací ► když to nejde - řeším ten správný problém? PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen Převedení na optimalizaci ► optimalizační metody (hledání minima) ► numericky a algoritmicky příjemnější ► zdánlivě složitější problém ► zkusíme řešení rovnic převést na optimalizaci ► definujeme funkci /(x) = |XFivx)2 ► pouze nezáporné hodnoty ► globální minima (0) v kořenech původního systému F PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen Převedení na optimalizaci ► optimalizační metody (hledání minima) ► numericky a algoritmicky příjemnější ► zdánlivě složitější problém ► zkusíme řešení rovnic převést na optimalizaci ► definujeme funkci /(x) = I^F^x)2 ► pouze nezáporné hodnoty ► globální minima (0) v kořenech původního systému F ► optimalizační metody hledají lokální minimum ► taková f(x) jich má zpravidla příliš mnoho ► nutný ještě přesnější odhad řešení ► využití v kombinovaných metodách Newtonova metoda s hledáním po přímce ► krok Newtonovy metody pro řešení F(x) =0 ► používá se k minimalizaci funkce f(x) = Xfí(x)2 ► v případě potřeby se krok zkrátí ► v reálných případech úspěšnější ► exaktní specifikace podmínek konvergence je komplikovaná a v praxi téměř neověřitelná ► metoda může uváznout v lokálním minimu / ► lze snadno detekovat ► je třeba zkusit jiný počáteční odhad Newtonova metoda s hledáním po přímce ► Newtonovský krok Ax je směrem dolů vůči / (V/)rAx = (FrJ)(-r1F) = -FrF < 0 PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen 17/25 Newtonova metoda s hledáním po přímce ► Newtonovský krok Ax je směrem dolů vůči / (V/)rAx = (F^M-J-iF) = -FrF < 0 ► směr Ax je správně ► při dostatečně krátkém kroku hodnota / poklesne ► celý krok může být příliš ► kvadratická aproximace může být nepřesná Newtonova metoda s hledáním po přímce ► Newtonovský krok Ax je směrem dolů vůči / (V/)rAx = (F^M-J-iF) = -FrF < 0 ► směr Ax je správně ► při dostatečně krátkém kroku hodnota / poklesne ► celý krok může být příliš ► kvadratická aproximace může být nepřesná ► volíme krok AAx pro vhodné 0 < A < 1 ► zdánlivě nejlepší je minimalizovat f(x + AAx) vůči A ► výpočetně náročné (počet vyhodnocení F) ► lze ukázat, že to není nezbytné k dosažení rychlé konvergence Newtonova metoda s hledáním po přímce ► Newtonovský krok Ax je směrem dolů vůči / (V/)rAx = (F^M-J-iF) = -FrF < 0 ► směr Ax je správně ► při dostatečně krátkém kroku hodnota / poklesne ► celý krok může být příliš ► kvadratická aproximace může být nepřesná ► volíme krok AAx pro vhodné 0 < A < 1 ► zdánlivě nejlepší je minimalizovat f(x + AAx) vůči A ► výpočetně náročné (počet vyhodnocení F) ► lze ukázat, že to není nezbytné k dosažení rychlé konvergence ► prosté kritérium f(x + AAx) < f(x) nestačí ► je třeba f(x + AAx) < f(x) + a\(V/)rAx ► stačí malé a, např. 10~4 Newtonova metoda s hledáním po přímce ► první pokus - plný krok Ax ► ideální by bylo minimalizovat g (A) = /(x + AAx) ► platí g'(0) = (V/)TAx PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen Newtonova metoda s hledáním po přímce ► první pokus - plný krok Ax ► ideální by bylo minimalizovat g (A) = f(x + AAx) ► platí gf(0) = (Vf)TAx ► aproximujeme g (A) polynomem ► nejprve kvadratický procházející g(0) a g(l) PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen a\2 + g'(0)\ + g(0) minimum je v Ai = -g'(0)/2a Newtonova metoda s hledáním po přímce ► první pokus - plný krok Ax ► ideální by bylo minimalizovat g (A) = f(x + AAx) ► platí gf(0) = (Vf)TAx ► aproximujeme g (A) polynomem ► nejprve kvadratický procházející g(0) a g(l) aA2 + #'(0)A + #(0) minimum je v Ai = -g'(0)/2a ► další kroky kubický polynom procházející g(0) a dvěma nej čerstvějšími A; aA3 + bA2 +gr(0)A + g(0) ► opakujeme do dosažení g (A) < g(0) + otAg'(0) PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen PA081: Newtonova metoda - shrnutí Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda ► rychlá konvergence, špatné globální vlastnosti ► vyžaduje derivace ► globální vlastnosti lze vylepšit řízeným zkrácením kroku ► i tak není metoda 100% spolehlivá ► neobejdeme se bez základní znalosti řešeného problému Vícerozměrná metoda sečen ► jednorozměrná metoda pracuje s aproximací derivace dí+i =- ► analogicky pro více rozměrů Bi+iAXi = AF; ► nedává jednoznačné řešení ► Broydenova metoda ► rekurentní formule Bí+i na základě B* ► minimální změna, která ještě vyhoví kritériu sečny jy _u , (AFj - BjAXj) 0 Axj (AXi)TAXi Vícerozměrná metoda sečen ► neexistuje jednoznačná analogie separace kořenů ► může zabloudit stejně jako Newtonova metoda ► lze aplikovat stejný princip hledání po přímce ► B není přesně Jakobián ► negarantuje, že Ax směřuje dolů vůči / = |FF ► základní předpoklad hledání po přímce ► v případě selhání náhrada aproximací J ► zpravidla dopředně diference PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen Vícerozměrná metoda sečen ► neexistuje jednoznačná analogie separace kořenů ► může zabloudit stejně jako Newtonova metoda ► lze aplikovat stejný princip hledání po přímce ► B není přesně Jakobián ► negarantuje, že Ax směřuje dolů vůči / = |FF ► základní předpoklad hledání po přímce ► v případě selhání náhrada aproximací J ► zpravidla dopředně diference PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen ► problém singularity B (resp. J v Newtonově metodě) ► nelze určit následující krok ► řeší se umělou modifikací B nebo „krokem stranou" ► de-facto restart z jiného, bližšího bodu Vícerozměrná metoda sečen ► neexistuje jednoznačná analogie separace kořenů ► může zabloudit stejně jako Newtonova metoda ► lze aplikovat stejný princip hledání po přímce ► B není přesně Jakobián ► negarantuje, že Ax směřuje dolů vůči / = |FF ► základní předpoklad hledání po přímce ► v případě selhání náhrada aproximací J ► zpravidla dopředně diference ► problém singularity B (resp. J v Newtonově metodě) ► nelze určit následující krok ► řeší se umělou modifikací B nebo „krokem stranou" ► de-facto restart z jiného, bližšího bodu ► existují i komplikovanější metody ► neomezují se jen na hledání po přímce Domácí úkol Návratový modul kosmické lodi Sojuz 28 se vrací z oběžné dráhy a blíží se k jaderné elektrárně Temelín. Určete vzdálenost mezi povrchy modulu a chladící věže v daném okamžiku. ► chladící věž je rotační hyperboloid o rovnici x2+y2 (z-C)2 A2 B2 = 1 ► modul Sojuz 28 je rotační paraboloid „špičkou dolu", tj (x - X)2 + (y - Y)2 =D2(z- Z) PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen ► A, B, C, D, X, Y, Z jsou v daném okamžiku konstanty Domácí úkol ► oba povrchy lze parametrizovat, např. chladící věž: Xt = r cos (x, y t = r sin a kde r = AA /1 + {zt-cy B2 ► podobně pro povrch Sojuzu parametry j8, zs ► vzdálenost dvou libovolných bodů na povrchu je PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen y](xt -xs)2 + (yt -ys)2 + (zt-zs)2 ► hledáme hodnoty a, j8, zt,zs minimalizující tuto vzdálenost ► derivace podle všech parametrů musí být nulová - systém 4 rovnic Domácí úkol ► v oblíbeném programovacím jazyce si najděte vhodnou knihovnu k řešení systémů nelineárních rovnic ► implementujte systém rovnic tak, aby bylo možno knihovnu použít ► vstupem programu jsou parametry A-Z z předchozích vzorců ► případně jednoduchou vizualizací znázorněte postup hledání řešení ► odevzdání 20.4., 1 bod ke zkoušce Shrnutí ► řešení systému nelineárních rovnic je jeden z nej nepříjemnějších problémů ► řešení se vyskytují „náhodně" ► nelze mechanicky separovat kořeny ► vždy je třeba mít základní představu o charakteristice problému ► metoda prosté iterace ► dokážeme-li splnit podmínky konvergence ► rychost konvergence neurčitá ► Newtonova a Broydenova metoda ► vylepšení globálních vlastností hledáním po po přímce ► nejsou zcela spolehlivé, mohou uváznout v lokálním minimu pomocné funkce ► detaily a implementace viz W.H.Press, Numerical Recipes in C ► netřeba se učit vzorce zpaměti ► důležité je znát podstatu problémů a principy metod PA081: Programování numerických výpočtů A. Křenek Metoda prosté iterace Newtonova metoda Hledání po přímce Metoda sečen 25/25