MN13
From Studia Informatyczne
Spis treści |
Wektory i wartości własne
<<< Powrót do strony głównej przedmiotu Metody numeryczne
Niech będzie dana rzeczywista kwadratowa macierz wymiaru
. Wektorem własnym
oraz
odpowiadającą mu wartością własną
nazwiemy taką parę, dla której

przy czym .
Zadanie wyznaczania wartości własnych i wektorów własnych macierzy ma bardzo szerokie zastosowania w tak odległych od siebie dziedzinach, jak np. analiza odporności konstrukcji mechanicznych (wieżowce, mosty, wagony kolejowe) na wibracje, czy też rankingowanie stron internetowych w wyszukiwarce Google.
Przykład: Odporność budynku na trzęsienie ziemi
Rozważmy prosty układ mechaniczny opisujący, naturalnie w pewnym
jedynie przybliżeniu, zachowanie się układu 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), 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 symetrycznej
macierzy
wymiaru ,
która powstaje ze współczynników równania różniczkowego opisującego dynamikę
tego układu.
Przykład: Macierz Google'a
Podstawowy algorytm rankingowania stron WWW w wyszukiwarce Google
sprowadza się do znalezienia rzeczywistego wektora własnego pewnej silnie
rozrzedzonej macierzy
(gigantycznego rozmiaru, równego liczbie indeksowanych
stron, czyli w chwili pisania tego tekstu około
stron), odpowiadającego wartości własnej równej 1:

Współrzędne wektora
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
gwarantują, że taki
wektor
zawsze istnieje i jest jedyny! Co więcej, wartość 1 jest
dominującą wartością własną
, a to z kolei ma ważne znaczenie dla tzw.
metody potęgowej numerycznego wyznaczania takiego wektora.
Przykład: Wyznaczanie miejsc zerowych wielomianu
Jak wiadomo, wartości własne to miejsca zerowe wielomianu charakterystycznego
macierzy . Zachodzi także fakt odwrotny, to
znaczy miejsca zerowe wielomianu są wartościami pewnej macierzy, np. miejsca
zerowe wielomianu

są wartościami własnymi m.in. macierzy stowarzyszonej,

Funkcja Octave'a compan(p)
wyznacza macierz stowarzyszoną dla zadanego
wielomianu o współczynnikach w wektorze . Z tej
macierzy korzysta następnie funkcja Octave'a
roots
, która właśnie w taki
sposób wyznacza pierwiastki wielomianów: jako wartości własne macierzy
stowarzyszonej.
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, gdy wartość własna jest wielokrotna?)
- 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łasnym
Każda macierz rzeczywista symetryczna wymiaru
ma rozkład

gdzie jest ortogonalna (tzn.
), a jej kolumnami są
wektory własne
, natomiast
jest diagonalna z
wartościami własnymi
na diagonali:

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 . Przy tym mogą być nam pomocne dwa fakty:
Fakt
Dowolna wartość własna macierzy
spełnia

gdzie jest dowolną normą macierzową indukowaną przez normę wektorową.
Dowód
Rzeczywiście, skoro istnieje wektor taki, że
, to stąd
, więc fakt powyższy wynika już z definicji normy
macierzy:


Drugie twierdzenie jest równie proste w dowodzie, ale daje trochę więcej informacji o lokalizacji widma.
Twierdzenie Gerszgorina, o lokalizacji widma macierzy
Wartości własne macierzy leżą w sumie (teoriomnogościowej) dysków
na
płaszczyźnie zespolonej,

Przykład: Koła Gerszgorina
Wyznaczanie pojedynczej pary własnej
Jak wiemy z algebry, nawet gdy jest macierzą rzeczywistą, jej
widmo może być zespolone! Analizując poniższe metody, będziemy zakładać, że poszukiwane wartości i wektory własne
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 spełniają

(to znaczy, istnieje dokładnie jedna dominująca wartość własna macierzy
.
Załóżmy także, że istnieje baza złożona z wektorów własnych tej
macierzy (tak jest np. dla macierzy symetrycznej na mocy
twierdzenia o własnościach symetrycznego zadania własnego).
Kierunek własny jakiejś macierzy
ma taką własność, że poddany działaniu przekształcenia
wydłuża się
razy, wobec tego dowolny wektor
poddany
działaniu
najbardziej wydłuży się w kierunku
. Iterując tę procedurę,
powinniśmy dostawać w wyniku wektory, w których coraz bardziej dominuje kierunek
. Formalnie, niech

wtedy

i w konsekwencji

Założenia, że istnieje dokładnie jedna dominująca wartość własna,
, wynika, że wyrażenie w nawiasie dąży do
i w konsekwencji wektory
dążą, gdy
, do kierunku wektora własnego
, to znaczy wektora
odpowiadającego dominującej wartości własnej
(o ile tylko
).
Szybkość zbieżności metody potęgowej jest liniowa, o współczynniku zależnym od
stosunku . W patologicznym przypadku, gdy
, może więc okazać się, że metoda praktycznie nie jest
zbieżna.
W praktyce nie wyznaczamy wzorem , lecz korzystamy z
metody iteracyjnej
Algorytm Metoda potęgowa
<math>\displaystyle x_0</math> = dowolny wektor startowy; k = 0; while( !stop ) { <math>\displaystyle y_k</math> = <math>\displaystyle Ax_{k-1}</math>; <math>\displaystyle x_k</math> = <math>\displaystyle y_k/||y_k||_\infty</math>; k++; }
Warunek normowania ma m.in. na celu zapobieżenie powstawania nadmiaru i
niedomiaru (gdy , to
, a gdy
, to
). Przy okazji, zauważ, że
, 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, lub warunek małego residuum,
, gdzie
jest przybliżeniem
dostępnym na
-tej iteracji.




Metoda potęgowa doskonale sprawdza się, gdy macierz jest macierzą
rozrzedzoną --- np. w przypadku macierzy Google'a.
Odwrotna metoda potęgowa
Zauważmy, że dla dowolnej macierzy kwadratowej o wartościach własnych
i odpowiadających im wektorach własnych
, mamy:
- Macierz
ma wartości własne
oraz wektory własne
,
- Jeśli dodatkowo
jest nieosobliwa, to macierz
ma wartości własne
oraz wektory własne
Z połączenia tych dwóch własności wynika, że
Stwierdzenie O transformacji widma macierzy
Macierz (o ile istnieje),
to ma wartości własne równe
i wektory własne
identyczne z
.
Skoro tak, to jeśli najbliższą wartością własną
jest
,
wówczas metoda potęgowa zastosowana do macierzy
zbiegnie do
. To prowadzi do następującego algorytmu, odwrotnej metody potęgowej:
Algorytm Odwrotna metoda potęgowa
<math>\displaystyle x_0</math> = dowolny wektor startowy; k = 0; while( !stop ) { Rozwiąż układ równań <math>\displaystyle (A-\sigma I)y_k = x_{k-1}</math>; <math>\displaystyle x_k</math> = <math>\displaystyle y_k/||y_k||_\infty</math>; k++; }
Metoda Rayleigh (RQI)
Z własności metody potęgowej wynika, że metoda odwrotna potęgowa jest zbieżna tym szybciej, im bliżej jest przesunięcie
(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
,
korzystając z dotychczas wyznaczonego wektora
i ilorazu Rayleigh:

Stąd nazwa metody, w skrócie RQI (Rayleigh Quotient Iteration).
Algorytm Metoda RQI
<math>\displaystyle x_0</math> = dowolny wektor startowy; <math>\displaystyle \sigma_0</math> = przybliżenie <math>\displaystyle \lambda_j</math>; k = 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> = <math>\displaystyle y_k/||y_k||_2</math>; <math>\displaystyle \sigma_{k+1}</math> = <math>\displaystyle x_k^TAx_k</math>; k++; }
(wybierając normowanie wektora w normie euklidesowej upraszczamy co nieco
algorytm).
Wielką zaletą metody RQI jest jej szybkość zbieżnoś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 i skończona precyzja arytmetyki pomagają...
Przez pewien czas numerycy odnosili się do tej metody z rezerwą,
twierdząc, i słusznie, że im lepszym przybliżeniem będzie
, tym
bardziej rośnie uwarunkowanie
, 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 rozwiązanie układu jest obarczone wielkim błędem --- to wektor błędu ma kierunek praktycznie zgodny z kierunkiem poszukiwanego wektora
, a tym samym złe uwarunkowanie macierzy i skończona precyzja arytmetyki pomagają 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 ).
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 operacji sprowadzić przez
ortogonalne podobieństwo zadanie z
macierzą gęstą
(w przypadku
niesymetrycznym) do zadania z macierzą Hessenberga, czyli macierzą, której element
jest zerowy gdy tylko
:

bądź wręcz trójdiagonalną, gdy była symetryczna.
Każdą macierz kwadratową da się sprowadzić do postaci Hessenberga sekwencją
przekształceń postaci



i oznaczmy . Możemy
wziąć na początek przekształcenie Householdera
takie, że
, gdzie
. Wtedy

To samo przekształcenie przyłożone z prawej strony zachowa pierwszą kolumnę i w efekcie nie zmieni struktury macierzy:

Dalej stosujemy tę samą metodę do podmacierzy wymiaru , itd. aż dochodzimy
do macierzy Hessenberga.
Gdy wyjściowa macierz jest symetryczna, to z definicji, macierz wynikowa
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
; wektory własne macierzy
także można
łatwo (jak?) odzyskać z wektorów własnych macierzy wynikowej.
Metoda "dziel i rządź" dla macierzy symetrycznej
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 --- DSYEVD
w LAPACKu).
Startując z symetrycznej macierzy już w postaci trójdiagonalnej, łatwo
widzieć, że "prawie" rozpada się ona na dwie mniejsze macierze trójdiagonalne:
dokładniej,

gdzie ,
są --- tak jak
--- macierzami
trójdiagonalnymi i symetrycznymi (jako podmacierze
), tylko o połowę mniejszego
wymiaru, gdy
Natomiast
, więc macierz
ma tylko cztery niezerowe elementy, każdy równy
.
Zgodnie ze swoją nazwą, metoda "dziel i rządź" sprowadza zadanie znajdowania par własnych macierzy wymiaru 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
--- 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 umiemy
rozwiązać zadanie własne tak, że znamy macierze:
--- ortogonalną oraz
--- diagonalną, takie, że

Wtedy łatwo widzieć, że dla łatwo wyznaczalnego wektora ,

W ten sposób zadanie własne dla oryginalnej macierzy wymiaru
jest
równoważne zadaniu własnemu macierzy diagonalnej zaburzonej o macierz rzędu 1.
Na szczęście łatwo pokazać, że jeśli nie jest wartością własną
macierzy diagonalnej
, to wartości własne
macierzy
spełniają równanie

gdzie są elementami na diagonali macierzy
.
W typowym przypadku będzie miała dokładnie
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 z
małą stałą.
Metoda QR
Dla zadania własnego z macierzą niesymetryczną najczęściej stosuje się metodę QR.
Jakkolwiek ostateczna wersja metody QR 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
<math>\displaystyle A_1 = A</math>; for k = 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>; }
Można sprawdzić, że mają te same wartości własne, bo
. Co więcej, powyższy algorytm (gdy
jest nieosobliwa) w
zasadzie jest równoważny hipotetycznemu algorytmowi iteracji prostej
zastosowanemu nie do pojedynczego wektora, ale do
wektorów naraz:
Algorytm Metoda potęgowa na przestrzeni
<math>\displaystyle V_1 = I</math>; for k = 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; }
Drugi krok w istocie ortogonalizuje kolumny . Gdyby nie ortogonalizować
zestawu wektorów
, oczywiście dostalibyśmy w efekcie zbieżność
wszystkich kolumn macierzy do tego samego wektora --- odpowiadającego
dominującej wartości własnej
. Zapewniając sobie ortogonalność
,
możemy liczyć na to, że kolejne kolumny macierzy
będą dążyć do wektorów
własnych odpowiadających kolejnym wartościom własnym
(przy stosownych
założeniach o
, m.in. że
wszystkie wartości własne
spełniają
dla
). 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
były dodatnie) mamy zależności
oraz
.
Tak więc, w sprzyjających warunkach, metoda QR, jako równoważna
metodzie potęgowej na przestrzeni, będzie zbieżna: ,
gdzie
jest macierzą trójkątną (bo wektory własne odpowiadające różnym
wartościom własnym są ortogonalne. Tym samym, wartościami własnymi
(a więc także
)
będą liczby na diagonali
.
Twierdzenie O zbieżności metody QR
Niech wartości własne spełniają
oraz macierz
o kolumnach
złożonych z
kolejnych wektorów własnych
ma taką własność, że
ma rozkład LU,
.
Wtedy w metodzie QR ciąg macierzy jest zbieżny do macierzy diagonalnej, a
ciąg
ma podciąg zbieżny do macierzy trójkątnej, której elementy diagonalne
są równe
dla
.
Powyższa wersja algorytmu QR jest mało praktyczna, m.in. jest zbieżna wolno i
przy poważnych ograniczaniach na . 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
<math>\displaystyle A_1 = A</math>; for k = 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>; }
Koszt wyznaczenia wszystkich wektorów i wartości własnych jest rzędu ze
stałą równą około 30. Omówienie sposobów wyboru "sprytnego przesunięcia" wykracza niestety poza ramy wykładu.
Metoda Jacobiego dla macierzy symetrycznej
Na zakończenie wspomnijmy o (bardzo starej) metodzie Jacobiego, która działa na oryginalnej macierzy symetrycznej (bez konieczności uprzedniego sprowadzenia do postaci trójdiagonalnej).
Nie będziemy szczegółowo jej tu omawiać z braku miejsca, wymienimy tylko jej dwie najważniejsze cechy odróżniające ją od metod omawianych wcześniej:
- jest znacznie wolniejsza
- (niekiedy) wyznacza dokładniejsze wartości i wektory własne niż inne metody
Dodatkową zaletą metody Jacobiego jest to, że łatwo ją zrównoleglić (w czym jest podobna do metody dziel i rządź).
Pomysł metody Jacobiego jest stosunkowo prosty: należy sekwencją przekształceń ortogonalnych sprowadzić wyjściową macierz symetryczną
do postaci (prawie) diagonalnej:

Tak więc iteracja Jacobiego ma postać:
Algorytm Metoda Jacobiego
for k = 0, 1, ... wybierz <math>\displaystyle J_k</math>; <math>\displaystyle A</math> = <math>\displaystyle J_k^T A J_k</math>;
Macierze wybieramy jako tzw. obroty Jacobiego, tzn. obroty Givensa dobrane tak, by w danym kroku iteracji wyzerować kolejną parę pozadiagonalnych elementów macierzy. W klasycznej wersji, zerowaniu podlega pozadiagonalna para o największym module --- w ten sposób najbardziej zredukujemy miarę niediagonalności macierzy wyrażoną jako suma kwadratów elementów pozadiagonalnych:

Twierdzenie O zbieżności metody Jacobiego
Klasyczna metoda Jacobiego jest zbieżna co najmniej liniowo, tzn.

gdzie oznacza miarę niediagonalności na
-tym kroku iteracji.
Można pokazać, że asymptotycznie (tzn. dostatecznie blisko granicy) zbieżność metody Jacobiego jest nawet kwadratowa.
Uwarunkowanie
Twierdzenie Bauera-Fike'a, o uwarunkowaniu wartości własnych
Niech będzie diagonalizowalna, to
znaczy dla pewnej macierzy
zachodzi

a więc (gdyż macierz po prawej stronie jest podobna do )
,
są
wartościami własnymi
. Rozważmy macierz zaburzoną
i jakąś jej
wartość własną
. Wtedy istnieje wartość własna
macierzy
taka, że

Ponieważ dla rzeczywistej macierzy symetrycznej macierz przejścia jest
ortogonalna,
, to mamy
i w konsekwencji zachodzi
Wniosek Wartości własne macierzy symetrycznej są doskonale uwarunkowane
Przy oznaczeniach jak twierdzeniu Bauera-Fike'a, jeśli
dodatkowo założymy, że macierz jest rzeczywista i symetryczna, to

Z drugiej strony, dla macierzy niediagonalizowalnych, uwarunkowanie wartości własnych może być dowolnie duże, co ilustruje poniższy
Przykład

Weźmy dla uproszczenia .
Wartości własne
to zera wielomianu
,
zatem
i w konsekwencji

gdy , 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
dla
ujemnego parametru
są zespolone!
Bardziej spektakularny przykład pochodzi od Wilkinsona:
Przykład: Perfidny wielomian Wilkinsona
Niech

Zmiana współczynnika przy o
skutkuje przesunię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 wyżej wymienionego zaburzenia mamy dodatkowo do czynienia z zaburzeniami powstałymi wskutek wyznaczenia współczynników wielomianu w arytmetyce zmiennoprzecinkowej.
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.
Biblioteki
LAPACK zawiera w sobie kolekcję doskonałych narzędzi do rozwiązywania różnych
wariantów zadania własnego, m.in. DGEEV
dla macierzy niesymetrycznych
oraz DSYEV
dla macierzy symetrycznych rozwiązują pełne zagadnienie własne, wyznaczając wszystkie wartości własne i wektory własne. Dla
macierzy symetrycznych mamy jeszcze m.in. funkcje DSYEVX
(dla wybranych
wartości własnych) i DSYEVD
(z algorytmem "dziel 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 eig
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 eigs
.
Literatura
W celu dogłębnego zapoznania się z omawianym na wykładzie materiałem, przeczytaj rozdział 5.1, 5.2 i 5.5 w
- D. Kincaid, W. Cheney Analiza numeryczna, Wydawnictwa Naukowo-Techniczne, Warszawa 2006, ISBN 83-204-3078-X.
Część wykładu oparto na materiałach zawartych w bardzo ciekawym podręczniku
- J. Demmel, Numerical linear algebra, SIAM, 1997.
Od dziesięcioleci, wspaniałym przeżyciem jest lektura książki ojca nowoczesnej analizy numerycznej,
- J. H. Wilkinson, The algebraic eigenvalue problem, Clarendon Press, 1965,
a także
- B. Parlett, The symmetric eigenvalue problem, Prentice-Hall, 1980.