Mariusz PYRZ
SIMR (PW), Instytut Pojazdów
Metody numeryczne w mechanice
Rozwiązywanie równań nieliniowych
2.
Rozwiązywanie równania nieliniowego
(z jednÄ… niewiadomÄ…)
Problem:
Dana funkcja f : R R, regularna w przedziale [a,b].
Wyznaczyć pierwiastek (pierwiastki) równania f(x)=0.
Iteracyjna metoda poszukiwania rozwiÄ…zania:
Iteracyjna metoda poszukiwania rozwiÄ…zania:
wyznaczanie i poprawianie kolejnych przybli\eń x1, x2, x3 & , zbiegających się do rozwiązania x*.
Numeryczne poszukiwanie rozwiÄ…zania:
Jak generować kolejne przybli\enia rozwiązania?
Od jakiego punktu rozpocząć poszukiwanie?
Kiedy zatrzymać algorytm?
Jak szybko zmierzamy do rozwiązania (czy mo\na ulepszyć procedurę)?
2
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Rys. 2.1
3
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Lokalizowanie pierwiastków
Omawiane metody numeryczne przeszukują przedział [a,b] zakładając, \e istnieje tam
jeden pierwiastek równania. Przed przystąpieniem do rozwiązywania problemu w
którym występuje wiele pierwiastków mo\na przeprowadzić próbę lokalizacji i
odseparowania pierwiastków, tzn. określenie przedziałów [a,b] do oddzielnego
przeszukiwania.
Do rozdzielenia obszarów występowania pierwiastków mo\na wykorzystać
twierdzenie Rolle a:
twierdzenie Rolle a:
Je\eli f(x) jest określona i ciągła w [a,b], i je\eli N oznacza
liczbę pierwiastków f(x)=0 w rozpatrywanym przedziale, wówczas:
Je\eli f(a)f(b)<0 to f(x) posiada co najmniej jeden pierwiastek w [a,b]
(N jest wówczas liczbą nieparzystą). Je\eli ponadto f (x)`"0 dla ka\dego x z
przedziału [a,b] , to istnieje w nim tylko jeden pierwiastek (N=1).
Je\eli f(a)f(b)>0 to f(x) nie posiada pierwiastków w [a,b] (N=0), lub w
rozpatrywanym przedziale występuje ich parzysta liczba.
4
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Rys. 2.2
5
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda kolejnych przybli\eń
Polega na budowaniu ciÄ…gu kolejnych
x1 = F ( x0 ,...)
przybli\eń rozwiązania, zaczynając od
x2 = F ( x1,...)
wartości początkowej x0
odpowiednio wybranej w [a,b].
ï"
xn = F ( xn -1,...)
Sposób rozwiązywania:
- określ granice przedziału w którym szukamy pierwiastka (f(a)f(b)<0)
- wybierz punkt poczÄ…tkowy x0
- oblicz kolejne przybli\enia rozwiÄ…zania xn+1=F(xn, & ) dla n=0,1,2,&
- zakończ iteracje po spełnieniu warunku zatrzymania procedury
6
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Klasy metod iteracyjnych
Metody rekurencyjne obliczenia wykonywane w danym kroku procesu
iteracyjnego uwzględniają jedynie rezultaty uzyskane na poprzednim kroku
przykład: metoda stycznych (Newtona)
xn+1 = F(xn )
Metody nierekurencyjne obliczenia na danym etapie procedury iteracyjnej
biorÄ… pod uwagÄ™ wyniki uzyskane na wszystkich poprzednich etapach
xn+1 = F(x0, x1,& , xn-1, xn) przykład: metoda połowienia (bisekcji)
7
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Kryteria zatrzymania obliczeń
xn - xn-1 < µ1 f (xn ) < µ2
Rys. 2.3 Rys. 2.4
8
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Kryteria zatrzymania obliczeń
xn - xn-1 < µ1 f (xn ) < µ2
Rys. 2.5
Rys. 2.6
9
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Inne kryteria zatrzymania obliczeń
xn - xn-1
n
< µ3
xn
N numer iteraji
N numer iteraji
WzglÄ™dna rÌ\nica
WzglÄ™dna rÌ\nica
N max maksymalna
kolejnych przybli\eń
dopuszczalna liczba iteracji
10
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
RzÄ…d metody
RzÄ…d metody informuje o tym, jak szybko ciÄ…g {xn} definiowany przez
kolejne przybli\enia xn = F(xn-1) zbiega siÄ™ do rozwiÄ…zania x* (informuje o
tym, jak szybko maleje błąd przybli\eń w kolejnych iteracjach).
Metoda iteracyjna rozwiązywania równania f(x)=0 jest rzędu p (pe"1) je\eli
µn+1
µk = x - xk
µk = x* - xk
lim = C, C `" 0
lim = C, C `" 0
p
n "
µn
C - stała nieujemna (stała błędu asymptotycznego)
p - wykładnik zbie\ności metody
Im wy\szy rząd metody p tym szybsza zbie\ność do x*.
µk H"10-2,µk +1 H"10-3,µk +2 H"10-4,µk +3 H"10-5,...
p=1 metoda liniowa
1µk H" 10-2,µk +1 H" 10-4,µk +2 H" 10-8,µk +3 H" 10-16,...
p=2 metoda kwadratowa
11
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda bisekcji (połowienia)
Punkt początkowy x0 poło\ony jest na środku przedziału początkowego [a,b].
W ka\dej iteracji przedział przeszukiwania jest dzielony na połowy.
W kolejnym kroku rozpatrywana jest połowa zawierająca rozwiązanie (tj.
mająca ró\ne znaki funkcji na granicach przedziału). W n-tej iteracji przedział
poczÄ…tkowy jest podzielony przez 2n a rozwiÄ…zanie jest ograniczone przez
a - b
" =
"n =
2n
Mo\na z góry określić maksymalną liczbę iteracji k, niezbędną do uzyskania
xk - xk -1 < µ1
dokładności
ëÅ‚ öÅ‚
b - a b - a
d" µ1
k = log2 ìÅ‚
2k
µ1 ÷Å‚
íÅ‚ Å‚Å‚
k zaokrąglone do liczby całkowitej (od góry)
12
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda bisekcji - algorytm
1) Oblicz f(a) i f(b) na granicach przedziału [a, b] powinno być f(a)f(b) < 0
2) Oblicz przybli\enie rozwiÄ…zania x=(a+b)/2 ;
Sprawdz kryterium stopu
Je\eli kryterium stopu spełnione zatrzymaj obliczenia
Je\eli nie, przejdz do punktu 3.
3) Oblicz f(x)
Je\eli f(a)f(x) > 0 , to rozwiÄ…zanie znajduje siÄ™ w przedziale [x, b];
zawę\amy granice przyjmując dla nowego przedziału
a=x, f(a)=f(x) i kontynuujemy algorytm przechodzÄ…c do kroku 2)
Je\eli f(a)f(x) < 0 , to rozwiÄ…zanie znajduje siÄ™ w przedziale [a, x];
zawę\amy granice przyjmując dla nowego przedziału
b=x, f(b)=f(x) i kontynuujemy algorytm przechodzÄ…c do kroku 2)
Je\eli f(x) = 0, to rozwiązanie zostało znalezione - zatrzymaj obliczenia.
13
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda bisekcji interpretacja geometryczna
Rys. 2.7
14
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda bisekcji przypadek patologiczny
Rys. 2.8
15
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda bisekcji zbie\ność
Zatrzymanie obliczeń następuje gdy
- znaleziono rozwiązanie dokładne ;
- spełniono kryterium zatrzymania (dotyczące np. wielkości przedziału
ograniczającego rozwiązanie lub wartości funkcji f(x))
Zbie\ność metody
Metoda bisekcji jest zbie\na liniowo z wykładnikiem p=1 i stałą C=1/2.
Metoda bisekcji jest zbie\na liniowo z wykładnikiem p=1 i stałą C=1/2.
Zalety: zbie\ność zawsze zapewniona, mo\liwość oszacowania liczby
iteracji potrzebnych do wygenerowania rozwiązania o określonej
dokładności.
Wady: wolna zbie\ność
16
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda cięciw (Regula Falsi)
Nowe przybli\enie pierwiastka jest wyznaczone przez punkt xp, w którym prosta
łącząca punkty o współrzędnych (a, f(a)) i (b, f(b)) przecina oś x
xp - a b - xp
af (b) - bf (a)
= , xp =
0 - f (a) f (b) - 0 f (b) - f (a)
Je\eli f(x ) = 0 to pierwiastek został znaleziony.
Je\eli f(xp) = 0 to pierwiastek został znaleziony.
Je\eli nie, to zawę\amy przedział ograniczający rozwiązanie w zale\ności
od zmiany znaku funkcji na granicach przedziału, wybierając
przedział [xp,b] je\eli f(xp)f(b) < 0 lub przedział [a,xp] jeśli f(xp)f(b) > 0 .
Rozwiązanie poszukiwane jest iteracyjnie, w ka\dej kolejnej iteracji przedział zawierający
pierwiastek jest zawę\any w zale\ności od znaku f(xp)f(b).
17
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda cięciw - algorytm
1) Wybierz przedział [a,b] (w którym funkcja zmienia znak); oblicz ya=f(a) i yb=f(b).
af (b) - bf (a) bya - ayb b - a
2) Oblicz
x = = = a - ya
w celu zmniejszenia
f (b) - f (a) ya - yb yb - ya
błędów zaokrąglenia
sprawdz kryterium zatrzymania algorytmu
3) Je\eli f(x)f(b) > 0 , pierwiastek znajduje siÄ™ w przedziale [a,x];
3) Je\eli f(x)f(b) > 0 , pierwiastek znajduje siÄ™ w przedziale [a,x];
podstaw nową granicę b=x, yb=f(x) i wróć do kroku 2)
Je\eli f(x)f(b) < 0 pierwiastek znajduje siÄ™ w przedziale [x,b];
podstaw nową granicę a=x, ya=f(x) i wróć do kroku 2)
Je\eli f(x) = 0, znaleziono rozwiÄ…zanie, zatrzymaj obliczenia
Obliczenia mogą być zatrzymane po spełnieniu
podobnych kryteriów stopu jak w metodzie bisekcji.
18
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda cięciw interpretacja geometryczna
Rys. 2.9
19
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda cięciw zbie\ność
Zbie\ność metody
1+ 5
Metoda cięciw jest rzędu
p = H" 1.62
2
Zalety: pewność, szybsza zbie\ność ni\ metoda bisekcji.
Zalety: pewność, szybsza zbie\ność ni\ metoda bisekcji.
Wady: zbie\ność mo\e być wolna gdy f jest silnie wypukła lub wklęsła.
20
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda siecznych
Niech [a,b] będzie przedziałem zawierającym pierwiastek równania f(x)=0 a x0 wartością
rzędnej [a,b] dobraną tak aby f(x0)`"0 . Wybierzmy inny punkt startowy o rzędnej x1. Sieczna
Å‚Ä…czÄ…ca punkty (x0,f(x0)) i (x1,f(x1)) przecina oÅ› x w punkcie x2. Iteracyjnie kontynuujemy
procedurÄ™ przyjmujÄ…c punkty (x1,f(x1)) et (x2,f(x2)) w celu wyznaczenia kolejnego przybli\enia
rozwiązania o odciętej x2, itd. Otrzymujemy w ten sposób ciąg {xi} zdefiniowany przez
f (xn) (xn - xn-1)
xn+1 = xn - := F. )
xn+1 = xn - := F. )
(xn
(xn
f (xn ) - f (xn-1)
-
Zbie\ność jest zagwarantowana jedynie je\eli początkowa wartość x0 wybrana jest w bliskim
sąsiedztwie poszukiwanego pierwiastka który spełnia warunek .
2
F (x) < 1
21
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda siecznych - algorytm
Wybierz dwa punkty startowe x0 i x1 z przedziału [a,b]
Oblicz kolejne przybli\enie rozwiÄ…zania (n=1,2,& )
f (xn)(xn - xn-1)
xn+1 = xn - ; f (xn+1); n = n +1
f (xn) - f (xn-1)
Przerwij obliczenia je\eli
- znaleziono rozwiÄ…zanie problemu f(xn+1)=0
- lub spełniono warunek zatrzymania algorytmu, np.:
xn+1 - xn d" µ1
Uwaga : Początkowe wzory metody siecznych są podobne do metody cięciw ale sieczne
buduje siÄ™ bez sprawdzania warunku f(x1)f(x2) < 0 i wykorzystuje siÄ™ zawsze dwa
ostatnio wyznaczone punkty xn i xn-1.
22
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda siecznych interpretacja geometryczna
Rys. 2.10
23
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda stycznych (Newtona-Raphsona)
Załó\my, ze f(x) jest klasy C2 w otoczeniu pierwiastka x*.
Rozwinięcie f(x) w szereg Taylora wokół oszacowania xn wartości x* :
(x* - xn)2
f (x*) = f (xn) + (x* - xn ) f 2 (xn) + f 2 2 (xn )+&
2!
p
µn
Dla xn wystarczająco bliskiego x*, mo\na pominąć składniki rzędu (pe"0)
Poniewa\ x* jest pierwiastkiem f(x*) =0 , otrzymujemy przybli\enie
Poniewa\ x* jest pierwiastkiem f(x*) =0 , otrzymujemy przybli\enie
f (xn ) + (x* - xn ) f 2 (xn ) H" 0
µn
Przybli\enie błędu związanego z zastąpieniem x* przez xn
f (xn )
.
µn = x* - xn H" -
f 2 (xn )
xn+1 = xn + µn
f (xn )
Mo\na zatem zbudować
xn+1 = xn - = F(xn)
f 2 (xn )
następującą procedurę iteracyjną
24
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda stycznych - algorytm
x0 "[a,b]
1) Wybierz punkt startowy , oblicz f(x0), f '(x0), n=0
2) Oblicz przybli\enie rozwiÄ…zania xn+1
f (xn )
xn+1 = xn -
2
f (xn )
3) Je\eli kryterium zatrzymania nie jest spełnione, podstaw n=n+1
i powróć do kroku 2)
i powróć do kroku 2)
Przerwij obliczenia je\eli
- znaleziono rozwiÄ…zanie problemu f(xn+1)=0
- lub spełniono warunek zatrzymania algorytmu, np.
x1 - x0 d" µ1
n d" Nmax
W praktyce dodaje siÄ™ te\ warunek
25
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda stycznych interpretacja geometryczna
Rys. 2.11
26
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metoda stycznych - zbie\ność
Metoda nie zawsze jest zbie\na
f (x) "C2([a,b]) f (a) f (b) < 0
Je\eli
Oraz f (x) i f (x) nie zmieniajÄ… znaku w [a,b], to przyjmujÄ…c x0 tak \eby f (x0)f (x0)>0
zapewnimy zbie\ność metody z wykładnikiem p=2 (do pierwiastka pojedynczego)
Zalety: bardzo szybka zbie\ność (kwadratowa p=2)
Wady: trzeba obliczać pierwszą pochodną funkcji f(x) (znajomość analityczna
lub przybli\enie numeryczne) i dobrze wybrać punkt startowy x0.
27
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Uwagi
Metoda siecznych
f (xn )(xn - xn-1) f (xn )
.
xn+1 = xn - = xn -
f (xn) - f (xn-1)
f (xn) - f (xn-1)
(xn - xn-1)
Przybli\enie pochodnej
f (xn)
xn+1 = xn -
Metoda stycznych
2
2
f (xn)
f (xn)
Strategia rozwiązywania trudnych przypadków:
najpierw wolno zbie\na ale pewna metoda bisekcji, następnie szybka metoda siecznych
.
28
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metody szukania punktu stałego
Równanie f(x)=0 mo\e być zapisane w równowa\nej postaci
x=F(x),
gdzie F : R R jest pewnÄ… funkcjÄ… zmiennej rzeczywistej x.
Poszukiwanie x spełniającego f(x)=0 jest równowa\ne szukaniu takiego x, aby na
przykład spełnić równanie x = x +a f(x) (stała a `" 0).
algorytm
Mo\na zatem przyjąć F(x)= x +a f(x).
x1 = F ( x0 )
Pierwiastek x* równania x=F(x) nazywany jest
x2 = F ( x1 )
punktem stałym" funkcji F(x).
ï"
F(x) mo\e posiadać wiele punktów stałych.
xn = F ( xn -1 )
29
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metody szukania punktu stałego -przykład
Rys. 2.12
30
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metody szukania punktu stałego - zbie\ność
Warunek wystarczający zbie\ności
Je\eli F(x) : R R jest ciągła w przedziale [a,b], oraz
2
F (x) <1 "x "[a,b]
to metoda punktu stałego jest zbie\na.
31
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Metody przyspieszania zbie\ności
Przyspieszanie zbie\ności metod pierwszego rzędu polega na budowaniu na
podstawie ciÄ…gu wartoÅ›ci przybli\onych {xi} nowego ciÄ…gu wartoÅ›ci {x°i}, który
zbiega siÄ™ do pierwiastka x* szybciej ni\ ciÄ…g poczÄ…tkowy.
Metoda Aitken a
Metoda Aitken a
Metoda Steffensen a
Rozwiazywanie układu równań nieliniowych
- patrz koniec rozdziału " Układy równań liniowych"
32
M.Pyrz Metody numeryczne w mechanice Równania nieliniowe 10.2011
Wyszukiwarka
Podobne podstrony:
MNM 3 2014
MNM 4 2014
MNM 5 2014
MNM 1 2014
MNM mgr 2014 przyklad obliczeniowy nr 4
MNM pytania do Wykladu 2014
MNM mgr 2014 przyklad obliczeniowy do lab 1
próbna 29 marca 2014
Biuletyn 01 12 2014
Audyt wewnętrzny 2014 86 95
2014 grudziadz zestaw 1
Darr @ The Mall (2014)
kol zal sem2 EiT 13 2014
WYTYCZNE TCCC 2014 WERSJA POLSKA
2014 xv smp final wyniki
więcej podobnych podstron