super metoda globalna do równań nieliniowych!

Mam taki super pomysł na udoskonalenie metody Newtona:

chodzi o to że metoda Newtona jest najlepsza w zasadzie, ale wysiada niekiedy,
więc wystarczy naprawić ten motyw i będzie super - nie?

M. Newtona polega na stycznych - w postaci prostych,
a jak to proste one źle celują po prostu... haha!

Uwaga - korekta!

Zamiast głupich bo prostych stycznych, wystarczy użyć krzywych stycznych - proste?! :)

A jak wiadomo najprostsze krzywe są parabole - nie?

No to do dzieła, wystarczy to zastosować: jak wygląda metoda parabol stycznych?

Co to będzie parabola styczna i jak ją efektywnie liczyć?

Istnieje nieskończenie wiele paraboli stycznych do krzywej gładkiej w każdym jej punkcie (dla przykładu, dla krzywej y = 0 w punkcie [0; 0] (więc i dla wszystkich funkcji których pochodna w zerze jest zerem i które przechodzą przez [0; 0]) styczne będą wszystkie ax^2).

Chcesz wybierać jakąś konkretną? Jak?

EDYCJA:
Sytuacja wygląda tak: mamy funkcję f(x), która ma na przedziale [m; n] dokładnie jedno miejsce zerowe i jest ona różniczkowalna na całym tym przedziale.
Wyznaczamy sobie jakiś punkt [λ; f(λ)], gdzie λ ∈ [m; n], i rozpatrujemy parabole o równaniu p(x) = ax² + bx + c (co nam wystarczy, jak widać poniżej).

Wówczas styczna będzie rodzina parabol spełniająca układ równań:
aλ² + bλ + c = f(λ)
2aλ + b = f′(λ)

I co dalej? Te parabole mają naprawdę różne zachowania. Którą chcesz wybrać, co chcesz z nią zrobić? Tak jak w metodzie Newtona, szukać miejsc zerowych? Mogą być dowolne, zależne głównie od wyboru paraboli…

Przykład: f(x) = ln(x), na przedziale [½; 4]. Strzelamy λ := 2.
4a + 2b + c = ln(2)
4a + b = ½

Rozwiązanie to rodzina parabol z b - ½ - 4a; c = 4a - 1 + ln(2). Jak nam to pomaga? Co robimy dalej? W zależności od wyboru a, mogą one mieć miejsca zerowe gdzie chcą…

Jest na forum LaTeX, tak w ogóle? Kiedyś był, teraz nie umiem znaleźć…

Wiem że jesteście bardzo przewidywalni, ale aż takiego badziewia się nie spodziewałem!

Wyprowadzenie metody prostych stycznych, czyli Newtona: robimy... no jak to się robi?

Z wielomianu liniowego, zgadza się?

y = f'(x-x0) + f(a), co jest równaniem prostej - tak?
stąd Newton - metoda stycznych!

Zatem aby zrobić metodę stycznych, ale parabol wystarczy użyć... no, no... czego? paraboli oczywiście, haha!

y = f''/2 (x-a)^2 + f'(x-a) + f;
gdzie: f, f' i f'' wyznaczamy w x=x0 = a.

Wystarczy to rozwiązać, znaczy robimy tak: y = 0 i wyliczamy za pomocą wzorków z gimnazjum - pamiętacie: dela, b/2a, itd.?

no to do dzieła!

Mhm. To ten sam przykład co wyżej: f(x) = ln(x), a := 2:
y = ((-1/(x-2)^2)^2)/2 + 1/(x-2) + ln(2), ups… Nie ma miejsc zerowych…

No zarąbista metoda…

Dla porównania, metoda Newtona jest tutaj prostsza obliczeniowo, a jej trzecia iteracja daje wynik 0.9961317034475227, przy szóstej dającym wynik zaokrąglający się w Pythonie (w nim naklepałem na szybko) do 1.0.

To najpierw trzeba rozwiązać - wyprowadzić wzór parabol, a nie żeby od razu wyliczać!

y''/2 (x-a)^2 + y'(x-a) + y = 0; uwaga: y lepiej wygląda od f!

i rozwiązujemy to zwyczajnie: delta = y'^2 - 2yy'', itd.

x = a - [y' - sqrt(delta) / y'';

i to jest już gotowa super metoda parabol III-go rzędu!
Takie trudne było? hihi!

Poprawmy to troszkę, bo tu poprzez proste przekształcenia można wyeliminować problem gdy f'' = 0:

y'/y'' * [1 - sqr([1 - 2yy''/y'^2)] = y'/y'' (1 - 1 + 2yy''/y'^2)/(1 + sqrt(1 - 2yy''/y'^2)] =
2y/[y' * (1 + sqrt[1 - sqr([1 - 2yy''/y'^2)] = 2y / [y' + sqrt(y'^2 - 2yy'']

ostateczny wzór metody:
x = x - 2y / [y' + sqrt(y'^2 - 2yy'')]

Ładny?

Co jest od dawna znane - pod nazwą: metoda Laguerra... haha!

Równania kwadratowe wylicza w 1 kroku, np.: y = x^2 - a, co ma zero: x0 = sqrt(a);

y' = 2x, y'' = 2

podstawiając do metody:
x = x - 2x/2 [1 - sqrt(1 - 2(x^2-a)2/4x^2)] = x - x (1 - sqrt(a)/x) = x - x + sqrt(a) = sqrt(a);
perfekt!

Ten przykład:
x0 = 2; y = ln(x); y' = 1/x; y'' = -1/x^2

wstawiamy:
x - 2ln(x)/(1/x + sqrt(1/x^2 + 2ln(x)/x^2)] = x - 2xln(x) / [1 + sqrt(1 + 2lnx)]

x1 = 0.91047294161718621219608839957073
x2 = 1.0002981779500368739614442174409
x3 = 0.99999999999116891569423241675851

Co jest od dawna znane - pod nazwą: metoda Laguerra... haha!

Nadal nie rozumiem WTF. Siedzisz sobie z jakąs książka do metod numeryczych i przepisujesz nam ją na forum. Po co? Znowu chcesz trafić do perełek? o_O Przecież to jest żałosne co teraz robisz :D Jutro przeczytasz o kolejnej metodzie i znówu przyjdziesz na forum z niby pytaniem żeby potem pokazać jaki jesteś mądry, bo przed chwilą coś przeczytałeś? :D

Wyglądasz mniej więcej jak ten ziomek w długich włosach:

I wszystko fajnie, tylko że:
a) to nie jest metoda do rozwiązywania równań nieliniowych, tylko wielomianowych (rozszerzenie zakresu zastosowań pozostawiamy jako ćwiczenie dla czytelników)
b) potrzebuje do szczęścia operowania na C, i to bardzo
c) jej związek z parabolami, a w szczególności parabolami stycznymi jest… nader niewielki

Więc, no… Co następne? Napiszesz nam, że wzory Viety są wystarczające do znalezienia wszystkich pierwiastków wielomianu, mając na myśli metodę Łobaczewskiego? Przeczytasz coś innego, napiszesz dwa zdania o nader luźnym powiązaniu z tym, czym chcesz zaszpanować, i każesz zgadywać o co chodzi?

Jaki jest cel? Co chcesz zrobić? Pokazać że znasz metody numeryczne? Spoko, przez jakiś rok po tym, jak to miałem na studiach, to też to umiałem… A potem zapomniałem, bo mi nie było potrzebne, a jakby było, to wiem gdzie znaleźć…

Co teraz?
Wiadomo co: ta metoda jest jednak trochę za słaba,
zatem teraz trzeba zrobić metodę stycznych kubicznych! :)

A żeby było zabawniej: opartą na dwóch punktach i.. dwóch pochodnych.

x = g(y(a),y(b), y'(a),y'(b))

na pewno będzie rakieta! hohi!

Althorion:

1. metoda nie wymaga zespolonych, no, chyba że szukamy zer zespolonych
2. nie jest to metoda specjalnie do wielomianów, lecz dobra dla dowolnych równań
3. jest to metoda parabol stycznych, co można sobie łatwo narysować
4. jeśli ktoś nie lubi zbyt trudnych wzorów, może sobie wyeliminować ten pierwiastek:

sqrt(1 - 2yy''/y'^2) = 1 - yy''/y'^2 - ...

wtedy otrzymamy taki wzorek:
x = x - 2y/[y' (1 + 1 - yy''/y'^2] = x - y/y' / (1 - 1/2 yy''/y'^2)

co jest znowu od dawna znane - metoda Halleya... pewnie komety tym liczył! haha!
https://en.wikipedia.org/wiki/Halley%27s_method

Dobra.
Czas na moje wyniki, mam nadzieję że tym razem będzie to naprawdę nowa i super metoda.

interpolacja odwrotna na 4 punktach: y(a), y(b), y'(a) i y'(b)

tradycyjny trójkąt do interpolacji x(y):

a p = y(a) 
a p 1/p' 
b q a-b/p-q ((1/p'-(a-b)/(p-q)/(p-q)
b q 1/q'    ((a-b)/(p-q)-1/q'))/(p-q) (1/p'+1/q'-2(a-b)/(p-q))/(p-q)^2

p = y(a); q = y(b);

z tego od razu wyliczamy: x(y=0), nie potrzeba zapisywać wzoru, który wyglądałby dość strasznie, o tak:

x = a - p/p' + ((1/p'-(a-b)/(p-q))/(p-q) pq - (1/p'+1/q'-2(a-b)/(p-q))/(p-q)^2 pqq;

no to niech teraz ktoś sprawdzi czy i jak to działa... bo mi się nie chce. :)

Ale to chyba słabe będzie.

Po prostu nie można liczyć dobrze w takim scenariuszu:

mamy dwa punkty x=a i b, i teraz wyliczymy - nawet nieźle, kolejne przybliżenie,
tylko że następnym kroku to nam zdechnie, niestety:

ponieważ modyfikujemy tu tylko jeden punkt,
no a ten drugi pozostanie - słaby, więc osłabi nam całokształt! :)

Metoda Newtona przybliża, krzywą prostą.
Jest na to uogólnienie i nazywa się, szereg Taylora.
Newton to po prostu przybliżenie używające dwa pierwsze wyrazy szeregu Taylora.
Spokojnie można kombinować z przybliżeniem z kolejnym wyrazem, bo liczenie miejsca zerowego dla trójmianu jest dość proste.
Co jednak jak przybliżenie nie ma miejsc zerowych? Trójmian tak ma.

Jestem pewien, że już dawno temu, ktoś w ten sposób uogólnił metodę Newtona.
Pytanie, czy gra jest warta świeczki?

Gdy wielomian, czy nawet dowolna funkcja, nie ma zer rzeczywistych,
wtedy ta metoda z pierwiastkiem znajdzie zespolone zero.

np.:
 y = x^2 + 1; zera: x0 = i lub -i

x0 = 1; wtedy: y = 2; y' = 2; y'' = 2
i wyliczamy:

x = 1 - 2y/y' 1/[1 + sqrt(1 - 2yy''/y'^2)] = 1 - 2 /[1 + sqrt(1 - 8/4)] = 1 - 2/(1+i) = 1 -2(1-i)/2 = 1 -1 + i = i

proponuję sprawdzić takie coś:
y = exp(x^2) + 1; jakie to ma zera?

Uważam, że dobrze kombinujesz, tylko niestety w złym kierunku. Użycie aproksymacji funkcją kwadratowa a nie liniową podwyższa rząd metody. Jest to o tyle fajne, że w takiej sytuacji otrzymasz metodę zbiegającą znacznie szybciej do np. podwójnych miejsc zerowych, gdzie klasyczny Newton jest dość powolny, tj. zbiega liniowo. Nie ma żadnego powodu aby to nie działało. To, że @Althorion podał przypadek, że nie działa nie dyskwalifikuje metody, bo dla klasycznej metody Newtona też można podać takie przypadki.

Problem jednak polega na tym, że NIKT NIE POTRZEBUJE "szybszego Newtona". W większości praktycznych problemów miejsca zerowe są jednokrotne, a tam Newton zbiega niesamowicie szybko - typowo wystarczy kilka iteracji, a te iteracje są też dość proste, więc możliwe, że będzie i tak szybciej niż przy aproksymacji nieliniowej.

Natomiast głównym problemem metody Newtona jest to, że ona jest w ogólności zbieżna jedynie lokalnie. Jak masz wredną funkcję i wystartujesz ze złego punktu, to nie zbiegnie nigdy. Podnosząc rząd aproksymacji tak naprawdę pogarszasz obszar zbieżności czyli powiększasz główną wadę tej metody. Chodzi o to, że aproksymacje wielomianem wyższego rzędu zwykle są mniej dokładne daleko od punktu wokół którego aproksymujemy niż aproksymacje niższego rzędu np. liniowa.

Jak chcesz ulepszać metodę Newtona to poczytaj sobie o metodach bazujących na homotopiach. Tzn. tam stosuje się taki trick, że przekształca się parametrycznie oryginalną funkcję do funkcji łatwiejszej np. liniowej, gdzie rozwiązanie jest trywialne a potem modyfikując odpowiednio parametr stopniowo wraca się do funkcji oryginalnej śledząc miejsca zerowe. W praktyce kiedyś tego próbowałem, ale jest trudniejsze niż się wydaje - miejsca zerowe mogą się rozdwajać albo niestety... znikać i pojawiać w innym miejscu.

Nieprawda.
Metody wyższego rzędu mają szerszy obszar zbieżności.

Przykład - wcześniej obliczany: y = ln(x)

już dla x0 = e; metoda Newtona wysiada:
x = e - eln(e) = e - e = 0; śmierć!

Natomiast metoda 3-go rzędu (z interpolacji) zadziała dla e!
y' = 1/x = 1/e; y'' = -1/x^2 = -1/e^2

i wsadzamy to do tej metody Lagguerra:
x = e - 2elne/(1 + sqrt(1 - 2lne (-1/e^2)/(1/e)^2) = e - 2e / (1 + sqrt(3)) = 0.72836142073579679474666427173803

co jak widać nadal pięknie jedzie do 1 = wynik.

Z tego co widzę po tych improwizacjach, najlepsze metody są te, które wyliczamy bezpośrednio z interpolacji odwrotnej, czyli takie:

x = x0 + x' (y-y0) + x''/2 (y-y0)^2 + ...

w szczególności:

1. x - y/y' -> Newton: to jest liniowa interpolacja
2. x - y/y' - 1/2 y''/y'^3 y^2; to jest kwadratowa
3. x - y/y' - 1/2 y''/y'^3 y^2 - (y'''/y'^4 - 3y''^2/y'^5)*y^3; to jest kubiczna - rząd 4

itd.

kolejne są coraz lepsze, ale i wymagają coraz więcej obliczeń.

Metoda rzędu n -> oo byłaby po prostu dokładną funkcją odwrotną do y, więc dawałaby wynik w jednym kroku: x(y=0) = x0;

Dobrze, a jak to będzie działać dla przypadku wielowymiarowego? Bo przypadek jednowymiarowy jest w praktyce dość mało interesujący. Np. ile będziesz musiał się naliczyć aby rozwiązać tak układ 10000 równań? Bo zwykły Newton daje radę, ale i tak już wtedy potrzebujesz macierzy pochodnych 10000x10000 (na szczęście często może być rzadka), której obliczenie i rozkład zajmuje ponad 90% czasu metody.

Podejrzewam że aż tak wielkie układy równań nieliniowych są nadal poza zasięgiem komputerów.
Liniowe - tak.

Generalnie przy liczeniu układów równań niewiele się chyba zmienia (w tych formułach):
zamiast y i y' mamy wektor Y i macierz Y', no i tyle.

wyższe pochodne: Y'', Y''' ?
pewnie są zbyt ciężkie... no ale można to interpolować niższymi, np.: y'' =~ y'(a) - y'(b) / a-b;
więc podobnie i macierz Y'' można załatwić.