Metoda Laguerre’a

W analizie numerycznej metoda Laguerre’a jest algorytmem znajdowania pierwiastków dostosowanym do wielomianów. Innymi słowy, metoda Laguerre’a może być użyta do numerycznego rozwiązywania równania ${\displaystyle p(x)=0}$ dla danego wielomianu ${\displaystyle p(x).}$ Jedną z najbardziej użytecznych właściwości tej metody jest to (na podstawie obszernych empirycznych badań) że jest bardzo bliska bycia całkiem niezawodną metodą, co znaczy, że jest prawie zawsze gwarantowana zbieżność do „któregoś” pierwiastka wielomianu, bez znaczenia jak wybrany jest punkt początkowy. Ta metoda jest nazwana na cześć Edmonda Laguerre’a, francuskiego matematyka.

Definicja

Metoda Laguerre’a znalezienia jednego pierwiastka wielomianu ${\displaystyle p(x)}$ stopnia n jest następująca:

• Wybierz punkt startowy ${\displaystyle x_{0}.}$
• Dla ${\displaystyle k=0,1,2,\dots }$
  • Jeśli ${\displaystyle p(x_{k})}$ jest bardzo małe, przerwij pętlę.
  • Wylicz ${\displaystyle G={\frac {p'(x_{k})}{p(x_{k})}}.}$
  • Wylicz ${\displaystyle H=G^{2}-{\frac {p''(x_{k})}{p(x_{k})}}.}$
  • Wylicz ${\displaystyle a={\frac {n}{G\pm {\sqrt {(n-1)(nH-G^{2})}}}},}$ gdzie znak wybierany jest w ten sposób aby mianownik miał większy moduł dla uniknięcia utraty cyfr znaczących.
  • Przypisz ${\displaystyle x_{k+1}=x_{k}-a.}$
• Powtarzaj dotąd aż ${\displaystyle a}$ stanie się dostatecznie małe lub będzie osiągnięta maksymalna ilość iteracji.

Uwaga: algorytm należy zaimplementować na liczbach zespolonych, ponieważ nawet gdy ostatecznie zostanie znaleziony pierwiastek z zerową częścią urojoną, to wcześniej bardzo często wartość pod pierwiastkiem jest ujemna.

Jeśli zostanie znaleziony pierwiastek rzeczywisty, ${\displaystyle p}$ może być podzielone przez odpowiedni czynnik liniowy. Jeśli zostanie znaleziony pierwiastek zespolony ${\displaystyle a+bi}$ to gdy mamy jedynie rzeczywiste współczynniki, musi istnieć też pierwiastek sprzężony ${\displaystyle a-bi,}$ więc możemy podzielić przez ${\displaystyle (x-(a+bi))*(x-(a-bi))=x^{2}-2ax+a^{2}+b^{2}}$ (bez części urojonej). Ten krok deflacji redukuje stopień wielomianu o jeden lub o dwa, tak że ostatecznie możemy znaleźć przybliżenia wszystkich pierwiastków ${\displaystyle p.}$

Uwaga: deflacja może prowadzić do aproksymacji czynników różniących się wyraźnie od dokładnych wartości. Ten błąd jest najmniejszy jeśli pierwiastki zostają znajdowane w kolejności rosnącego modułu. Sposobem, który pozwala na uniknięcie utraty dokładności po kolejnych dzieleniach wielomianu jest kolejne wyszukanie pierwiastka zaczynając od znalezionej przybliżonej wartości, ale tym razem dla pełnego wielomianu^[1].

Wyprowadzenie

Zasadnicze twierdzenie algebry mówi że każdy wielomian ${\displaystyle p}$ stopnia ${\displaystyle n}$ może być napisany w formie

${\displaystyle p(x)=C(x-x_{1})(x-x_{2})\dots (x-x_{n}),}$

tak że ${\displaystyle x_{k}}$ gdzie ${\displaystyle (k=1,2,\dots ,n)}$ są pierwiastkami wielomianu. Jeśli weźmiemy logarytm naturalny obu stron, możemy znaleźć, że

${\displaystyle \ln |p(x)|=\ln |C|+\ln |x-x_{1}|+\ln |x-x_{2}|+\ldots +\ln |x-x_{n}|.}$

Oznaczmy pochodną przez

${\displaystyle G={\frac {d}{dx}}\ln |p(x)|={\frac {1}{x-x_{1}}}+{\frac {1}{x-x_{2}}}+\ldots +{\frac {1}{x-x_{n}}},}$

a zanegowaną drugą pochodną przez

${\displaystyle H=-{\frac {d^{2}}{dx^{2}}}\ln |p(x)|={\frac {1}{(x-x_{1})^{2}}}+{\frac {1}{(x-x_{2})^{2}}}+\ldots +{\frac {1}{(x-x_{n})^{2}}}.}$

Następnie zrobimy to co Acton zwał a ‘drastycznym zbiorem założeń’, że pierwiastek, który szukamy, powiedzmy ${\displaystyle x_{1}}$ jest w pewnej odległości od naszego przybliżenia ${\displaystyle x}$ i wszystkie inne pierwiastki są zgrupowane razem w pewnej odległości. Jeśli oznaczymy te odległości przez

${\displaystyle a=x-x_{1}}$

i

${\displaystyle b=x-x_{i},\quad i=2,3,\dots ,n}$

wtedy nasze równanie dla ${\displaystyle G}$ może być zapisane jako

${\displaystyle G={\frac {1}{a}}+{\frac {n-1}{b}}}$

i wyrażenie dla ${\displaystyle H}$ staje się

${\displaystyle H={\frac {1}{a^{2}}}+{\frac {n-1}{b^{2}}}.}$

Rozwiązując te równania dla ${\displaystyle a}$ znajdujemy że

${\displaystyle a={\frac {n}{G\pm {\sqrt {(n-1)(nH-G^{2})}}}},}$

gdzie pierwiastek kwadratowy liczby zespolonej jest wybrany taki, by powstał większy moduł mianownika, lub równoważnie, by spełnić:

${\displaystyle \operatorname {Re} \,({\overline {G}}{\sqrt {(n-1)(nH-G^{2})}}\,)>0,}$

gdzie Re oznacza rzeczywistą część liczby zespolonej, a ${\displaystyle {\overline {G}}}$ jest zespolonym sprzężeniem G; lub

${\displaystyle a={\frac {p(x)}{p'(x)}}\cdot \left({\frac {1}{n}}+{\frac {n-1}{n}}\,{\sqrt {1-{\frac {n}{n-1}}\,{\frac {p(x)p''(x)}{p'(x)^{2}}}}}\right)^{-1},}$

gdzie pierwiastek kwadratowy liczby zespolonej jest wybrany by mieć nieujemną część rzeczywistą.

Dla małych wartości ${\displaystyle p(x)}$ ten wzór różni of offsetu trzeciego rzędu metody Halleya błędem ${\displaystyle O(p(x)^{3})}$ tak że zbieżność do pierwiastka jest sześcienna.

Uwaga, nawet jeśli ‘drastyczny zbiór założeń’ nie działa dla pewnego szczególnego wielomianu ${\displaystyle p(x),}$ ${\displaystyle p(x)}$ może być transformowane do powiązanego wielomianu ${\displaystyle r}$ dla którego te założenia są prawidłowe, tj. przez przesunięcie początku w kierunku odpowiedniej liczby zespolonej ${\displaystyle q(x)=p(x-w),}$ dając odrębne pierwiastki odrębnych wielkości jeśli trzeba (co będzie jeśli niektóre pierwiastki są sprzężone), a następnie dostając ${\displaystyle r}$ z ${\displaystyle q(x)}$ przez powtarzające się stosowanie transformacji pierwiastkiem kwadratowym użytej w metodzie Graeffe’a wystarczająco razy by uczynić mniejsze pierwiastki znacząco mniejsze niż największy pierwiastek (i zgrupowane w porównaniu); przybliżony pierwiastek metody Graeffe’a może być następnie użyty jako start nowej iteracji metody Laguerre’a na ${\displaystyle r.}$ Przybliżony pierwiastek ${\displaystyle p(x)}$ może być następnie otrzymany prosto z tego dla ${\displaystyle r.}$

Jeśli zrobimy silniejsze założenie że człony w G odpowiadające pierwiastkom ${\displaystyle x_{i},\ i=2,3,\dots ,n}$ są pomijalnie małe w porównaniu to członu odpowiedzialnego za pierwiastek ${\displaystyle x_{1},}$ to prowadzi do metody Newtona.

Właściwości

Jeśli ${\displaystyle x}$ jest pojedynczym pierwiastkiem wielomianu ${\displaystyle p(x)}$ wtedy metoda Laguerre’a zbiega sześciennie niezależnie czy początkowe przybliżenie ${\displaystyle x_{0}}$ jest wystarczająco blisko pierwiastka ${\displaystyle x.}$ Jeśli natomiast ${\displaystyle x}$ jest pierwiastkiem wielokrotnym, wtedy zbieżność jest jedynie liniowa. To otrzymujemy wraz z karą obliczania wartości wielomianu wraz z jego pierwszą i drugą pochodną w każdej iteracji.

Dużą zaletą metody Laguerre’a jest to, że jest niemal zagwarantowana zbieżność do „pewnego” pierwiastka wielomianu bez względu jaki punkt początkowy został wybrany. Tak jest w odróżnieniu od innych metod takich jak metoda Newtona-Raphsona, które mogą skończyć się niepowodzeniem zbieżności dla źle wybranego początkowego przybliżenia. Może nawet zbiegać do zespolonego pierwiastka wielomianu ponieważ pierwiastek kwadratowy brany do obliczenia a powyżej może być liczbą ujemną. To może być uważane za zaletę lub wadę w zależności od zastosowania metody. Empirycznie pokazano, że brak zbieżności zdarza się ekstremalnie rzadko, czyniąc tę metodę dobrym kandydatem na algorytm znajdowania pierwiastków wielomianu ogólnego przeznaczenia. Jakkolwiek, ze względu na dość ograniczone teoretyczne zrozumienie algorytmu, wiele osób zajmujących się numeryczna analizą wahają się jej używać i preferują lepiej zrozumiałe metody jak algorytm Jenkinsa-Trauba, dla którego została opracowana solidniejsza teoria. Niemniej jednak, algorytm jest dość prosty do użycia w porównaniu do innych „murowanych” metod, łatwy dość do użycia ręcznie lub z pomocą kalkulatora kiedy komputer jest niedostępny. Prędkość zbieżności tej metody oznacza że bardzo rzadko potrzeba więcej niż kilku iteracji do otrzymania pełnej dokładności.

Kod

Metoda wymaga wyliczenia wartości wielomianu (najlepiej z użyciem schematu Hornera) dla liczb zespolonych choć w najczęstszym przypadku współczynniki wielomianu pozostają rzeczywiste. Poza tym wymaga liczenia pierwszej i drugiej pochodnej.

Wikibooks:pl:Kody źródłowe/Metoda Laguerre'a

O ile pierwiastek jednokrotny metoda osiąga szybko i z pełną dokładnością rzędu maszynowego Epsilon, to wylicza pierwiastek n-krotny jedynie z dokładnością ${\displaystyle {\sqrt[{n}]{epsilon}}}$ ponieważ w dość dużym pobliżu wartość wielomianu jest porównywalna z epsilon. Aby zwiększyć tę dokładność, należy, gdy wartość wielomianu będzie rzędu epsilon zacząć szukać pierwiastka pochodnej wielomianu i brać znaleziony punkt, jeśli nadal w tym punkcie wartość wielomianu jest wystarczająco mała. Dla więcej niż dwukrotnych pierwiastków procedura dla pochodnej jest kontynuowana, gdzie szuka się drugiej pochodnej itd.

Wikibooks:pl:Kody źródłowe/Metoda Laguerre'a dla pierwiastków wielokrotnych

Przypisy

Bibliografia

• Forman S. Acton: Numerical Methods that Work. Harper & Row, 1970. ISBN 0-88385-450-3. Brak numerów stron w książce
• S. Goedecker. Remark on Algorithms to Find Roots of Polynomials. „SIAM J. Sci. Comput.”. 15 (5), s. 1059–1063, 1994. DOI: 10.1137/0915064. 
• Wankere R. Mekwi: Iterative Methods for Roots of Polynomials. [w:] Master’s thesis, University of Oxford [on-line]. 2001. [dostęp 2016-11-04]. [zarchiwizowane z tego adresu (2012-12-23)].
• V.Y. Pan. Solving a Polynomial Equation: Some History and Recent Progress. „SIAM Rev.”. 39 (2), s. 187–220, 1997. DOI: 10.1137/S0036144595288554. 
• Section 9.5.3. Laguerre’s Method. W: WH Press, SA Teukolsky, WT Vetterling, BP Flannery: Numerical Recipes: The Art of Scientific Computing. Wyd. 3rd. New York: Cambridge University Press, 2007. ISBN 978-0-521-88068-8. Brak numerów stron w książce
• Anthony Ralston, Philip Rabinowitz: A First Course in Numerical Analysis. McGraw-Hill, 1978. ISBN 0-07-051158-6. Brak numerów stron w książce