|
|
Linia 1: |
Linia 1: |
| ''Uwaga: przekonwertowane latex2mediawiki;
| |
| prawdopodobnie trzeba wprowadzi� poprawki!''
| |
|
| |
|
| ==Wyznaczanie wektorów i wartości własnych==
| |
|
| |
| Niech będzie dana rzeczywista kwadratowa macierz <math>\displaystyle A</math> wymiaru <math>\displaystyle N</math>. Wektorem własnym <math>\displaystyle x\in C^N</math> oraz
| |
| odpowiadającą mu wartością własną <math>\displaystyle \lambda \in C</math> nazwiemy taką parę, dla której
| |
|
| |
| <center><math>\displaystyle A x = \lambda x,
| |
| </math></center>
| |
|
| |
| przy czym <math>\displaystyle x\neq 0</math>.
| |
|
| |
| Zadanie wyznaczania wartości własnych i wektorów własnych macierzy ma bardzo
| |
| szerokie zastosowania w tak odległych do siebie dziedzinach jak np. analiza
| |
| odporności konstrukcji mechanicznych (wieżowce, mosty, wagony kolejowe) na
| |
| wibracje, czy też rankingowanie stron internetowych w wyszukiwarce Google.
| |
|
| |
| {{przyklad|Odporność budynku na trzęsienie ziemi||
| |
|
| |
| Rozważmy prosty układ mechaniczny opisujący, naturalnie w pewnym
| |
| jedynie przybliżeniu, zachowanie się układu <math>\displaystyle N</math> ciężkich płyt połączonych ze
| |
| sobą relatywnie elatycznymi dźwigarami --- co może np. modelować konstrukcję
| |
| wieżowca.
| |
|
| |
| Wiadomo, że jeśli częstotliwości drgań własnych tego wieżowca będą bliskie
| |
| częstotliwości siły wymuszającej (o niewielkiej amplitudzie), to konstrukcja
| |
| wpadnie w rezonans i w końcu rozpadnie się wskutek zbyt wielkich przemieszczeń.
| |
| Wychylenia naszych płyt z położenia równowagi są opisywane układem pewnych
| |
| równań różniczkowych.
| |
| Teoria matematyczna takich równań różniczkowych pokazuje, że częstotliwości
| |
| drgań własnych to nic innego jak ''wartości własne'' pewnej
| |
| macierzy
| |
| wymiaru <math>\displaystyle 2N</math>,
| |
| która powstaje ze współczynników równania różniczkowego opisującego dynamikę
| |
| tego układu.
| |
| }}
| |
|
| |
| {{przyklad|Macierz Google'a||
| |
|
| |
| Podstawowy algorytm rankingowania stron WWW w [http://www.wikipedia.org/pagerank wyszukiwarce Google]
| |
| sprowadza się do znalezienia rzeczywistego ''wektora własnego'' <math>\displaystyle \pi</math> pewnej silnie
| |
| rozrzedzonej macierzy <math>\displaystyle A</math> (gigantycznego rozmiaru, równego liczbie indeksowanych
| |
| stron, czyli w chwili pisania tego tekstu około <math>\displaystyle 2.5\cdot 10^{10}</math> stron), odpowiadającego wartości własnej równej 1:
| |
|
| |
| <center><math>\displaystyle A \pi = \pi.
| |
| </math></center>
| |
|
| |
| Współrzędne wektora <math>\displaystyle \pi</math>
| |
| interpretuje się jako wartość rankingową kolejnych stron WWW. Aby wszystko miało
| |
| sens, współrzędne wektora muszą być z przedziału [0,1]. Pewne
| |
| twierdzenia matematyczne i subtelny dobór macierzy <math>\displaystyle A</math> gwarantują, że taki
| |
| wektor <math>\displaystyle \pi</math> zawsze istnieje i jest jedyny! Co więcej, wartość 1 jest
| |
| dominującą wartością własną <math>\displaystyle A</math>, a to z kolei ma ważne znaczenie dla tzw.
| |
| [[sec:metoda-potegowa|Uzupe�nij: metody potęgowej]] numerycznego wyznaczania takiego wektora.
| |
| }}
| |
|
| |
| {{przyklad|Wyznaczanie miejsc zerowych wielomianu||
| |
|
| |
| Jak wiadomo, wartości własne to miejsca zerowe wielomianu charakterystycznego
| |
| macierzy <math>\displaystyle P(\lambda) = \det(A - \lambda I)</math>. Zachodzi także fakt odwrotny, to
| |
| znaczy miejsca zerowe wielomianu są wartościami pewnej macierzy, np. miejsca
| |
| zerowe wielomianu
| |
|
| |
| <center><math>\displaystyle p(\lambda) = p_1 \lambda^N + \ldots + p_N \lambda + p_{N+1}
| |
| </math></center>
| |
|
| |
| są wartościami własnymi m.in. macierzy stowarzyszonej,
| |
|
| |
| <center><math>\displaystyle A = \begin{pmatrix}
| |
| -p_2/p_1 & -p_3/p_1 & \cdots & -p_{N+1}/p_1\\
| |
| 1 & & & \\
| |
| & 1 & & \\
| |
| & & \ddots & \\
| |
| & & & 1
| |
| \end{pmatrix}
| |
| </math></center>
| |
|
| |
| Funkcja Octave'a <code>compan(p)</code> wyznacza macierz stowarzyszoną dla zadanego
| |
| wielomianu o współczynnikach w wektorze <math>\displaystyle p = [p_1,\ldots,p_N, p_{N+1}]^T</math>. Z tej
| |
| macierzy korzysta następnie funkcja Octave'a <code>roots</code>, która właśnie w taki
| |
| sposób wyznacza pierwiastki wielomianów: jako wartości własne macierzy
| |
| stowarzyszonej.
| |
| }}
| |
|
| |
| {{przyklad|||
| |
| Praktyczne zadanie z macierzą symetryczną
| |
| }}
| |
|
| |
| W praktyce obliczeniowej spotyka się zazwyczaj kilka typów zagadnień:
| |
| * Wyznaczenie dominującej wartości własnej (to znaczy: największej co do
| |
| modułu) i odpowiadającego jej wektora własnego (a może kilku wektorów?)
| |
| * Wyznaczenie najmniejszej co do modułu wartości własnej i wektorów jej
| |
| odpowiadających (zauważmy, że to jest np. zadanie wyznaczenia ''jądra
| |
| macierzy osobliwej'' --- wtedy wiemy a priori, że szukana najmniejsza co do modułu
| |
| wartość własna to zero)
| |
| * Wyznaczenie wartości własnej najbliższej zadanej liczbie (to jest właśnie
| |
| odpowiedź na pytanie jak blisko częstości wymuszającej są częstości drgań
| |
| własnych budynku)
| |
| * Wyznaczenie wszystkich wartości własnych (na przykład, w celu znalezienia
| |
| wszystkich pierwiastków zadanego wielomianu)
| |
| * Wyznaczenie wszystkich wartości i wektorów własnych (tzw. pełne
| |
| zagadnienie własne)
| |
|
| |
| Jak domyślamy się, dla macierzy rozrzedzonych dużego wymiaru pełne zagadnienie
| |
| własne jest zbyt kosztowne, gdyż najczęściej macierz wektorów własnych --- nawet
| |
| dla macierzy rzadkiej --- jest gęsta.
| |
|
| |
| Ponieważ w zastosowaniach bardzo często pojawiają się macierze rzeczywiste
| |
| symetryczne (powyższe przykłady pokazują, że nie tylko!) szczegółową analizę
| |
| metod numerycznych ograniczymy do tego przypadku, gdyż wtedy zachodzi
| |
|
| |
| {{twierdzenie|o symetrycznym zadaniu włanym||
| |
|
| |
| Każda macierz rzeczywista symetryczna <math>\displaystyle A</math> wymiaru <math>\displaystyle N</math> ma rozkład
| |
|
| |
| <center><math>\displaystyle A = Q\Lambda Q^T,
| |
| </math></center>
| |
|
| |
| gdzie <math>\displaystyle Q\in R^{N\times N}</math> jest ortogonalna (tzn. <math>\displaystyle Q^TQ = I</math>), a jej kolumnami są
| |
| wektory własne <math>\displaystyle A</math>, natomiast <math>\displaystyle \Lambda\in
| |
| R^N</math> jest diagonalna z
| |
| wartościami własnymi <math>\displaystyle A</math> na diagonali:
| |
|
| |
| <center><math>\displaystyle \Lambda = \begin{pmatrix} \lambda_1 & & \\ & \ddots & \\ & &
| |
| \lambda_N\end{pmatrix} .
| |
| </math></center>
| |
|
| |
| }}
| |
|
| |
| ===Uwarunkowanie zadania===
| |
|
| |
| {{twierdzenie|Bauer-Fike||
| |
|
| |
| Niech <math>\displaystyle A\in R^{N\times N}</math> będzie diagonalizowalna, to
| |
| znaczy dla pewnej macierzy <math>\displaystyle X</math> zachodzi
| |
|
| |
| <center><math>\displaystyle X^{-1} A X = \begin{pmatrix} \lambda_1 & & \\ & \ddots & \\ & &
| |
| \lambda_N\end{pmatrix} ,
| |
| </math></center>
| |
|
| |
| a więc (gdyż macierz po prawej stronie jest podobna do <math>\displaystyle A</math>) <math>\displaystyle \lambda_i\in C</math>,
| |
| <math>\displaystyle i=1,\ldots,N</math> są
| |
| wartościami własnymi <math>\displaystyle A</math>. Rozważmy macierz zaburzoną <math>\displaystyle \widetilde{A}</math> i jakąś jej
| |
| wartość własną <math>\displaystyle \widetilde{\lambda}</math>. Wtedy istnieje wartość własna <math>\displaystyle \lambda_j</math>
| |
| macierzy <math>\displaystyle A</math> taka, że
| |
|
| |
| <center><math>\displaystyle |\lambda_j - \widetilde{\lambda}| \leq \mbox{cond} _2(X) ||A - \widetilde{A}||_2.
| |
| </math></center>
| |
|
| |
| }}
| |
|
| |
| Ponieważ dla rzeczywistej macierzy symetrycznej macierz przejścia <math>\displaystyle X</math> jest
| |
| ortogonalna,
| |
| <math>\displaystyle X^{-1} = X^T</math>, to mamy <math>\displaystyle \mbox{cond} _2(X) = 1</math> i w konsekwencji zachodzi
| |
|
| |
| {{wniosek|Wartości własne macierzy symetrycznej są doskonale uwarunkowane||
| |
|
| |
| Przy oznaczeniach jak [[thm:Bauer-Fike|Uzupe�nij: twierdzeniu Bauera-Fike'a]], jeśli
| |
| dodatkowo założymy, że macierz <math>\displaystyle A</math> jest rzeczywista i symetryczna, to
| |
|
| |
| <center><math>\displaystyle \min_{j=1,\ldots,N}|\lambda_j - \widetilde{\lambda}| \leq ||A - \widetilde{A}||_2.
| |
| </math></center>
| |
|
| |
| }}
| |
|
| |
| Z drugiej strony, dla macierzy niediagonalizowalnych, uwarunkowanie wartości
| |
| własnych może być
| |
| dowolnie duże, co ilustruje poniższy
| |
|
| |
| {{przyklad|||
| |
|
| |
| <center><math>\displaystyle A_\epsilon = \begin{pmatrix} a & 1 \\ \epsilon & a \end{pmatrix}
| |
| </math></center>
| |
|
| |
| Weźmy dla uproszczenia <math>\displaystyle a=0</math>.
| |
| Wartości własne <math>\displaystyle A_\epsilon</math> to zera wielomianu <math>\displaystyle p_\epsilon(\lambda) = \lambda^2 - \epsilon</math>,
| |
| zatem <math>\displaystyle \lambda_\epsilon = \pm \sqrt{\epsilon}</math> i w konsekwencji
| |
|
| |
| <center><math>\displaystyle |\lambda_\epsilon - \lambda_0| / ||A_\epsilon - A_0|| = \sqrt{\epsilon}/\epsilon
| |
| \rightarrow \infty,
| |
| </math></center>
| |
|
| |
| gdy <math>\displaystyle \epsilon \rightarrow 0^+</math>, a więc uwarunkowanie takiego zadania jest
| |
| nieskończone: dowolnie mała zmiana macierzy powoduje zaburzenie wartości
| |
| własnych niewspółmiernie wielkie wobec zaburzenia danych. Dodatkowo, wartości własne i wektory własne macierzy <math>\displaystyle A</math> dla
| |
| ujemnego parametru <math>\displaystyle \epsilon</math> są zespolone!
| |
|
| |
| [[Image:MNeigencond.png|frame|200px|Zachowanie się wartości własnych macierzy <math>\displaystyle A</math> (z
| |
| parametrem <math>\displaystyle a=1</math>) w otoczeniu <math>\displaystyle \delta = 0</math>]]
| |
|
| |
| }}
| |
|
| |
| Bardziej spektakularny przykład pochodzi od Wilkinsona:
| |
|
| |
| {{przyklad|Perfidny wielomian Wilkinsona||
| |
|
| |
| Niech
| |
|
| |
| <center><math>\displaystyle p(\lambda) = (\lambda -1)(\lambda - 2) \cdots (\lambda - 20).
| |
| </math></center>
| |
|
| |
| Zmiana współczynnika przy <math>\displaystyle \lambda^{19}</math> o <math>\displaystyle 10^{-7}</math> skutkuje presunięciem niektórych
| |
| miejsc zerowych nawet o kilka jednostek na płaszczyźnie zespolonej! Poniżej
| |
| pokazujemy to na numerycznym przykładzie, gdzie prócz w/w zaburzenia mamy
| |
| dodatkowo z zaburzeniami powstałymi wskutek wyznaczenia współczynników
| |
| wielomianu w arytmetyce zmiennoprzecinkowej.
| |
|
| |
| [[Image:MNwilkinson.png|frame|200px|Zera oryginalnego i lekko zaburzonego perfidnego wielomianu
| |
| Wilkinsona.]]
| |
|
| |
| Jak widzimy, zera bardzo mało zaburzonego wielomianu mogą stać się wyraźnie nie-rzeczywiste!
| |
|
| |
| }}
| |
|
| |
| Jeśli chodzi o wektory własne, ich wrażliwość na zaburzenia macierzy jest
| |
| bardziej skomplikowana i zależy m.in. od uwarunkowania wartości własnych (czego
| |
| łatwo się domyślić) oraz od tego, jak blisko siebie leżą wartości własne.
| |
|
| |
| ===Lokalizacja wartości własnych===
| |
|
| |
| Jak okaże się za chwilę, czasem warto mieć ogólne rozeznanie o tym, gdzie ''z
| |
| grubsza'' leżą wartości własne danej macierzy <math>\displaystyle A</math>. W tym celu mogą być nam
| |
| pomocne dwa fakty:
| |
|
| |
| {{fakt|||
| |
| Dowolna wartość własna <math>\displaystyle \lambda\in C</math> macierzy <math>\displaystyle A</math> spełnia
| |
|
| |
| <center><math>\displaystyle |\lambda| \leq ||A||,
| |
| </math></center>
| |
|
| |
| gdzie <math>\displaystyle ||A||</math> jest dowolną normą macierzową indukowaną przez normę wektorową.
| |
| }}
| |
|
| |
| Rzeczywiście, skoro istnieje wektor <math>\displaystyle x\neq 0</math> taki, że <math>\displaystyle Ax = \lambda x</math>, to stąd
| |
| <math>\displaystyle ||Ax||/||x|| = |\lambda|</math>, więc fakt powyższy wynika już z definicji normy
| |
| macierzy:
| |
|
| |
| <center><math>\displaystyle ||A|| = \max_{y\neq 0}\frac{||Ay||}{||y||} \geq ||Ax||/||x||.
| |
| </math></center>
| |
|
| |
| Drugie twierdzenie jest równie proste w dowodzie, ale daje trochę więcej
| |
| informacji o lokalizacji widma.
| |
|
| |
| {{twierdzenie|Gerszgorina||
| |
|
| |
| Wartości własne macierzy <math>\displaystyle A</math> leżą w sumie (teoriomnogościowej) dysków <math>\displaystyle K_i</math> na
| |
| płaszczyźnie zespolonej,
| |
|
| |
| <center><math>\displaystyle K_i = \{z \in C: |z - a_{ii}| \leq \sum_{j\neq i} |a_{ij}| \}, \qquad i =
| |
| 1,\ldots N.
| |
| </math></center>
| |
|
| |
| }}
| |
|
| |
| {{przyklad|Koła Gerszgorina||
| |
|
| |
| Niech
| |
|
| |
| <center><math>\displaystyle A = \begin{pmatrix}
| |
| 1.08930 & 1.38209 & -1.00037 & 0.69355 & 2.32178 \\
| |
| 0.14211 & 1.74696 & 1.68440 & 0.30664 & 1.26718 \\
| |
| -0.74620 & 2.02686 & -0.68293 & 0.19684 & 0.35854 \\
| |
| 0.83517 & 0.74987 & 1.71331 & 1.09765 & -0.44321 \\
| |
| 1.02132 & -2.62155 & 0.79247 & 1.11408 & 0.48076 \\
| |
| \end{pmatrix}
| |
| </math></center>
| |
|
| |
| [[Image:MNgershgorindisks.png|frame|200px|Lokalizacja wartości własnych macierzy <math>\displaystyle A</math> kołami Gerszgorina oraz zgrubna
| |
| lokalizacja wewnątrz okręgu
| |
| o promieniu równym <math>\displaystyle ||A||_1</math>. Dokładne wartości własne zaznaczone trójkącikami.]]
| |
|
| |
| }}
| |
|
| |
| ===Metoda potęgowa, odwrotna potęgowa, RQI===
| |
|
| |
| Jak wiemy z algebry, nawet gdy <math>\displaystyle A</math> jest macierzą rzeczywistą, jej
| |
| widmo może być zespolone! Analizując poniższe metody, będziemy zakładać, że poszykiwane wartości i wektory
| |
| własne <math>\displaystyle A</math> są
| |
| rzeczywiste. Iterując na liczbach rzeczywistych nie mamy wszak szansy, by
| |
| dotrzeć do liczb zespolonych!...
| |
|
| |
| ====Metoda potęgowa====
| |
|
| |
| Przypuśćmy, że wartości własne macierzy <math>\displaystyle A\in R^{N\times N}</math> spełniają
| |
|
| |
| <center><math>\displaystyle |\lambda_1| > |\lambda_2| \geq \ldots \geq |\lambda_N|,
| |
| </math></center>
| |
|
| |
| (to znaczy, istnieje dokładnie jedna ''dominująca'' wartość własna macierzy
| |
| <math>\displaystyle A</math>.
| |
|
| |
| Załóżmy także, że istnieje baza złożona z wektorów własnych <math>\displaystyle q_1,\ldots,q_N</math> tej
| |
| macierzy (tak jest np. dla macierzy symetrycznej na mocy
| |
| [[thm:symetric-eig|Uzupe�nij: twierdzenia o własnościach symetrycznego zadania
| |
| własnego]]).
| |
|
| |
| Kierunek własny <math>\displaystyle q_k</math> jakiejś macierzy <math>\displaystyle A</math> ma taką własność, że poddany działaniu przekształcenia
| |
| <math>\displaystyle A</math> wydłuża się <math>\displaystyle \lambda_k</math> razy, wobec tego, dowolny wektor <math>\displaystyle x\in R^N</math> poddany
| |
| działaniu <math>\displaystyle A</math> najbardziej wydłuży się w kierunku <math>\displaystyle q_1</math>. Iterując tę procedurę,
| |
| powinniśmy dostawać w wyniku wektory, w których coraz bardziej dominuje kierunek
| |
| <math>\displaystyle q_1</math>. Formalnie, niech
| |
|
| |
| <center><math>\displaystyle x = \alpha_1q_1 + \ldots + \alpha_Nq_N,
| |
| </math></center>
| |
|
| |
| wtedy
| |
|
| |
| <center><math>\displaystyle Ax = A \left( \sum_i \alpha_iq_i \right) = \sum_i \alpha_i A q_i
| |
| = \sum_i \alpha_i \lambda_i q_i
| |
| </math></center>
| |
|
| |
| i w konsekwencji
| |
|
| |
| <center><math>\displaystyle A^kx = \sum_i \alpha_i \lambda_i^k q_i = \lambda_1^k\left(\alpha_1q_1 +
| |
| \alpha_2\left(\frac{\lambda_2}{\lambda_1}\right)^kq_2 + \ldots +
| |
| \alpha_N\left(\frac{\lambda_N}{\lambda_1}\right)^kq_N \right).
| |
| </math></center>
| |
|
| |
| Ponieważ z założenia, że istnieje dokładnie jedna dominująca wartość własna,
| |
| <math>\displaystyle \left|\frac{\lambda_N}{\lambda_1}\right| < 1</math>, to wyrażenie w nawiasie dąży do
| |
| <math>\displaystyle \alpha_1q_1</math> i w konsekwencji wektory <math>\displaystyle x_k = A^kx</math> dążą, gdy
| |
| <math>\displaystyle k\rightarrow\infty</math>, do kierunku wektora własnego <math>\displaystyle q_1</math>, to znaczy wektora
| |
| odpowiadającego dominującej wartości własnej <math>\displaystyle A</math> (o ile tylko <math>\displaystyle \alpha_1
| |
| \neq 0</math>).
| |
|
| |
| Szybkość zbieżności metody potęgowej jest liniowa, o współczynniku zależnym od
| |
| stosunku <math>\displaystyle \lambda_2/\lambda_1|</math>. W patologicznym przypadku, gdy <math>\displaystyle |\lambda_1|
| |
| \approx |\lambda_2|</math>, może więc okazać się, że metoda praktycznie nie jest
| |
| zbieżna.
| |
|
| |
| W praktyce nie wyznaczamy wzorem <math>\displaystyle x_k = (A^k)\cdot x</math>, lecz raczej korzystamy z
| |
| metody iteracyjnej
| |
|
| |
| {{algorytm|Metoda potęgowa||
| |
| <pre>
| |
|
| |
| <math>\displaystyle x_0</math> <nowiki>=</nowiki> dowolny wektor startowy; k <nowiki>=</nowiki> 0;
| |
| while( !stop )
| |
| {
| |
| <math>\displaystyle y_k</math> <nowiki>=</nowiki> <math>\displaystyle Ax_{k-1}</math>;
| |
| <math>\displaystyle x_k</math> <nowiki>=</nowiki> <math>\displaystyle y_k/||y_k||_\infty</math>;
| |
| k++;
| |
| }
| |
| </pre>}}
| |
|
| |
| Warunek normowania ma m.in. na celu zapobieżenie powstawania nadmiaru i
| |
| niedomiaru (gdy <math>\displaystyle |\lambda_1| < 1</math>, to <math>\displaystyle ||A^kx|| \rightarrow 0</math>, a gdy
| |
| <math>\displaystyle |\lambda_1| > 1</math>, to <math>\displaystyle ||A^kx|| \rightarrow \infty</math>). Przy okazji,
| |
| <math>\displaystyle ||y_k||_\infty \rightarrow |\lambda_1|</math>, a więc mamy także sposób na wyznaczenie
| |
| przybliżenia dominującej wartości własnej.
| |
|
| |
| Zazwyczaj jako warunek stopu wybiera się kryterium małej poprawki, <math>\displaystyle ||x_k -
| |
| x_{k-1}|| \leq \epsilon</math>, lub warunek małego residuum, <math>\displaystyle ||Ax_k - \lambda_{1,k}
| |
| x_k||\leq \epsilon</math>, gdzie <math>\displaystyle \lambda_{1,k}</math> jest przybliżeniem <math>\displaystyle \lambda_1</math>
| |
| dostępnym na <math>\displaystyle k</math>-tej iteracji.
| |
|
| |
| [[Image:MNXXX.png|frame|200px|Zasada działania metody potęgowej]]
| |
|
| |
| Metoda potęgowa doskonale sprawdza się, gdy macierz <math>\displaystyle A</math> jest macierzą
| |
| rozrzedzoną --- np. w przypadku macierzy Google'a.
| |
|
| |
| ====Odwrotna metoda potęgowa====
| |
|
| |
| Zauważmy, że dla dowolnej macierzy kwadratowej <math>\displaystyle A</math> o wartościach własnych
| |
| <math>\displaystyle \lambda_k</math> i odpowiadających im wektorach własnych <math>\displaystyle q_k</math>, mamy:
| |
| * Macierz <math>\displaystyle A-\sigma I</math> ma wartości własne <math>\displaystyle \lambda_k - \sigma</math> oraz wektory
| |
| własne <math>\displaystyle q_k</math>,
| |
| * Jeśli dodatkowo <math>\displaystyle A</math> jest nieosobliwa, to macierz <math>\displaystyle A^{-1}</math> ma wartości
| |
| własne <math>\displaystyle 1/\lambda_k</math> oraz wektory własne <math>\displaystyle q_k</math>
| |
|
| |
| Łącząc te dwie własności mamy, że
| |
|
| |
| {{stwierdzenie|Transformacja widma macierzy||
| |
|
| |
| Macierz <math>\displaystyle (A-\sigma I)^{-1}</math> (o ile istnieje),
| |
| to ma wartości własne równe <math>\displaystyle \frac{1}{\lambda_k - \sigma}</math> i wektory własne
| |
| identyczne z <math>\displaystyle A</math>.
| |
| }}
| |
|
| |
| Skoro tak, to jeśli najbliższą <math>\displaystyle \sigma</math> wartością własną <math>\displaystyle A</math> jest <math>\displaystyle \lambda_j</math>,
| |
| wówczas metoda potęgowa zastosowana do macierzy <math>\displaystyle (A-\sigma I)^{-1}</math> zbiegnie do
| |
| <math>\displaystyle q_j</math>. To prowadzi do następującego algorytmu, odwrotnej metody potęgowej:
| |
|
| |
| {{algorytm|Odwrotna metoda potęgowa||
| |
| <pre>
| |
|
| |
| <math>\displaystyle x_0</math> <nowiki>=</nowiki> dowolny wektor startowy; k <nowiki>=</nowiki> 0;
| |
| while( !stop )
| |
| {
| |
| Rozwiąż układ równań <math>\displaystyle (A-\sigma I)y_k = x_{k-1}</math>;
| |
| <math>\displaystyle x_k</math> <nowiki>=</nowiki> <math>\displaystyle y_k/||y_k||_\infty</math>;
| |
| k++;
| |
| }
| |
| </pre>}}
| |
|
| |
| ====Metoda Rayleigh====
| |
|
| |
| Z własności metody potęgowej, metoda odwrotna potęgowa jest zbieżna tym
| |
| szybciej, im bliżej <math>\displaystyle \lambda_j</math> jest przesunięcie <math>\displaystyle \sigma</math> (w stosunku do
| |
| pozostałych wartości własnych). Dlatego dobrze byłoby --- dla zwiększenia
| |
| szybkości zbieżności iteracji --- poprawiać wartość przesunięcia <math>\displaystyle \sigma</math>,
| |
| korzystając z dotychczas wyznaczonego wektora <math>\displaystyle x_k \approx q_j</math> i ilorazu
| |
| Rayleigh:
| |
|
| |
| <center><math>\displaystyle \lambda_j = \frac{q_j^TAq_j}{q_j^Tq_j} \approx \frac{x_k^TAx_k}{x_k^Tx_k}
| |
| </math></center>
| |
|
| |
| {{algorytm|Metoda RQI (Rayleigh Quotient Iteration)||
| |
| <pre>
| |
|
| |
| <math>\displaystyle x_0</math> <nowiki>=</nowiki> dowolny wektor startowy; <math>\displaystyle \sigma_0</math> <nowiki>=</nowiki> przybliżenie <math>\displaystyle \lambda_j</math>; k <nowiki>=</nowiki> 0;
| |
| while( !stop )
| |
| {
| |
| Rozwiąż układ równań <math>\displaystyle (A-\sigma_k I)y_k = x_{k-1}</math>;
| |
| <math>\displaystyle x_k</math> <nowiki>=</nowiki> <math>\displaystyle y_k/||y_k||_2</math>;
| |
| <math>\displaystyle \sigma_{k+1}</math> <nowiki>=</nowiki> <math>\displaystyle x_k^TAx_k</math>;
| |
| k++;
| |
| }
| |
| </pre>}}
| |
|
| |
| (wybierając normowanie wektora <math>\displaystyle x</math> w normie euklidesowej upraszczamy co nieco
| |
| algorytm).
| |
|
| |
| Wielką zaletą metody RQI jest jej szybkość zbiezności: kwadratowa gdy wartość
| |
| własna jest pojedyncza, a nawet sześcienna w przypadku macierzy symetrycznej.
| |
|
| |
| Wadą metody RQI jest to, że na każdym jej kroku należy rozwiązywać układ równań
| |
| z ''inną'' macierzą.
| |
|
| |
| {{uwaga|Gdy złe uwarunkowanie pomaga...||
| |
|
| |
| Przez pewien czas numerycy odnosili się do tej metody z rezerwą,
| |
| twierdząc, i słusznie, że im lepszym przybliżeniem <math>\displaystyle q_j</math> będzie <math>\displaystyle \sigma_k</math>, tym
| |
| bardziej rośnie uwarunkowanie <math>\displaystyle A-\sigma_k I</math>, a tym samym --- błąd numerycznego
| |
| rozwiązywania układu z tą macierzą będzie coraz większy i metoda będzie tracić
| |
| stabilność. Tymczasem okazuje się, że --- choć rzeczywiście tak jest ---
| |
| wektor błędu ma kierunek praktycznie zgodny z kierunkiem poszukiwanego wektora
| |
| <math>\displaystyle q_j</math>, a tym samym tylko ''pomaga'' w zbieżności metody!
| |
| }}
| |
|
| |
| ===Metody rozwiązywania pełnego zadania własnego===
| |
|
| |
| Najszybszą obecnie znaną metodą rozwiązywania pełnego zadania własnego (to
| |
| znaczy znajdowania wszystkich wartości i wektorów własnych) macierzy ''symetrycznej'' jest metoda '''"dziel i rządź"'''.
| |
|
| |
| Dla macierzy niesymetrycznych, najbardziej dopracowanym i przetestowanym, a więc
| |
| godnym zaufania algorytmem jest '''metoda QR z przesunięciami'''
| |
| (wykorzystująca, jak łatwo się domyślić, rozkład QR macierzy). Metoda QR
| |
| przewyższa także metodę dziel i rządź w przypadku symetrycznym, gdy wymiar
| |
| macierzy jest mały (mniej więcej <math>\displaystyle N \leq 25</math>).
| |
|
| |
| Obie metody są oczywiście metodami iteracyjnymi, jednak przyjęło się nazywać je
| |
| metodami bezpośrednimi, gdyż praktycznie zawsze potrzebują z góry ograniczonej
| |
| liczby iteracji do tego, by zbiec do wyniku o (niemal) maksymalnej rozsądnej
| |
| dokładności.
| |
|
| |
| Dla efektywności obu metod kluczowy jest ''preprocessing'' macierzy,
| |
| pozwalający niezbyt wygórowanym kosztem <math>\displaystyle O(N^3)</math> operacji sprowadzić przez
| |
| ortogonalne podobieństwo zadanie z
| |
| macierzą gęstą <math>\displaystyle A</math> do zadania z macierzą Hessenberga (w przypadku
| |
| niesymetrycznym)
| |
|
| |
| <center><math>\displaystyle \begin{pmatrix}
| |
| * & * & * & * & \cdots & * \\
| |
| * & * & * & * & \cdots & * \\
| |
| & * & * & * & \cdots & * \\
| |
| & & \ddots & \ddots & \ldots & \vdots \\
| |
| & & & \ddots & \ddots & * \\
| |
| & & & & * & *
| |
| \end{pmatrix}
| |
| </math></center>
| |
|
| |
| bądź wręcz trójdiagonalną, gdy <math>\displaystyle A</math> była symetryczna.
| |
|
| |
| Każdą macierz kwadratową <math>\displaystyle A</math> da się sprowadzić do postaci Hessenberga sekwencją
| |
| przekształceń postaci
| |
|
| |
| <center><math>\displaystyle
| |
| A := Q_k A Q_k^T,
| |
| </math></center>
| |
|
| |
| gdzie <math>\displaystyle Q_k</math> jest pewnym przekształceniem Householdera. Niech <center><math>\displaystyle
| |
| A = \begin{pmatrix}
| |
| d_1 & * & * & * & \cdots & * \\
| |
| a_1 & * & * & * & \cdots & * \\
| |
| a_2 & * & * & * & \cdots & * \\
| |
| \vdots & * & \ddots & \ddots & \ldots & \vdots \\
| |
| \vdots & * & * & \ddots & \ddots & * \\
| |
| a_{N-1} & * & * & * & * & *
| |
| \end{pmatrix}
| |
| </math></center>
| |
|
| |
| i oznaczmy <math>\displaystyle a = [a_1,\ldots,a_{N-1}]^T</math>. Możemy
| |
| wziąć na początek przekształcenie Householdera <math>\displaystyle \widetilde{Q}_1</math> takie, że
| |
| <math>\displaystyle \widetilde{Q}_1a = c\cdot
| |
| e_1</math>, gdzie <math>\displaystyle e_1 = [1,0,\ldots, 0]^T</math>. Wtedy
| |
|
| |
| <center><math>\displaystyle
| |
| \begin{pmatrix}
| |
| 1 & \\
| |
| & \widetilde{Q}_1
| |
| \end{pmatrix}
| |
| \cdot A = \begin{pmatrix}
| |
| d_1 & * & * & * & \cdots & * \\
| |
| c & * & * & * & \cdots & * \\
| |
| & * & * & * & \cdots & * \\
| |
| & * & \ddots & \ddots & \ldots & \vdots \\
| |
| & * & * & \ddots & \ddots & * \\
| |
| & * & * & * & * & *
| |
| \end{pmatrix}
| |
| </math></center>
| |
|
| |
| To samo przekształcenie przyłożone z prawej strony zachowa pierwszą kolumnę i
| |
| w efekcie nie zmieni struktury macierzy:
| |
| <center><math>\displaystyle
| |
| \begin{pmatrix}
| |
| 1 & \\
| |
| & \widetilde{Q}_1
| |
| \end{pmatrix}
| |
| \cdot A
| |
| \cdot
| |
| \begin{pmatrix}
| |
| 1 & \\
| |
| & \widetilde{Q}_1
| |
| \end{pmatrix}
| |
| = \begin{pmatrix}
| |
| d_1 & * & * & * & \cdots & * \\
| |
| c & * & * & * & \cdots & * \\
| |
| & * & * & * & \cdots & * \\
| |
| & * & \ddots & \ddots & \ldots & \vdots \\
| |
| & * & * & \ddots & \ddots & * \\
| |
| & * & * & * & * & *
| |
| \end{pmatrix} .
| |
| </math></center>
| |
|
| |
| Dalej stosujemy tę samą metodę do podmacierzy wymiaru <math>\displaystyle N-1</math>, itd. aż dochodzimy
| |
| do macierzy Hessenberga.
| |
|
| |
| Gdy wyjściowa macierz <math>\displaystyle A</math> jest symetryczna, to z definicji, macierz wynikowa
| |
| <math>\displaystyle \begin{pmatrix}
| |
| I & \\
| |
| & \widetilde{Q}_{N-2}\end{pmatrix}
| |
| \cdots \begin{pmatrix}
| |
| 1 & \\
| |
| & \widetilde{Q}_1
| |
| \end{pmatrix}
| |
| A \begin{pmatrix}
| |
| 1 & \\
| |
| & \widetilde{Q}_1\end{pmatrix}
| |
| \cdots \begin{pmatrix}
| |
| I & \\
| |
| & \widetilde{Q}_{N-2}
| |
| \end{pmatrix} </math> też jest symetryczna i jednocześnie
| |
| Hessenberga --- a więc musi być trójdiagonalna! Ponadto, macierz wynikowa będzie
| |
| miała te same wartości własne, co <math>\displaystyle A</math>; wektory własne macierzy <math>\displaystyle A</math> także można
| |
| łatwo (jak?) odzyskać z wektorów własnych macierzy wynikowej.
| |
|
| |
| ====Metoda dziel i rządź====
| |
|
| |
| Jest to obecnie najefektywniejsza metoda rozwiązywania zagadnienia własnego
| |
| macierzy symetrycznej wymiaru powyżej kilkudziesięciu. Omówimy w zarysie jej
| |
| najprostszy wariant (obarczony pewnymi wadami, usuniętymi w wersji
| |
| bibliotecznej --- <code>DSYEVD</code> w LAPACKu).
| |
|
| |
| Startując z symetrycznej macierzy <math>\displaystyle A</math> już w postaci trójdiagonalnej, łatwo
| |
| widzieć, że "prawie" rozpada się ona na dwie mniejsze macierze trójdiagonalne:
| |
| dokładniej,
| |
|
| |
| <center><math>\displaystyle
| |
| \begin{pmatrix}
| |
| a_1 & b_1 & & \\
| |
| b_1 & a_2 & \ddots & \\
| |
| & \ddots & \ddots & b_{N-1} \\
| |
| & & b_{N-1} & a_N
| |
| \end{pmatrix}
| |
| =
| |
| \begin{pmatrix}
| |
| T_1 & \\
| |
| & T_2
| |
| \end{pmatrix}
| |
| + b_{m} uu^T,
| |
| </math></center>
| |
|
| |
| gdzie <math>\displaystyle T_1 = \begin{pmatrix}
| |
| a_1 & b_1 & & \\
| |
| b_1 & a_2 & \ddots & \\
| |
| & \ddots & \ddots & b_{m-1} \\
| |
| & & b_{m-1} & a_m - b_m
| |
| \end{pmatrix}
| |
| </math>, <math>\displaystyle T_2 = \begin{pmatrix}
| |
| a_{m+1} - b_m & b_{m+1} & & \\
| |
| b_{m+1} & a_{m+2} & \ddots & \\
| |
| & \ddots & \ddots & b_{N-1} \\
| |
| & & b_{N-1} & a_N
| |
| \end{pmatrix}
| |
| </math> są --- tak jak <math>\displaystyle A</math> --- macierzami
| |
| trójdiagonalnymi i symetrycznymi (jako podmacierze <math>\displaystyle A</math>), tylko o połowę mniejszego
| |
| wymiaru, gdy <math>\displaystyle m \approx N/2</math>. Natomiast <math>\displaystyle u = e_{m} + e_{m+1}</math>, tak więc macierz
| |
| <math>\displaystyle b_{m} uu^T</math> ma tylko cztery niezerowe elementy, każdy równy <math>\displaystyle b_m</math>.
| |
|
| |
| Zgodnie ze swoją nazwą, metoda dziel i rządź sprowadza zadanie znajdowania par
| |
| własnych macierzy wymiaru <math>\displaystyle N</math> do dwóch takich zadań dla macierzy dwa razy
| |
| mniejszych. Te z kolei można potraktować w taki sam sposób i iteracyjnie
| |
| zmniejszyć wymiar macierzy do tak małego (około 25), by opłacało się zastosować
| |
| metodę QR (teoretycznie, można byłoby oczywiście doprowadzić podział do momentu,
| |
| gdy macierze trójdiagonalne są rozmiaru <math>\displaystyle 1\times 1</math> --- dla których rozwiązanie
| |
| zadania włanego jest trywialne --- ale taki algorytm byłby bardziej kosztowny od
| |
| wariantu z udziałem QR).
| |
|
| |
| Rzeczywiście, przypuśćmy, że dla obu macierzy trójdiagonalnych <math>\displaystyle T_1,T_2</math> umiemy
| |
| rozwiązać zadanie własne tak, że znamy macierze: <math>\displaystyle Q_i</math> --- ortogonalną oraz
| |
| <math>\displaystyle D_i</math> --- diagonalną, takie, że
| |
| <center><math>\displaystyle
| |
| Q_i^T T_i Q_i = D_i \qquad i=1,2.
| |
| </math></center>
| |
|
| |
| Wtedy łatwo widzieć, że dla łatwo wyznaczalnego wektora <math>\displaystyle v</math>,
| |
|
| |
| <center><math>\displaystyle
| |
| \begin{pmatrix}
| |
| Q_1^T & \\
| |
| & Q_2^T
| |
| \end{pmatrix}
| |
| \left(
| |
| \begin{pmatrix}
| |
| T_1 & \\
| |
| & T_2
| |
| \end{pmatrix}
| |
| + b_{m} uu^T
| |
| \right)
| |
| \begin{pmatrix}
| |
| Q_1 & \\
| |
| & Q_2
| |
| \end{pmatrix}
| |
| =
| |
| \begin{pmatrix}
| |
| D_1 & \\
| |
| & D_2
| |
| \end{pmatrix}
| |
| + b_{m} vv^T.
| |
| </math></center>
| |
|
| |
| W ten sposób zadanie własne dla oryginalnej macierzy <math>\displaystyle T</math> wymiaru <math>\displaystyle N</math> jest
| |
| równoważne zadaniu własnemu macierzy diagonalnej zaburzonej o macierz rzędu 1.
| |
|
| |
| Na szczęście łatwo pokazać, że jeśli <math>\displaystyle \lambda</math> nie jest wartością własną
| |
| macierzy diagonalnej <math>\displaystyle D = \begin{pmatrix}
| |
| D_1 & \\
| |
| & D_2
| |
| \end{pmatrix} </math>, to wartości własne <math>\displaystyle \lambda</math> macierzy <math>\displaystyle D+ b_{m} vv^T</math>
| |
|
| |
| spełniają równanie
| |
|
| |
| <center><math>\displaystyle
| |
| f(\lambda) \equiv 1 + b_{m} \sum_{j=1}^N\frac{v_j^2}{d_j - \lambda} = 0,
| |
| </math></center>
| |
|
| |
| gdzie <math>\displaystyle d_j</math> są elementami na diagonali macierzy <math>\displaystyle D</math>.
| |
|
| |
| [[Image:MNsecular.png|frame|50|Wykres <math>\displaystyle f(\lambda)</math> dla macierzy jednowymiarowego
| |
| laplasjanu rozmiaru 10. Zwróć uwagę na asymptoty pionowe tej funkcji oraz jej
| |
| monotoniczność.]]
| |
|
| |
| W typowym przypadku <math>\displaystyle f</math> będzie miała dokładnie <math>\displaystyle N</math> pojedynczych miejsc zerowych
| |
| i wykres zachęcający do stosowania do niej metody Newtona. Okazuje się, że
| |
| ogólny przypadek nie jest istotnie trudniejszy, choć wymaga ważnych modyfikacji,
| |
| zarówno w celu szybszego rozwiązywania powyższego równania nieliniowego, jak i w
| |
| celu zapewnienia lepszej stabilności algorytmu.
| |
|
| |
| Ostateczny koszt wyznaczenia wszystkich wektorów i wartości własnych jest rzędu <math>\displaystyle O(N^3)</math> z
| |
| małą stałą.
| |
|
| |
| ====Metoda QR====
| |
|
| |
| Dla zadania własnego z macierzą niesymetryczną najczęściej stosuje się metodę
| |
| QR.
| |
|
| |
| Jakkolwiek ostateczna wersja metody <math>\displaystyle QR</math> działa dla macierzy niesymetrycznych,
| |
| wygodnie będzie nam założyć dla przejrzystości ekspozycji, że macierz jest
| |
| symetryczna i w konsekwencji ma rzeczywiste widmo.
| |
|
| |
| W najprostszym wariancie (bez przesunięć), algorytm QR ma postać:
| |
|
| |
| {{algorytm|Metoda QR||
| |
| <pre>
| |
|
| |
| <math>\displaystyle A_1 = A</math>;
| |
| for k <nowiki>=</nowiki> 1, 2, ...
| |
| {
| |
| wykonaj rozkład <math>\displaystyle A_k = Q_{k}R_{k}</math>;
| |
| <math>\displaystyle A_{k+1} = R_{k}\cdot Q_{k}</math>;
| |
| }
| |
| </pre>}}
| |
|
| |
| Można sprawdzić, że <math>\displaystyle A, A_1, A_2,\ldots</math> mają te same wartości własne, bo <math>\displaystyle A_{k+1} =
| |
| Q_{k+1}^TA_kQ_{k+1}</math>. Co więcej, powyższy algorytm (gdy <math>\displaystyle A</math> jest nieosobliwa) w
| |
| zasadzie jest równoważny teoretycznemu algorytmowi iteracji prostej
| |
| zastosowanemu nie do pojedynczego wektora, ale do <math>\displaystyle N</math> wektorów naraz:
| |
|
| |
| {{algorytm|Iteracja prosta na przestrzeni||
| |
| <pre>
| |
|
| |
| <math>\displaystyle V_1 = I</math>;
| |
| for k <nowiki>=</nowiki> 1, 2, ...
| |
| {
| |
| <math>\displaystyle W_{k+1} = A\cdot V_k</math>;
| |
| wyznacz rozkład QR <math>\displaystyle W_{k+1} = V_{k+1} R_{k+1}</math>, gdzie <math>\displaystyle V_{k+1}</math> jest ortogonalna;
| |
| }
| |
| </pre>}}
| |
|
| |
| Drugi krok w istocie ortogonalizuje kolumny <math>\displaystyle W_{k+1}</math>. Gdyby nie ortogonalizować
| |
| zestawu wektorów <math>\displaystyle W_{k+1}</math>, oczywiście dostalibyśmy w efekcie zbieżność
| |
| wszystkich kolumn macierzy do tego samego wektora --- odpowiadającego
| |
| dominującej wartości własnej <math>\displaystyle A</math>. Zapewniając sobie ortogonalność <math>\displaystyle V_{k+1}</math>,
| |
| możemy liczyć na to, że kolejne kolumny macierzy <math>\displaystyle V_k</math> będą dążyć do wektorów
| |
| własnych odpowiadających kolejnym wartościom własnym <math>\displaystyle A</math> (przy stosownych
| |
| założeniach o <math>\displaystyle A</math>, m.in. że
| |
| wszystkie wartości własne <math>\displaystyle A</math> spełniają <math>\displaystyle |\lambda_i| \neq |\lambda_j|</math> dla <math>\displaystyle i
| |
| \neq j</math>). Jeśli założyć dla uproszczenia, że oba używane rozkłady QR mają
| |
| jednoznacznie określone czynniki rozkładu (na przykład, wymuszając, by
| |
| diagonalne elementy macierzy <math>\displaystyle R</math> były dodatnie) mamy zależności <math>\displaystyle V_{k+1} =
| |
| Q_1\cdots Q_k</math> oraz <math>\displaystyle A_{k+1} = V_{k+1}^TAV_{k+1}</math>.
| |
|
| |
| Tak więc, w sprzyjających warunkach, metoda QR, jako równoważna
| |
| iteracji prostej na podprzestrzeni, będzie zbieżna: <math>\displaystyle A_k \rightarrow A_\infty</math>,
| |
| gdzie <math>\displaystyle A_\infty</math> jest macierzą trójkątną (bo wektory własne odpowiadające różnym
| |
| wartościom własnym są ortogonalne), a tym samym wartościami własnymi <math>\displaystyle A_\infty</math>
| |
| (a więc także <math>\displaystyle A</math>)
| |
| będą liczby na diagonali <math>\displaystyle A_\infty</math>.
| |
|
| |
| {{twierdzenie|Zbieżność metody QR w szczególnym przypadku||
| |
|
| |
| Niech wartości własne <math>\displaystyle A\in R^{N\times N}</math> spełniają <math>\displaystyle |\lambda_1|,\ldots,
| |
| |\lambda_N| > 0</math> oraz macierz <math>\displaystyle T = [x_1,\ldots,x_N]</math> o kolumnach <math>\displaystyle x_i</math> złożonych z
| |
| kolejnych wektorów własnych <math>\displaystyle A</math> ma taką własność, że <math>\displaystyle T^{-1}</math> ma rozkład LU,
| |
| <math>\displaystyle T^{-1} = LU</math>.
| |
|
| |
| Wtedy w metodzie QR, ciąg macierzy <math>\displaystyle Q_k</math> jest zbieżny do macierzy diagonalnej, a
| |
| ciąg <math>\displaystyle A_k</math> ma podciąg zbieżny do macierzy trójkątnej, której elementy diagonalne
| |
| <math>\displaystyle u_{ii}</math> są równe <math>\displaystyle \lambda_i</math> dla <math>\displaystyle i = 1,\ldots, N</math>.
| |
| }}
| |
|
| |
| Powyższa wersja algorytmu QR jest mało praktyczna, m.in. jest zbieżna wolno i
| |
| przy poważnych ograniczaniach na <math>\displaystyle A</math>. Sprytna modyfikacja algorytmu wyjściowego
| |
| daje w wyniku tzw. metodę QR z przesunięciami, która jest praktycznie
| |
| niezawodna dla dowolnej macierzy.
| |
|
| |
| {{algorytm|Metoda QR z przesunięciami||
| |
| <pre>
| |
|
| |
| <math>\displaystyle A_1 = A</math>;
| |
| for k <nowiki>=</nowiki> 1, 2, ...
| |
| {
| |
| wybierz sprytnie przesunięcie <math>\displaystyle \sigma_k</math>;
| |
| wykonaj rozkład <math>\displaystyle A_k - \sigma_kI = Q_{k}R_{k}</math>;
| |
| <math>\displaystyle A_{k+1} = R_{k}\cdot Q_{k} + \sigma_kI</math>;
| |
| }
| |
| </pre>}}
| |
|
| |
| Koszt wyznaczenia wszystkich wektorów i wartości własnych jest rzędu <math>\displaystyle O(N^3)</math> ze
| |
| stałą równą około 30.
| |
|
| |
| ===Biblioteki===
| |
|
| |
| LAPACK zawiera w sobie kolekcję doskonałych narzędzi do rozwiązywania różnych
| |
| wariantów zadania własnego, m.in. <code>DGEEV</code> dla macierzy niesymetrycznych
| |
| oraz <code>DSYEV</code> dla macierzy symetrycznych rozwiązują pełne zagadnienie własne
| |
| (wyznaczając wszystkie wartości własne i, opcjonalnie, wektory własne). Dla
| |
| macierzy symetrycznych mamy jeszcze m.in. funkcje <code>DSYEVX</code> (dla wybranych
| |
| wartości własnych) i <code>DSYEVD</code> (z algorytmem dzieli i rządź)
| |
|
| |
| Fortranowska biblioteka ARPACK rozwiązuje zadanie własne dla macierzy
| |
| rozrzedzonych, znajdując kilka wybranych (np. największych co do modułu)
| |
| wartości i wektorów własnych.
| |
|
| |
| Funkcja <code>eig</code> w Octave i MATLABie wyznacza wszystkie wartości własne (i
| |
| opcjonalnie wektory własne) zadaniej gęstej macierzy --- oczywiście korzystając
| |
| z LAPACKa. Jak dotąd, tylko MATLAB potrafi skorzystać z ARPACKa dla wyznaczenia
| |
| fragmentów widma macierzy rzadkiej, za pomocą funkcji <code>eigs</code>.
| |