Numerické metody řešení nelineárních rovnic

MASARYKOVA UNIVERZITA Přírodovědecká fakulta

Bakalářská práce

Numerické metody řešení nelineárních rovnic

Studijní program: Aplikovaná matematika Studijní obor: Matematika - ekonomie Brno 2011

Lukáš Jagoš

ii

Prohlašuji, že jsem svou bakalářskou práci napsal samostatně a výhradně s použitím citovaných pramenů. Souhlasím se zapůjčováním práce a jejím zveřejňováním.

Jméno: Lukáš Jagoš Datum: 20.12. 2010 Podpis:

Vedoucí práce: Mgr. Jiří Zelinka, Dr.

iv

v

Abstract

Title: Numerical Methods for solving Nonlinear Equations Author: Lukáš Jagoš Department of Mathematics and Statistics, Faculty of Science, MU Supervisor: Mgr. Jiří Zelinka, Dr. Abstract: This Thesis is a Tool for Students of a Subject Numerical methods. It is divided into four Chapters including Introductory Terms, Solving Nonlinear Equations, Polynoms and Systems of Nonlinear Equations. Chapters cover Rehearsal of Essential Theory for particular Method, Sample Solved Excercise and one Solved Excercise for Comparison of listed Methods. Keywords: numerical, methods, solving, nonlinear, equations

vi

Abstrakt

Název práce: Numerické metody řešení nelineárních rovnic Autor: Lukáš Jagoš Ústav matematiky a statistiky Přírodovědecké fakulty, MU Vedoucí bakalářské práce: Mgr. Jiří Zelinka, Dr. Abstrakt: Bakalářská práce Numerické řešení nelineárních rovnic je pomůckou pro studenty předmětu Numerické metody. Práce je rozdělena do čtyř základních kapitol zahrnujících seznámení se základními pojmy, řešení nelineárních rovnic, polynomů a soustav nelineárních rovnic. Kapitoly obsahují výčet teorie nutné k aplikaci dané metody, ukázkový řešený příklad a na závěr vždy jeden společný příklad pro srovnání uvedených metod. Klíčová slova: vybrané, metody, řešení, nelineárních, rovnic

vii

Obsah

Čestné prohlášení

iv

Abstract

vi

Abstrakt

vii 2 5 5 6

1 Základní pojmy 1.1 Volba počáteční aproximace 1.2 Zastavení výpočtu 1.3 Řád metody 2 Řešení nelineárních rovnic 2.1 Metoda bisekce 2.2 Metoda prosté iterace 2.2.1 Hledání vhodného tvaru iterační funkce 2.3 Newtonova metoda 2.3.1 Iterační metody vyšších řádů 2.4 Metoda sečen 2.5 Shrnutí

7 7 10 12 14 19 19 22

3 Polynomy 3.1 Odhad polohy a počtu kořenů 3.2 Laguerrova iterační metoda 3.3 Graeffova-Lobačovského iterační metoda

26 26 28 30

4 Soustavy nelineárních rovnic

34 viii

Obsah

4.1 4.2

Metoda prosté iterace Newtonova metoda pro systémy nelineárních rovnic

35 37

Závěr

40

Literatura

42

ix

Úvod

Řešení rovnic a jejich soustav patří mezi tradiční součásti středoškolské matematiky. Ta si ale vybírá převážně takové třídy rovnic, které je možné řešit elementárními vzorci. Obecně platí, že kořeny nelineární rovnice f (x) = 0 nelze nalézt explicitním vzorcem. K hledání kořenů tedy používáme metody iteračního charakteru. Vycházíme ze známé počáteční aproximace řešení a z ní stanovujeme novou hodnotu aproximace (iteraci) postupem, kterým se jednotlivé metody odlišují. Přesný popis kroků realizujících numerickou metodu označujeme jako algoritmus numerické metody. Jeden krok algoritmu nazýváme iterací. Obvyklým postupem pro řešení nelineárních rovnic můžeme narazit na některé problémy: • jak poznat vhodný typ numerické metody • jak určit hodnotu počáteční aproximace • jak poznat, zda bude daná metoda konvergovat ke kořenu Protože pro řešení nelineárních rovnic neexistuje univerzální metoda, je vhodné volit metodu, která nejvíce odpovída povaze a dostupné informaci o rozložení a vlastnostech daného řešení.

1

Kapitola 1 Základní pojmy

V rámci textu této práce budeme předpokládat hledání kořenů rovnice f (x) = 0,

(1.1)

kde x ∈ R. Dále budeme předpokládat, že funkce f (x) je na konkrétním intervalu spojitá a má na něm tolik spojitých derivací, kolik jich bude v daném případě potřeba. Následující kapitola představuje seznámení se základními pojmy numerických metod. Definice 1.1. Nechť X je libovolná neprázdná množina a nechť zobrazení % : X × X → R+ splňuje pro všechna x, y, z ∈ X následující axiomy: 1. %(x, y) = 0 ⇔ x = y 2. %(x, y) = %(y, x) 3. %(x, y) + %(y, z) = ρ(x, z) Potom zobrazení % nazveme metrikou a dvojici (X, %) metrickým prostorem. 2


Formální definice metrického prostoru nám dovoluje pracovat v textu s pojmem vzdálenost. Příklad 1.1.

Typy některých metrik X ≡ R,

%(x, y) = |x − y|

X ≡ Rn ,

%(x, y) =

v u n uX t (x

i

− yi )2

i=1

X ≡ Rn ,

%(x, y) = max |xi − yi | i=1,...,n

X ≡ Rn ,

%(x, y) =

n X

|xi − yi |

i=1

Definice 1.2. Bod x ∈ X je limitou posloupnosti {xn }∞ n=1 v metrickém prostoru X, jestliže platí ∀ε > 0 ∃n0 ∈ N takové, že ∀n > n0 : %(xn , x) < ε. Poznámka 1. limitu.

Za konvergentní budeme považovat tu posloupnost, která má

Věta 1.1. (Bolzanova) Nechť funkce f ∈ C[a, b] a nechť f nabývá v koncových bodech intervalu hodnot s opačnými znaménky, tj. f (a)f (b) < 0. Potom uvnitř tohoto intervalu existuje alespoň jeden bod c takový, že f (c) = 0. V případě, že první derivace funkce f má na tomto intervalu konstantní znaménko, pak se zde nachází právě jeden takový bod. Definice 1.3. Mějme metrický prostor (X, %). Zobrazení ϕ : X → X nazveme kontrakcí, jestliže existuje konstanta 0 ≤ k < 1 taková, že platí %(ϕ(x), ϕ(y)) ≤ k.ϕ(x, y) Číslo k nazýváme koeficientem kontrakce. Definice 1.4. Posloupnost {xn }∞ n=1 se nazývá cauchyovská, pokud pro ∀ε > 0 ∃n0 ∈ N takové, že pro ∀n > n0 a pro ∀t ∈ N platí %(xn , xn+t ) < ε.

3


Definice 1.5. Řekneme, že metrický prostor (X, %) je úplný, pokud v něm má každá cauchyovská posloupnost limitu. Věta 1.2. (Banachova) Nechť A : X → X je kontrakce na X a (X, %) je neprázdný úplný metrický prostor. Pak existuje právě jeden prvek x ∈ X takový, že Ax = x. Věta 1.3. Nechť g ∈ C[a, b], g : [a, b] → [a, b], potom g má na tomto intervalu pevný bod ξ. Je-li navíc pro q ∈ [0, 1] | g(x) − g(y) |≤ q | x − y | ,

∀x, y ∈ ha, bi

pak g má na [a, b] jednoduchý pevný bod. Definice 1.6.

Pevný bod ξ funkce g ∈ C [a, b] se nazývá

• přitahující, pokud existuje takové okolí U bodu ξ, že pro ∀x0 ∈ U posloupnost n o∞ xk konverguje k ξ, k=0

• odpuzující, pokud existuje takové okolí V bodu ξ, že pro ∀x0 ∈ V, x0 6= ξ existuje k takové, že xk ∈ / V. Věta 1.4. nazýváme

Nechť g ∈ C[a, b], g : [a, b] → [a, b] a nechť ξ je pevný bod. Potom ξ

• přitahujícím pevným bodem, jestliže g(x) − g(ξ) x−ξ

< 1.

• odpuzujícím pevným bodem, jestliže g(x) − g(ξ) x−ξ

Definice 1.7.

> 1.

Chybou aproximace ε nazýváme odchylku přibližné (aproximo4


vané) hodnoty od skutečné hodnoty. Je-li x přibližná hodnota veličiny a x¯ její skutečná hodnota, pak absolutní chybu aproximace vyjadřuje vztah ε = |x − x¯| . Relativní chyba aproximace δ je určena vztahem δ=

x¯ − x . x¯

1.1 Volba počáteční aproximace Protože se jedná o nelineární rovnice, je možné, že počet řešení a odhad intervalů není na první pohled zřejmý. Pokud je to z hlediska složitosti výpočtu možné, postupujeme buď vyšetřením průběhu funkce užitím znalostí matematické analýzy nebo převedením rovnice f (x) = 0 na ekvivalentní tvar f1 (x) = f2 (x) a hledáním x-ových souřadnic průsečíků grafů funkcí f1 a f2 . Méně praktickou, ale spolehlivou možností by byl výpočet hodnot funkce a vynesení takového množství souřadnic, aby tvar vynášeného grafu odpovídal charakteristice dané funkce.

1.2 Zastavení výpočtu Předpokládejme, že zvolená numerická metoda konverguje a generovaná posloupnost iterací se blíži ke kořenu. Potom pro posouzení přesnosti výpočtu můžeme použít některý ze vztahů

k+1 x

− xk < ε

f (xk )

<ε

xk+1 − xk k x

kde ε je požadovaná přesnost. 5

<ε


1.3 Řád metody n

o∞

generovaná iterační metodou konverguje k x¯. Řekneme, Nechť posloupnost xk k=0 že tato iterační metoda konverguje k x¯ lineárně, pokud existuje σ ∈ (0, 1) takové, že k+1 x − x¯ lim = σ. k→∞ |xk − x ¯| Číslo σ nazýváme řád konvergence. Řekneme, že daná metoda je řádu p (p ∈ R, p > 1), pokud platí lim

k→∞

k+1 x

− x¯

|xk − x¯|p

= L, kde L > 0.

Konvergenci, která je řádu 2 nazýváme kvadratickou, řádu 3 kubickou atd. Poznámka 2. Řád p = 1 nám říká, že chyba se v každém dalším kroku zmenšuje na konkrétní zlomek předchozí chyby, řád p = 2 nám říká, že chyba se v každém dalším kroku zmenšuje na konkrétní zlomek druhé mocniny předchozí chyby atd.

6

Kapitola 2 Řešení nelineárních rovnic

V následující kapitole rozebereme různé typy metod, kterými lze řešit nelineární rovnice f (x) = 0.

2.1 Metoda bisekce Metoda bisekce je jednou z nejstarších metod pro hledání kořenů nelineárních rovnic. Jedná se o velmi jednoduchou a univerzální metodu, nevýhodou je však její relativně pomalá konvergence. Používá se tedy spíše k hledání přibližného řešení. Metoda využívá věty 1.1. • funkce f musí být spojitá na intervalu [a, b] • f (a)f (b) < 0, tedy na intervalu [a, b] se musí nacházet alespoň jeden kořen f . Pokud jsou tyto podmínky splněny, pak je metoda bisekce generuje posloupnost, která konverguje ke kořenu x¯.

7


Nyní můžeme sestrojit posloupnost intervalů (a1 , b1 ) ⊃ (a2 , b2 ) ⊃ . . . , které obsahují hledaný kořen. Intervaly ak+1 , bk+1 , k = 0, 1, . . . určímé tímto způsobem: 1. Jako počáteční interval bereme [a, b]

2. Středem intervalu ak , bk je bod xk+1 , xk+1 =

3. Pokud f (xk+1 ) 6= 0, ak+1 , bk+1 =

ak + b k 2

(2.1)

   ak , xk+1 ,   

když

f ak f (xk+1 ) < 0

     xk+1 , bk ,

když

f ak f (xk+1 ) > 0

4. Pokud f (xk+1 ) = 0, pak xk+1 = x¯. Věta 2.1. Nechť jsou splněny podmínky věty 1.1 a nechť má funkce f na intervalu [a, b] jediný kořen. Potom posloupnost určená vztahem (2.1) konverguje ke kořenu x¯. Důkaz. Pro ohraničené posloupnosti {ak }, {bk } platí, že jsou monotonní, takže konvergentní. Dále platí b 0 − a0 b k − ak = 2k a proto lim (bk − ak ) = 0. Potom k→∞

lim ak = lim bk = lim xk+1 .

k→∞

k→∞

k→∞

Z podmínek pro konvergenci funkce f (spojitost, opačná znaménka funkčních hodnot v koncových bodech intervalu) vyplývá lim f (ak ) = lim f (bk ) = 0 .

k→∞

k→∞

8


Protože platí f (¯ x) = 0, dokázali jsme konvergenci metody bisekce.

Jednou z výhod metody bisekce je to, že jsme při její aplikaci schopni zjistit, jak dobře daná iterace aproximuje hledaný kořen. Bez důkazu uvedeme následující větu. Věta 2.2. Nechť f ∈ C[a, b] a f (a)f (b) < 0 s tím, že funkce f má na [a, b] pouze jeden kořen x¯. Potom posloupnost iterací generovaná metodou bisekce aproximuje hledaný kořen podle vztahu k+1 x

− x¯ ≤

b−a 2n+1

Příklad 2.1. Metodou bisekce najděte odhad kořene rovnice f (x) = 0 s přesností −4 ε < 10 , pokud f (x) = x3 + 5x − 3

4

y

2

0

-2

-4

-2

-1

0

1

2

x

Obr 2.1: Graf funkce f (x) = x3 + 5x − 3

Řešení. Z obr. 2.1 ihned poznáváme, že se kořen bude určitě nacházet v intervalu [0, 1]. Jako počáteční interval tedy zvolíme [a0 , b0 ] = [0, 1]. Připomeňme, že f (0) < 0, f (1) > 0. Posloupnost intervalů zaznamenáme do tabulky.

9


k

ak

bk

xk+1

f (xk+1 )

0

0.00000

1.00000

0.50000

1

0.50000

1.00000

0.750000 1.171875

2

0.50000

0.75000

0.62500

3

0.50000

0.62500

0.562500 −0.00952

4

0, 562500

0.62500

0.593750 0.178070

5

0.562500

0.59375

0.578130 0.083890

6

0.562500

0.57813

0.570315 0.037075

7

0.562500

0.570315 0.566408 0.013750

8

0.562500

0.566408 0.564454 0.00210

9

0.562500

0.564454 0.563478 −0.0037

10

0.563480

0.564454 0.563966 −7.963 × 10−4

−0.37500 0.369141

Po deseti iteracích má interval [a10 , b10 ] délku 0,000974 a 11.iterace x11 aproximuje kořen x¯ s chybou nepřesahující |f (x11 ) − f (x10 )| = 0, 0029. Na tomto příkladu jsme mohli vidět, že metoda bisekce opravdu konverguje ke kořenu x¯ velmi pomalu.

2.2 Metoda prosté iterace Metoda prosté iterace patří mezi základní iterační metody. Je založena na převodu původního tvaru rovnice (1.1) na ekvivalentní tvar x = g(x), který získáme vhodnými úpravami. V tomto tvaru již nehledáme průsečík funkce s osou x, ale pevný bod funkce. Je důležité si uvědomit, že úpravy funkce nejsou jednoznačné a většinou je můžeme provést více způsoby (viz podkapitola 2.2.1). Věta 2.3. Předpokládejme, že existuje funkce g, která pro ∀x ∈ I = [a, b] splňuje Lipschitzovy podmínky uvedené ve větě 1.3. • g(x) ∈ [a, b]

∀x ∈ [a, b]

(2.2)

• |g 0 (x)| ≤ q < 1

∀x ∈ [a, b]

(2.3) 10


Potom rovnice xk+1 = g(xk ) má na intervalu [a, b] jediné řešení a posloupnost n o∞ určená předpisem xk k=0

xk+1 = g(xk )

(2.4)

k němu konverguje pro libovolný počáteční bod x0 ∈ I. Důkaz. Věta je důsledkem Banachovy věty o pevném bodě. Metrikou je %(x, y) = |x − y| a vzhledem k ní je (I, %) úplným metrickým prostorem. Jsou tedy splněny všechny podmínky věty 1.2. Příklad 2.2. Funkce x = g(x), 0 ≤ x ≤ 1 má pevný bod pro každé x ∈ [0, 1]. Funkce g(x) = x − sin(πx) má na intervalu [0, 1] 2 pevné body, v bodě x = 0 a dále v bodě x = 1. Geometrická interpretace této iterační metody je zobrazena na obr. 2.2. Hledat pevný bod funkce geometricky znamená hledat průsečík funkcí y = x a y = g(x).

1.5

Y

1

0.5

0

-0.5

-1

x0

0

B

x Obr 2.2: Metoda prosté iterace xk+1 = sin(xk )

Definice 2.1.

Máme-li funkci g ∈ C [a, b], g : [a, b] → [a, b], potom metoda

11


n

prosté iterace generuje posloupnost xk

o∞ k=0

xk+1 = g(xk )

definovanou pro x0 ∈ [a, b]:

k = 0, 1, 2 . . .

Platí tedy, že x1 = g(x0 ), x2 = g(x1 ) = g(g(x0 )), atd. Pracnost výpočtu metodou prosté iterace závisí na vhodné volbě iterační funkce a počáteční aproximace. Metoda prosté iterace patří mezi jednokrokové iterační metody, neboť výpočet xk+1 závisí pouze na jedné předchozí aproximaci xk .

2.2.1 Hledání vhodného tvaru iterační funkce Při hledání vhodného tvaru iterační funkce je důležité, aby byly vždy splněny předpoklady nutné pro konvergenci (2.2), (2.3). Platí, že různé tvary iteračních funkcí mají různé rychlosti konvergence. Příklad 2.3. 10 = 0. Řešení.

Použitím metody prosté iterace najděte řešení rovnice x3 + 4x2 −

Zadaná rovnice má na intervalu [1, 2] jediný kořen. Algebraickými úpravami sestavíme několik ekvivalentních zápisů této rovnice a porovnáme jejich konvergenční vlastnosti metodou prosté iterace. 1. x = g1 (x) = −x3 − 4x2 + 10 + x 2. x = g2 (x) = 3. x = g3 (x) =

1 2

√ 10 − x3

q

10−4x2 x

4. x = g4 (x) = x −

x3 +4x2 −10 3x2 +8x

Za startovací bod bereme x0 = 1.5. Výsledky opět shrneme do tabulky.

12


k 0

g1 xk

g2 xk

g3 xk

g4 xk

1.5000

1.500000

1.5000

1.50000

1 -0.8750

1.286954

0.8165

1.37333

2

1.402542

2.9969

1.36526

1.345459

∅

1.365230014

6.7324

3 -469.72 4

1.028 × 10

8

1.375171

5

1.360094

6

1.367853

7

1.363887

8

1.365917

9

1.364878

10

1.365414

20

1.365230236

30

1.365230013

1.365230013

Kořenem rovnice x3 + 4x2 − 10 = 0 je číslo x¯ = 1, 365230013. Ke správnému kořenu tedy konvergovala funkce g2 a g4 , g4 navíc velmi rychle. Funkce g1 a g3 divergovala. Důvod je jednoduchý. Pro funkci g1 a g3 nejsou splněny podmínky konvergence. Platí • g 01 (1, 5) = −17, 75 • g 0 2 (1, 5) = −0, 65561 • g 0 3 (1, 5) = −5, 17 • g 0 4 (1, 5) = 0, 1148 Bod x0 = 1, 5 byl tedy pro funkce g1 a g3 odpuzující. Geometrickou interpretaci iteračního procesu jednotlivých funkcí můžeme pozorovat na následujícím obrázku.

13


10 3

8

2

y

y

6

4

1

2

x1

0

x0

0

x1 x2 x0 -2 -1

-10

-8

-6

-4

-2

0

2

4

6

8

-3

10

-2

-1

x

0

1

2

3

x

a)

b) 4

4 3 3

y

y

2 2

1 1

0

x1

0

x0

-1

x1 x0

-1 -4

-3

-2

-1

0

1

2

3

4

-3

-2

-1

0

x

x

c)

d)

1

2

3

Obr 2.3: Metoda prosté iterace pro iterační funkci: a) g1 (x), b) g2 (x) c) g3 (x) d) g4 (x).

2.3 Newtonova metoda Někdy se této metodě říká také metoda tečen. Její princip je odvozen z nahrazení h i křivky v okolí kořene tečnou procházející bodem xk , f (xk ) .

xk+1

Jako počáteční aproximaci volíme x0 . Obecný předpis pro výpočet iterace h i získame tak, že bodem xk , f (xk ) vedeme tečnu ke křivce y = f (x). Průsečíh

i

kem tečny s osou x je pak bod xk+1 , 0 .

14


0.5

0.4

0.3

Y

0.2

0.1

0

x2

x1

x0

-0.1

-0.2 1

2

3

4

5

6

x

Obr. 2.4: Newtonova metoda pro funkci f (x) = x2 − 5.1x + 6.2

Do rovnice tečny y = f (xk ) + f 0 (xk )(x − xk ) dosadíme y := 0, vyjádříme x a položíme x := xk+1 . 0 = f (xk ) + f 0 (xk )x − f 0 (xk )xk x=

f 0 (xk )xk − f (xk ) f 0 (xk )

xk+1 = xk −

f (xk ) f 0 (xk )

(2.5)

Ke stejnému předpisu dospějeme také odvozením z Taylorova rozvoje. Předpokládejme jednoduchý reálný kořen na intervalu I = [a, b]. Pokud na I existují nenulové derivace funkce f , můžeme ji v okolí libovolného bodu x ∈ I rozvinout do Taylorovy řady: f (x) = f (x0 ) + f 0 (x0 )(x − x0 ) +

1 00 f (x0 )(x − x0 )2 + . . . 2!

Nyní proveďme linearizaci prvních dvou členů rozvoje, tzn. v rovnici f (x) = 0 na-

15


hradíme f (x) prvními dvěma členy Taylorova rozvoje a určíme kořen x1 . f (x0 ) + f 0 (x0 )(x − x0 ) = 0 x1 = x 0 −

f (x0 ) f 0 (x0 )

Stejně můžeme postupovat i v okolí bodu xi , i = 1, . . . Obecný předpis bude opět odpovídat vztahu (2.5). Definice 2.2. Za předpokladu, že f 0 (xk ) 6= 0 pro k = 0, 1, . . ., nazýváme metodu určenou vztahem (2.5) Newtonovou metodou. Výpočet ukončíme a xk+1 bereme za dostatečně přesnou aproximaci kořene pokud jsou splněny podmínky uvedené v podkapitole 1.2. Věta 2.4. Nechť f ∈ C 2 [a, b]. Nechť x¯ ∈ [a, b] je kořenem rovnice f (x) = 0 a n o∞ f 0 (¯ x) 6= 0. Pak existuje δ > 0 tak, že posloupnost xk+1 generovaná Newtonovou k=0 metodou konverguje k bodu x¯ pro každou počáteční aproximaci x0 ∈ [¯ x − δ, x¯ + δ] ⊆ [a, b] . Řád metody je 2. Důkaz je uveden v [4]. Monotonní konvergenci Newtonovy metody zajišťují Fourierovy podmínky, které jsou shrnuty v následující větě. Věta 2.5. Nechť funkce f ∈ C 2 [a, b] a nechť rovnice f (x) = 0 má na intervalu [a, b] jediné řešení x¯. Nechť f 0 (x) a f 00 (x) jsou spojité a nemění znaménko na intervalu [a, b], přičemž f 0 (x) 6= 0, ∀x ∈ [a, b]. Nechť počáteční aproximace x0 je ten z krajních bodů [a, b], v němž je znaménko funkce f stejné jako znaménko funkce druhé n o∞ derivace f 00 . Pak posloupnost xk+1 generovaná Newtonovou metodou konverk=0 guje monotonně k bodu x¯ . Důkaz. Pro f (a) < 0 a f 0 , f 00 > 0 platí podle Fourierových podmínek x0 = b. Dokážeme, že potom x¯ < xk+1 < xk .

16


Platí, že 1 x − xk )2 f 00 (y) 0 = f (¯ x) = f (xk ) + (¯ x − xk )f 0 (xk ) + (¯ 2 Vzhledem k předpokladu nezápornosti f 00 (y) musí platit f (xk ) + (¯ x − xk )f 0 (xk ) < 0 Z čehož plyne f (xk ) f 0 (xk ) < xk .

x¯ < xk − xk+1 Indukcí se dokazují ostatní případy.

Newtonovou metodou nalezněte kořen rovnice f (x) = x3 −

Příklad 2.4.

√

6.

8 6 4

y

2 0 -2 -4 -6 -8 -3

-2

-1

0

1

2

3

x

Obr. 2.5. Funkce f (x) = x3 −

√

6

Dosazením funkce f do předpisu (2.5) dostáváme Newtonovu metodu ve tvaru

xk+1

√ k 3 (x ) − 6 = xk − . 2 k 3 (x )

Jako počáteční aproximaci volíme bod x0 = 2. Posloupnost iterací sepíšeme do tabulky.

17


k

xk

f (xk )

0

2.000000

5.550510

1

1.537457

1.184715

2

1.370392

0.124072

3

1.348370

0.001983

4

1.348006

−8.428 × 10−7

Je patrné, že již u 4. iterace získáváme velmi dobrý odhad kořene x¯. Newtonova metoda tedy konvergovala dostatečně rychle. Na tom samém příkladu ještě otestujeme, jak bude Newtonova metoda konvergovat v případě, že počáteční aproximaci zvolíme x0 = 13 , x0 = 9, x0 = 15 . x(k)

f (xk )

xk

-2.412453

4.000000

61.550510

1 5

7.570691

431.468485

2.717698

17.623102

20.545748

8670.480000

2

5.061373

127.210240

1.922347

4.654383

13.699094

2568.400000

3

3.406121

37.067186

1.502513

0.942501

9.137083

760.370000

4

2.341125

10.381907

1.363350

0.084600

6.101169

224.660000

5

1.709722

2.548283

1.348178

0.000938

4.089381

65.940000

2.775078

18.920000

1.956076

5.035000

1.517445

1.044640

9

1.366222

0.100647

10

1.348248

0.001318

k

x(k)

0

1 3

1

f (xk )

6

1.419136

0.408575

7

1.351512

0.019161

8

−7

1.348006 −8.428 × 10 −7

1.348006 −8.428 × 10

18

f (xk ) -2.441490


2.3.1 Iterační metody vyšších řádů Nyní uvedeme dvě iterační metody, jejichž řád konvergence je v případě jednoduchého kořene větší než 2. Uvažujme nejprve metodu určenou předpisem xk+1 = xk −

[f (xk )]2 f 00 (xk ) f (xk ) − f 0 (xk ) 2[f 0 (xk )]3

(2.6)

Řád této metody je 3. Další z používaných metod je tzv. Halleyova metoda xk+1 = xk −

f (xk )f 0 (xk ) , [f 0 (xk )]2 − 12 f (xk )f 00 (xk )

(2.7)

která konverguje rovněž kubicky. Obě tyto metody jsou příkladem toho, že jednobodová iterační metoda vyžaduje k výpočtu jedné iterace hodnoty tím vyšších derivací, čím vyššího je řádu. Protože výpočty vyšších derivací mohou být v některých případech náročné na použití, můžeme dosáhnout zjednodušení aproximací pomocí prvních derivací nebo funkčních hodnot (viz metoda sečen). Zjednodušení metody (2.6) dosáhneme následující úpravou fˆ00 (xk ) = −

(xk

6 2 [f (xk ) − f (xk−1 )] + k [2f 0 (xk ) + f 0 (xk−1 ] k−1 2 −x ) x − xk−1 h

k

xk+1 = xk −

i2

f (xk )

f (x ) − fˆ00 (xk ). f 0 (xk ) 2[f 0 (xk )]3

Tento vzorec závisí, podobně jako v případě Newtonovy metody, na f a f 0 a pro řád konvergence platí, že je přibližne roven číslu 2,732. Podmínkou je ovšem zmiňovaná existence jednoduchého kořene.

2.4 Metoda sečen Z důvodu výpočtu první derivace u Newtonovy metody je možné, že výpočet bude mnohem obtížnější, než tomu bylo v případě příkladu 2.4. Je-li zadaná rovnice na výpočet derivace náročná, můžeme dosáhnout zjednodušení vhodnou úpravou New19


tonovy metody. Uvažujme následující aproximaci:

f 0 (xk ) ≈

f (xk+1 ) − f (xk ) xk+1 − xk

k = 0, 1, ...

Po dosazení f 0 (xk ) do vztahu (2.5) získáváme vzorec pro výpočet iterací metodou sečen. Při výpočtu vycházíme ze dvou počátečích aproximací x0 , x1 s tím, že každou další aproximaci počítáme podle předpisu

xk+1 = xk −

xk − xk−1 f (xk ) . f (xk ) − f (xk−1 )

(2.8)

Definice 2.3. Iterační metodu vyjádřenou vztahem (2.8) nazýváme metodou sečen. Název této metody pochází z její geometrické interpretace: xk+1 je označení xové souřadnice průsečíku osy x a přímky vedené body[xk−1 , f (xk−1 )] a [xk , f (xk )]. y = f (xk ) +

f (xk ) − f (xk−1 ) (x − xk ) = 0 xk − xk−1 ⇒ x = xk+1

Konvergence metody sečen je ilustrována na obr. 2.6.

20


1.5

1.0

y

0.5

x0

x1 x4

0.0

x3

x2

-0.5

-1.0

-1.5 -2

-1

0

1

2

x

Obr. 2.6. Metoda sečen pro funkci y = arctan(x), x0 = −1.7, x1 = −0.4.

Věta 2.6. Nechť rovnice f (x) má na intervalu [a, b] kořen x¯ a nechť derivace f 0 a f 00 jsou v okolí bodu x¯ spojité, přičemž f 0 (x) 6= 0 v tomto okolí. Potom posloupnost určená metodou sečen konverguje pro dostatečně blízké počáteční aproximace x0 a x1 ke kořenu x¯. Řád metody je přibližně 1,618. Důkaz je uveden v [6]. Příklad 2.5. Metodou sečen najděte kořen funkce f (x) = xsin(x) − x3 + 4 ležící v intervalu I = [1, 5; 2, 5]. Řešení. x0 = 1, 5 x1 = 2, 5 Dosazováním do (2.8) získáváme následující iterace

21


f (xk )

k

x

0

1.500000

2.121243

1

2.500000

−10.128820

2

1.673162

0.980442

3

1.746133

0.395431

4

1.795458

−0.037631

5

1.791173

0.001237

6

1.791309 3.66332 × 10−6

Přicházíme k závěru, že metoda sečen konverguje také dostatečně rychle.

2.5 Shrnutí Metoda bisekce byla první z uvedených metod. Její výhodou jsou především nízké nároky na vlastnosti vyšetřované funkce a dále fakt, že při splnění podmínek konverguje vždy. Je to metoda jednoduchá a efektivní. V případě, že se na počátečním intervalu nachází více kořenů, nalezneme metodou bisekce vždy jen jeden. Hlavní nevýhodou je její pomalá konvergence, jak jsme se mohli přesvědčit na příkladu 2.1. Další z uváděných metod byla metoda prosté iterace. Rychlost její konvergence záleží na vhodně zvoleném tvaru iterační funkce. Výhodou je zejména nenáročnost na výpočet. Newtonova metoda je ze všech uváděných metod patrně nejznámější. Její podmínkou je jednoduchý kořen, ke kterému po splnění daných podmínek konverguje dostatečně rychle. Jak jsme mohli vidět, rychlost konvergence je závislá na vhodné volbě počáteční aproximace. Pokud derivaci ze vztahu pro výpočet iterací Newtonovou metodou nahradíme diferencí, získáme metodu sečen. Ta konverguje vždy, pokud zvolíme počáteční hodnoty dostatečně blízko kořene. Výhodou této metody je vysoká efektivnost. Při výpočtu není potřeba pracovat s derivací, ve srovnání s Newtonovou metodou je tedy jednodušší na výpočet. Mezi její nevýhody patří opět podmínka, že hledaný kořen dané rovnice musí být jeden a prostý. Srovnání metod nyní předvedeme na příkladu.

22


Příklad 2.6. Najděte odhad kořene rovnice −6x4 − 2x2 − x + 1 = 0. ležícího v intervalu srovnání. Řešení.

h

1 3 , 4 4

i

. K výpočtu použijte probírané metody a proveďte jejich

Metoda bisekce.

k

ak

bk

xk+1

f (xk+1 )

1

0.250000

0.75000

0.500000 −0.375000

2

0.250000

0.50000

0.375000

3

0.375000

0.50000

0.437500 −0.040131

4

0.375000 0.437500 0.406250

0.100245

5

0.406250 0.437500 0.421875

0.032110

6

0.421875 0.437500 0.429688 −0.003485

7

0.421875 0.429688 0.425782

0.014440

8

0.425782 0.429688 0.427735

0.005511

9

0.427735 0.429688 0.428712

0.001019

0.225098

10 0.428712 0.429688 0.429200 −0.001232

Metoda prosté iterace. Iterační funkci zvolíme ve tvaru s

g(x) =

−x + 1 6x2 + 2

23


k

g(xk )

x(k)

f (xk )

1

0.500000 0.377965 −0.375000

2

0.377965 0.466596

3

0.466596 0.401660 −0.186410

4

0.401660 0.448997

5

0.448997 0.414336 −0.096045

6

0.414336 0.439643

7

0.439643 0.421122 −0.050372

8

0.421122 0.434655

9

0.434655 0.424755 −0.026661

0.213870 0.119513 0.065483 0.035486

10 0.424755 0.431991

0.019111

Newtonova metoda.

x

k+1

k

=x −

−6 xk

4

− 2 xk

2

− xk + 1

−24 (xk )3 − 4xk − 1

k

xk+1

f (xk+1 )

0

0.250000

0.601563

1

0.503289

−0.394853

2

0.438269

−0.043796

3

0.429094 −7.42 × 10−4

4

0.428933

3.53 × 10−7

Metoda sečen. x0 = 0.25,

24

x1 = 0.75


k

xk+1

f (xk+1 )

1

0.339120

0.351486

2

0.464410

−0.174841

3

0.422789

0.002800

4

0.428533

0.001843

5

0.428937 −1.809×10−5

6

0.428933

3.526×10−7

7

0.428933

3.526×10−7

Vyhodnocení. • Nejrychleji konvergovala jednoznačne Newtonova metoda, pro tento typ úloh se tedy hodí nejvíce. V případě, že by zadaná rovnice byla složitější na výpočet derivace, tak bychom mohli úspěšně využít metodu sečen. Ta konvergovala také velmi rychle. • Rychlost Newtonovy metoda i metody sečen je podmíněna vhodnou volbou počátečních aproximací. V praxi se může stát, že při jejich nesprávné volbě budou dané metody konvergovat pomaleji. • Metoda bisekce byla jednoduchá na výpočet, konvergovala ale velmi pomalu. Je ideální pro nalezení dobré počáteční aproximace. • Metoda prosté iterace konvergovala nejpomaleji. Zrychlení bychom mohli dosáhnout jinou volbou iterační funkce.

25

Kapitola 3 Polynomy

Polynomem nazýváme výraz tvaru Pn (x), kde Pn (x) = an xn + an−1 xn−1 + ... + a1 x + a0 ,

n ∈ N.

Platí, že pro odhad reálných kořenů polynomu můžeme použít všechny dosud uvedené iterační metody. Polynomy se ale často vyskytují v nejrůznějších aplikacích a mají své vlastní specifické vlastnosti, věnujeme jim tedy tuto samostatnou kapitolu.

3.1 Odhad polohy a počtu kořenů Označme x¯1 , x¯2 , ..., x¯n kořeny polynomu Pn (x). Nechť A = max (|an−1 | , ..., |a0 |) , B = max (|an | , ..., |a1 |) ,

26

Kapitola 3 Polynomy

kde ai , i = 0, ..., n, a a0 an 6= 0 jsou koeficienty polynomu P . Potom pro všechny kořeny x¯i polynomu P platí\gg A 1 xi | ≤ 1 + . B ≤ |¯ |an | 1 + |a0 | Příklad 3.1.

Odhadněte polohu kořenů pro rovnici P5 (x) = 0, kde P5 (x) = x5 + 101x4 + 420x3 − 425x2 + 111x − 16.

Pro P5 (x) máme A = 425, B = 425 a po dosazení získáme 1 ≤ |¯ xi | ≤ 1 + 425 1 + 425 16 0, 03628 ≤ |¯ xi | ≤ 426. Absolutní hodnota všech kořenů polynomu P5 (x) bude tedy patřit do intervalu [0, 03628; 426]. Bez důkazu uvedeme ještě větu 4.1, které se říká také Cauchyova věta o poloze kořenů. Věta 3.1. vině, kde

Všechny kořeny polynomu Pn (x) náleží do kruhu τ v komplexní ro ak κk = max . 0≤k≤n−1 an

τ = {z ∈ C : |z| ≤ 1 + κk }

Nyní budeme věnovat pozornost určení počtu kořenů polynomu Pn (x). Definice 3.1. Označme M (x) = Pn (x), M1 (x) = M 0 (x), M2 (x) zbytek po dělení M (x) : M1 (x) násobený číslem (-1), M3 (x) zbytek po dělení M1 (x) : M2 (x) násobený číslem (-1),...,Mr (x) zbytek po dělení Mr−2 (x):Mr−1 (x) násobený číslem (-1), přičemž zbytek Mr+1 (x) po dělení Mr−1 (x):Mr (x) je nulový. Potom posloupnost M (x), M1 (x), ..., Mr (x) nazýváme Sturmovou posloupností.

27

Kapitola 3 Polynomy

Věta 3.2. Pokud platí, že polynom Pn (x) nemá násobné kořeny a Pn (a), Pn (b) 6= 0, potom počet reálných kořenů daného polynomu ležících v intervalu (a, b) je roven N (a) − N (b), kde N (α) značí počet znaménkových změn ve Sturmově posloupnosti pro x = α při vynechání nulových prvků. Důkaz je uveden v [6]. Poznámka 3. Pokud nebude v textu uvedeno jinak, budeme v následující části předpokládat reálné koeficienty polynomu P .

3.2 Laguerrova iterační metoda Předpokládejme prosté, reálné kořeny polynomu Pn (x), označme je λ1 , λ2 , ..., λn . Tyto kořeny si seřadíme vzestupně podle velikosti: λ1 < λ2 < ... < λn . Nyní definujme interval Ii = [λi , λi+1 ] , i = 0, ..., n, kde koeficienty λ0 = −∞, λn+1 = +∞. Potom pro každé λ 6= λj , j = 0, ..., n + 1 existuje právě jedno i takové, že platí λ ∈ Ii . Nyní sestrojíme takovou kvadratickou funkci, jejíž oba kořeny budou ležet v intervalu Ii a v bodě λ bude nabývat záporné hodnoty1 . Princip Laguerrovy metody spočívá v sestrojení takové křivky kvadratické funkce, která by osu x protla co nejblíže krajním bodů intervalu Ii . Tímto odvozením dostáváme předpis Laguerrovy iterační metody ve tvaru λk+1 = λk −

nPn (λk ) Pn0 (λk ) ±

q

Γ(λk )

,

(3.1)

kde Γ(λ) = (n − 1)((n − 1)(P 0n (λ))2 − nPn (λ)Pn00 (λ)). Poznámka 4. Znaménko před odmocninou volíme rovné znaménku P 0 (λk ). Platí-li pro první iteraci P 0 (λk ) = 0, volíme znaménko libovolně. 1

Platí, že takovýchto funkcí můžeme sestrojit nekonečně mnoho.

28

Kapitola 3 Polynomy

Věta 3.3. Mějme reálný polynom Pn (x) stupně n ≥ 1, jehož kořeny jsou prosté a reálné. Nechť λ0 je počáteční aproximace a P (λ0 ) 6= 0. Potom posloupnost iterací generovaná předpisem (4.1) konverguje k některému z kořenů Pn (x). Důkaz je uveden v [3]. Příklad 3.2.

Určete kořeny rovnice x4 − 2x3 + 6x2 + 2x − 3 = 0.

Řešení: Zádán máme polynom 4. stupně. Začneme výpočtem Laguerrovy metody. Za počáteční aproximaci zvolíme λ0 =0,9. λk+1 = λk −

nPn (λk ) Pn0 (λk ) ±

q

Γ(λk )

,

Γ(λk ) = (n − 1)[((n − 1)(P 0n (λk ))2 − nPn (λk )Pn00 (λk )].

Posloupnost iterací zapíšeme do tabulky. k

λk

P (λk )

0

0.9000

2.85810

1

0.7209

1.08077

2

0.6417

0.39516

3

0.6106

0.14190

4

0.5989

0.04891

5

0.5949

0.01741

6

0.5934

0.00563

7

0.5930

0.00250

8

0.5928 9.2744 × 10−4

29

Kapitola 3 Polynomy

3.3 Graeffova-Lobačovského iterační metoda Tato iterační metoda slouží k přibližnému určení polohy všech kořenů polynomu. Budeme předpokládat a0 6= 0. Upravme polynom P na rozklad kořenových činitelů an xn + an−1 xn−1 + ... + a1 x + a0 = an (x − x1 )(x − x2 )...(x − xn ). Roznásobením pravé strany a následným srovnáním koeficientů získáme tzv. Vietovy vzorce an−1 = −(x1 + x2 + ... + xn ), an an−2 = (−1)2 (x1 x2 + x1 x3 + ... + xn−1 xn ), an ··· a0 = (−1)n (x1 x2 ...xn ). an Princip metody předvedeme na příkladu. Uvažujme následující rovnici a2 x2 + a1 x + a0 = 0. Podle Vietových vzorců platí x1 + x2 = −

a1 , a2

x1 x2 =

a0 . a2

(3.2)

V případě, že se kořeny x1 a x2 od sebe výrazně liší (tj. x1 x2 ), platí x2

x1

potom

1,

a1 , a2 x2 a1 x1 1 + =− . x1 a2 x1 + x 2 = −

Zanedbáním členu xx21 získáme x1 ≈ − aa21 a následným dosazením do vztahu x1 x2 = dostáváme hodnotu a0 x2 ≈ − , a1 30

a0 a2

Kapitola 3 Polynomy

kterou můžeme považovat za první přiblížení. Nyní umocníme vztah (4.2), úpravou dostaneme a2 − 2a0 a2 a2 x21 + x22 = 1 2 , x21 x22 = 02 . a2 a2 Protože součet i součin druhých mocnin kořenů známe, můžeme sestavit rovnici b2 y 2 + b1 y + b0 = 0, která má kořeny y1 = x21 a y2 = x22 , kde b2 = a22 , b1 = −(a21 − 2a0 a2 ), b0 = a20 . Analogicky určíme i jejich aproximace y1 ≈ −

b1 b0 , y2 ≈ − , b2 b1

z čehož vyplývá |x1 | ≈

Poznámka 5. ještě přesnější.

v u u b1 t , b2

|x2 | ≈

v u u b0 t . b1

x2 Pro xx21 1 platí, že zlomek x22 bude ještě menší a aproximace

1

Nyní jsme popsali základní princip Graeffovy-Lobačovského iterační metody a můžeme tak shrnout obecný postup. Hledáme aproximace kořenů x¯1 , ..., x¯n rovnice Pn (x) = 0. Předpokládáme, že kořeny jsou dobře rozlišitelné, tzn. jsou reálné a prosté. Vytvoříme posloupnost polynomů P 0 , P 1 , ..., P k , kde P 0 = Pn (x), P k = akn xn + akn−1 xn−1 + ... + ak1 x + ak0 2 akn = (ak−1 n ) , 2 k−1 k−1 −akn−1 = (ak−1 n−1 ) − 2an an−2 , 2 k−1 k−1 k−1 k−1 akn−2 = (ak−1 n−2 ) − 2an−1 an−3 + 2an an−4 ,

... 2 k−1 k−1 (−1)n−1 ak1 = (ak−1 1 ) − 2a2 a0 ,

31

Kapitola 3 Polynomy

2 (−1)n ak0 = (ak−1 0 ) .

Pro kořeny platí |¯ xj | =

v u k u k an−j 2t k , an−j+1

j = 1, 2, ..., n.

)2 ∀i = 0, 1, 2, ... Výpočet zastavíme, pokud aki ≈ (ak−1 i Příklad 3.3. Použitím Graeffovy-Lobačovského metody určete kořeny rovnice x4 − 4x3 + 3x2 + 2x − 6 = 0. Řešení.

Nejprve sestavíme odpovídající posloupnost polynomů s koeficienty a0 , ..., a4

k

a4

a3

a2

0 1 2 3

1 1 1 1

−4 −10 −74 −6594

a1

3 2 13 −40 −559 −664 218801 −1889824

a0 −6 36 1296 1679616

A dále pokračujeme výpočtem s

|α1 | ≈

8

6594 . = 3.0019 1

dosazením do rovnice polynomu ze zadání dostáváme P (α1 ) ≈ 0.00381, P (−α1 ) ≈ 204.44 Za první odhad tedy bereme x1 = 3.0019. Spočítáme další iteraci s

|α2 | ≈

8

1679616 . = 0.9854 1889824

P (0.9854) ≈ −4.00639, P4 (−0.9854) ≈ −0.28755

32

Kapitola 3 Polynomy

Vidíme, že v dalších sloupcích se mění znaménka a je tedy třeba hledat komplexní kořeny. Pro náročnost výpočet nyní ukončíme. Jak je vidět, aplikace této metoda je poměrně komplikovaná. Její výhodou je to, že nevyžaduje odhad počáteční aproximace. Jejím použitím můžeme navíc hledat více než 1 kořen. V praxi bychom tedy mohli v této části výpočtu vzít hodnoty α1 a −α2 a použít je jako dobré startovací aproximace pro některou z rychleji konvergujících metod.

33

Kapitola 4 Soustavy nelineárních rovnic

Nyní se budeme zabývat soustavami nelineárních rovnic o n neznámých f1 (x1 , ..., xn ) = 0 .. .

(4.1)

fn (x1 , ..., xn ) = 0 Soustavu (4.1) můžeme zapisovat také ve vektorovém tvaru: f (x) = 0, kde



 

f (x , ..., xn ) f (x)  1 1   1 .. ..    = f (x) =  . .    fn (x1 , ..., xn ) fn (x)

    



a



0    ..  0 =  .  ∈Rn 

0



x1 , ..., x¯n )T , pro který platí Řekneme, že řešením soustavy (4.1) je každý vektor x =(¯ f (x) = 0. V textu budeme dále předpokládat, že funkce (4.1) jsou spojité a mají v daném případě tolik spojitých derivací, kolik jich je potřeba. Při hledání kořenů soustav nelineárních rovnic je získání řešení podmíněno splněním konvergenčních podmínek, proto je důležité umět polohu vektoru x buď vhodně 34


odhadnout nebo se k ní umět přiblížit. Pro soustavy nelineárních rovnic neexistuje univerzální metoda, která by se dala použít ke spolehlivému odhadu počáteční aproximace. Více je o této problematice pojednáno v [2]. Definice 4.1. Metriku v Rn definujeme jako %(x, y) = max |xi − yi |. 1≤i≤n

Prostor s touto metrikou se nazývá úplný metrický prostor. Poznámka 6. O ukončení výpočtu se rozhodujeme na základě některého z kriterií

k+1

x − xk

≤ ε

k+1

x

− xk ≤ ε xk

f (xk+1 )

≤ ε,

kde ε je požadovaná přesnost.

4.1 Metoda prosté iterace Věta 4.1. Nechť (M, %) je úplný metrický prostor s metrikou %(x, y) a zobrazení G : M → M má vlastnost %(G(x), G(y)) ≤ k%(x, y) pro ∀ x, y ∈ X ⊂ M, 0 ≤ k < 1. Pak soustava rovnic x = G(x) má na X jediné řešení x ¯. Důkaz plyne z aplikace Banachovy věty o pevném bodě. Požadované řešení x ¯ lze získat následujícím algoritmem: 1. Zvolíme vhodnou počáteční aproximaci x0 ∈ X

35


2. xk+1 = G(xk ) pro n = 0, 1, 2, . . . Jestliže je funkce G diferencovatelná, pak můžeme podmínku 2 nahradit podmínkou ∃q : kG0 (X)k ≤ q < 1 ∀X ∈ D kde G0 (X) je matice s prvky ai,j = vzorce

∂ϕ(x1 ,...,xn ) . ∂xj

V případě dvou rovnic můžeme použít





d −b  1  G0 (X) = , ad − bc −c a kde





a b  G(X) =  . c d Příklad 4.1.

Metodou prosté iterace určete řešení soustavy rovnic x = 0, 2 + 0, 1(−xy 2 + 3x) y = 0, 6 + 0, 1(−x2 y 3 − 2y)

v oblasti Ω = {[x, y]; 0 ≤ x ≤ 1, 0 ≤ y ≤ 1}. Řešení. Platí, že pro [x, y] ∈ Ω je [g1 (x, y), g2 (x, y)] ∈ Ω, 1. podmínka je tedy splněna. Dále je nutno sestavit si Jacobiovu matici 

G0 (X) = 

∂g1 (x,y) ∂x ∂g2 (x,y) ∂x

∂g1 (x,y) ∂y ∂g2 (x,y) ∂y







−0, 1y 2 + 0, 3 −0, 2xy =  −0, 2xy 3 −0, 3x2 y 2 − 0, 2

a pro ni provést odhad např. v k·k∞ normě.

kG0 (X)k∞ = max{ −0, 1y 2 + 0, 3 + |−0, 2xy| ; −0, 2xy 3 + −0, 3x2 y 2 − 0, 2 } Zřejmě platí kG0 (X)k∞ ≤ max{|0, 3 + 0, 2; 0, 2 + 0, 5|} = 0, 7

pro X = [x, y] ∈ Ω.

Pro q = 0, 7 je tedy splněna i druhá podmínka. Počáteční aproximace zvolíme

36


[x0 , y 0 ] = [0, 0]. x1 = g1 (x0 , y 0 )

y 1 = g2 (x0 , y 0 )

x2 = g1 (x1 , y 1 )

y 2 = g2 (x1 , y 1 )

x3 =

y3 =

...

...

Posloupnosti iterací zaznamenáme do tabulky.

k

xk

yk

xk − xk−1

y k − y k−1

0

0.00000

0.00000

1

0.20000

0.60000

0.200000

0.600000

2

0.25280

0.48086

0.052800

-0.119136

3

0.26999

0.50454

0.017195

0.023638

4

0.27413

0.50003

0.004131

-0.004509

5

0.27538

0.50094

0.001258

0.009050

6

...

...

...

...

15 0.275832 0.500796

−9

2.936 × 10

16 0.275832 0.500796 8.041 × 10−10

1.099 × 10−10 −1.012 × 10−12

4.2 Newtonova metoda pro systémy nelineárních rovnic Odvození Newtonovy metody v n dimenzích se dá opět provést pomocí Taylorova rozvoje. Uvažujme soustavu nelineárních rovnic (4.1) s předpokládaným řešením x =(¯ x1 , ..., x¯n )T a počáteční aproximací x0 = (x01 , .., x0n )T . Linearizací (analogicky jako v jednorozměrném případě) dospějeme ke vztahu f 0 (xk )(xk+1 −xk )+f (xk ) = 0

37

(4.2)


kde f 0 (x) je Jacobiova matice funkce f (x) 

∂f1 (x) ∂x1

··· .. .

∂f1 (x) ∂xn

∂fn (x) ∂x1

...

∂fn (x) ∂xn

.. .

  

f 0 (x) = 

.. .

    

Vztah (4.2) nyní upravíme do tvaru, s nímž můžeme počítat jednotlivé iterace:

xk+1 = xk − f 0 (xk )f (xk ),

k = 0, 1, . . .

(4.3)

¯ , nechť f 0 (x) je regulární matice Věta 4.2. Nechť soustava f (x) = 0 má řešení x se spojitými prvky v okolí O(¯ x), přičemž

−1

0

(f (x))

∞

≤ k,

k = konst.

pro všechna x z tohoto okolí. Dále nechť funkce fi , i = 1, ..., n mají spojité druhé parciální derivace v okolí O(¯ x). Potom při dostatečně blízké počáteční aproximaci x0 ¯. bude posloupnost generovaná Newtonovou metodou konvergovat ke kořenu x Důkaz je uveden v [6]. Příklad 4.1.

Newtonovou metodou určete řešení soustavy rovnic g1 (x, y) = 0, g2 (x) = 0

ležící ve čtverci Ω = {[x, y]; −2 ≤ x − 1, 1 ≤ y ≤ 2}, kde

g1 (x, y) = x3 − xy 2 − 1 g2 (x, y) = y 3 − 2x2 y + 2 Řešení.

38


Iterace budeme počítat podle předpisu (4.3). Položme 



x3 − xy 2 − 1  f (x, y) =  3 y − 2x2 y + 2 a





3x2 − y 2 −2xy  f 0 (x, y) =  . −4x 3y 2 − 2x2 Budeme začínat s počáteční aproximací [x0 , y 0 ] = [−1, 1].

k

xk

yk

g1 (xk , y k )

g2 (xk , y k )

0

-1.000000 1.000000

-1.000000

1.000000

1

-1.500000 2.000000

1.625000

1.000000

2

-1.379562 1.673966

0.240186

0.318969

3

-1.392137 1.629876

0.000181

0.012211

4

-1.394071 1.631182

-0.000004

-0.000014

5

−6

-1.394069 1.631182 2.411 × 10

4.428 × 10−6

Jak již bylo zmíněno, problematika stanovení počáteční aproximace pro tuto metodu je mnohem náročnější, než v případě nelineárních rovnic. Při aplikaci Newtonovy metody pro systémy nelineárních rovnic je navíc potřeba na každém kroku iterace řešit soustavu (4.2), z čehož vyplývá náročnost této metody na praktické použití.

39

40

Závěr

Hlavním cílem této práce bylo seznámit čtenáře s běžně používanými metody pro řešení nelineárních rovnic. Existuje spousta dalších metod, které v textu nebyly zmiňovány, protože by tím byl výrazně překročen plánovaný rozsah práce. V úvodní části jsme se věnovali základním definicím a pojmům. Protože uváděné věty nebyly cílem této práce, byly vypsány bez důkazů. Pro bližší seznámení se doporučuje literatura [5], [1].

41

Literatura

[1] Horský, Z., Vektorové prostory, SNTL: MVšT-II, 1980. [2] Keller, H.B., Isaacson, E., Analysis of numerical methods, New York: Dover publications inc., 1994. [3] Mekwi, W.R., Iterative methods for roots of polynomials, Trinity: University of Oxford, 2001. [4] Stoer, J., Bulirsch,R., Introduction to numerical analysis, New York: Heidelberg, Berlin, 1980. [5] Šalát, T., Metrické priestory, Bratislava: Alfa, 1981. [6] Zelinka, J., Horová, I., Numerické metody, Brno: Masarykova univerzita, 2008.

42

Numerické metody řešení nelineárních rovnic

Recommend Documents