background image

Wyznaczanie wartości i wektorów własnych macierzy 

(problem własny)

Plan wykładu:

1. Pojęcia podstawowe, definicje

2. Metoda Kryłowa poszukiwania pierwiastków równania 

charakterystycznego

3. Lokalizacja (szacowanie) wartości własnych – tw. Gerschgorina
4. Proste metody iteracyjne 

a)

Metoda potęgowa dla wartości własnej o największym module

b)

przyśpieszanie zbieżności

c) redukcja macierzy metodą Hottelinga
d) redukcja macierzy metodą Wielandta

5. Sprowadzanie macierzy symetrycznych/hermitowskich do postaci 

trójdiagonalnej

a) metoda Hauseholdera

b) Metoda Lanczosa

6. Sprowadzanie macierzy symetrycznych/hermitowskich do postaci 

macierzy Hessenberga

7. Wyznaczanie wartości i wektorów własnych macierzy trójdiagonalnych 

i macierzy Hessenberga

a) Metoda bisekcji
b) Rozkład LR i QR, rozkład QR  metodą Hauseholdera  

     8.Uogólniony problem własny

background image

2

Pojęcia podstawowe

Często przy tworzeniu modeli matematycznych  
wykorzystywanych do symulacji zjawisk 
fizycznych czy zachowania się układu, zachodzi 
potrzeba rozwiązania tzw. 

problemu własnego  

(np. rów. Schrodingera):

- A jest macierzą kwadartową o nxn
- x

k

 jest 

wektorem własnym

 macierzy 

odpowiadającym 

wartości własnej

 

k

Nie zawsze układ równań, którego chcemy 
znaleźć rozwiązanie, przyjmuje tak prostą 
postać. Nierzadko mamy do czynienia z tzw. 

uogólnionym problemem własnym

Jeśli macierz B jest nieosobliwa to problem 
uogólniony można przekształcić 
do postaci:

Ax

k

= ¸

k

x

k

a

ij

; ¸

k

; x

(k)

m

2 C

A

¢ x = ¸B ¢ x

B

¡1

Ax = ¸x

Liczbę 

 nazywamy  wartością własną macierzy jeśli 

istnieje taki niezerowy wektor x dla którego zachodzi: 

Wektor x nazywamy (prawostronnym) wektorem 
własnym przynależnym do wartości własnej 

. Ciąg 

wszystkich wartości własnych nazywamy widmem 
macierzy 
A i oznaczamy: Sp(A).
Z powyższej definicji wynika: 

Macierz (A-

I) jest osobliwa, więc:

Wyznacznik ten jest wielomianem stopnia n zmiennej 
:

Każda wartość własna 

k

 jest pierwiastkiem 

wielomianu charakterystycznego

 macierzy A.

Ax = ¸x

(A

¡ ¸I)x = 0

det(A

¡ ¸I) = 0

w(¸)

=

det(A

¡ ¸I) =

=

(

¡1)

n

n

+ a

n

¡1

¸

n

¡1

+ : : : + a

0

)

A = [a

ii

]

background image

3

Def. Wartości i wektory własne macierzy 
transponowanej A

T

 nazywamy

 lewostronnymi 

wartościami

 i 

lewostronnymi wektorami 

własnymi

 macierzy A.

Wyznacznik macierzy po jej transponowaniu nie 
ulega zmianie. Dlatego widmo macierzy A jest 
równe widmu lewostronnemu. 

Tw. Jeżeli 

jest wartością własną macierzy, a 

jest jej lewostronną wartością własną
oraz gdy

Wówczas wektor własny x

l

 jest ortogonalny do 

lewostronnego wektora własnego x

p

.

Dla macierzy symetrycznej A=A

T

 wektory 

własne są zarazem wektorami 

lewostronnymi. 

Jeżeli więc wektory własne przynależą do różnych 
wartości własnych to są do siebie 

ortogonalne

.

¸

p

6= ¸

l

p

6= ¸

l

)

Def. Macierze A i B są podobne jeśli istnieje 
nieosobliwa macierz podobieństwa P, że:

Tw. Jeżeli macierze A i B są podobne to mają 
identyczne widmo wartości własnych.

Tw. Macierz Q

mxn

 (m

≥ n) nazywamy ortogonalną 

jeśli:

Tw. Jeżeli macierz Q

nxn

 jest ortogonalna to:

Tw. Macierz symetryczna A jest ortogonalnie 
podobna do macierzy diagonalnej D:

Tw. Wartości własne macierzy symetrycznej są 
rzeczywiste.

Def. Macierz o elementach zespolonych i własności

nazywamy macierzą hermitowską.
Wartości własne macierzy hermitowskiej są 
rzeczywiste a wektory własne są do siebie 
ortogonalne.

P

¡1

AP = B

Q

T

Q = I

n

£n

QQ

T

= I

n

£n

Q

T

AQ = D

x

x

x

T

l

x

x

x

p

= 0

x

x

x

T

l

Ax

x

x

p

=

¡

A

T

x

x

x

l

¢

T

x

x

x

p

= ¸

l

x

x

x

T

l

x

x

x

p

l

¡ ¸

p

)x

x

x

l

x

x

x

p

= 0

x

x

x

T

l

Ax

x

x

p

= ¸

p

x

x

x

T

l

x

x

x

p

A = (A

T

)

¤

) a

ii

2 R

background image

4

Tw. Dla dowolnej macierzy A istnieje macierz 
nieosobliwa P, która  może mieć elementy 
zespolone i zachodzi pomiędzy nimi związek

Powyższa macierz definiuje 

postać kanoniczną 

Jordana

J

k

 jest macierzą zdefiniowaną następująco

i

 jest wartością własną macierzy A i może wystąpić 

w wielu macierzach J

k

.

 

k = 1; 2; : : : ; K

· n

J

k

=

2
6

6

6

6

4

¸

i

1

0

: : :

0

0

0

¸

i

1

: : :

0

0

: : :

: : :

: : :

: : :

: : :

: : :

0

0

0

: : :

¸

i

1

0

0

0

: : :

0

¸

i

3
7

7

7

7

5

Wyznaczniki

są 

dzielnikami elementarnymi

 macierzy A. 

Liczba m jest stopniem macierzy J

k

. Dla m=1 

macierz J

k

 stanowi 

dzielnik liniowy

.

Jeśli wszystkie wartości własne macierzy A są różne 
to  wszystkie dzielniki elementarne są liniowe.

Iloczyn skalarny (wewnętrzny)

Iloczyn zewnętrzny

det(J

k

¡ ¸I) = (¸

i

¡ ¸)

m

P

¡1

AP =

2

6

6

4

J

1

0

: : :

0

0

J

2

: : :

0

: : :

: : :

: : :

: : :

0

0

: : :

J

k

3

7

7

5

haaa; bbbi = aaa

T

bbb = [®

1

; ®

2

; : : : ; ®

n

]

2
6

6

6

4

¯

1

¯

2

..

.

¯

n

3
7

7

7

5

a

a

a

ihbbb = aaabbb

T

=

2

6

6

6

4

®

1

®

2

..

.

®

n

3

7

7

7

5

1

; ¯

2

; : : : ; ¯

n

]

=

2
6

6

4

®

1

¯

1

®

1

¯

2

: : :

®

1

¯

n

®

2

¯

1

®

2

¯

2

: : :

®

2

¯

n

: : :

: : :

: : :

: : :

®

n

¯

1

®

n

¯

2

: : :

®

n

¯

n

3
7

7

5

background image

5

Tw. (Schura) Suma kwadratów modułów wartości 
własnych jest ograniczona od góry przez kwadrat 
normy euklidesowej:

Tw. Widmo macierzy ulega przesunięciu po 
dodaniu do niej macierzy jednostkowej 
pomnożonej przez liczbę:

widmo:

zostaje zastąpione przez:

Tw. (Cayleya-Hamiltona)
Jeśli

jest równaniem charakterystycznym macierzy A 
to

n

X

i=1

j¸j

2

· kAk

2

E

Sp(A + cI) = Sp(A) + c

Sp(A) =

1

; ¸

2

; : : :

g

Sp(A + cI) =

1

+ c; ¸

2

+ c; : : :

g

Metoda Kryłowa poszukiwania zer równania 
charakterystycznego

Korzystając z tw. CH 

Co dla dowolnego wektora y daje

układ n równań na n niewiadomych

Do jego utworzenia potrzeba jednak 
n

3

 obliczeń, oraz n

3

/3 aby go rozwiązać.

Uwaga:
Wyznaczenie zer wielomianu charakterystycznego 
może być źle uwarunkowane.

w(¸) = det(A

¡ ¸I) = 0

w(A) = 0

w(¸) = ¸

n

+

n

¡1

X

i=0

b

i

¸

i

= 0

w(A) = A

n

+

n

¡1

X

i=0

b

i

A

i

= 0

b

0

; b

1

; b

2

; : : : ; b

n

¡1

A

n

yyy +

n

¡1

X

i=0

b

i

A

i

yyy = 0

background image

6

Lokalizacja wartości własnych

Tw. (Gershgorina) Niech C

i

 oznaczają koła 

domknięte na płaszczyźnie zespolonej o 
środkach w punktach a

ii

 (elementy diagonalne 

macierzy A) i promieniach równych sumie 
modułów elementów z danego (i-tego) wiersza: 

Wówczas:

Jeżeli k kół C

i

 tworzy zbiór rozłączny z 

pozostałymi kołami, to w zbiorze tym leży 
dokładnie k wartości wlasnych macierzy A.

Wnioski:
1. Jeżeli macierz  jest symetryczna i diagonalnie 
dominująca o nieujemnych elementach na 
diagonali, to jest nieujemnie określona, a jeśli 
jest ona dodatkowo nieosobliwa to jest dodatnio 
określona. Macierz symetryczna silnie 
diagonalnie dominująca jest dodatnio określona 
wtedy i tylko wtedy gdy elementy na diagonali 
są dodatnie.
2. W każdym kole rozłącznym z pozostałymi 
leży dokładnie jedna wartość własna.

C

i

=

fz : jz ¡ a

ii

j ·

n

X

j=1

j

6=i

ja

ij

jg

Sp(A)

½

n

[

i=1

C

i

Wartości własne macierzy A leżą na płaszczyźnie 
zespolonej i zawarte są w kole o środku w 0 i promieniu 
równym promieniowi spektralnemu tej macierzy. 
Ponieważ:

więc można przyjać że:

Widma wartości własnych i lewostronnych wartości 
własnych są identyczne. 

Aby otrzymać lepsze 

oszacowanie położenia wartości własnych  można więc 
zastosować twierdzenie Geshgorina dla A

T

. Koła 

zawierające wartości własne mają środki w a

ii

 i promienie 

równe sumie modułów pozostałych elementów

w i-tych kulmnach.

Przykład. Podać lokalizację wartości własnych macierzy

½(A)

· kAk

p

;

p = 1; 2;

1; E

¸

i

2 fz : jzj · kAk

p

g

A =

2

4

¡2 ¡1 0

2

0

0

0

0

2

3

5

background image

7

a) najgorsze oszacowanie – lokalizacja w kole o promieniu 4
b) tw. Gershgorina – lokalizacja w kole o promieniu 3
c) tw. Gershgorina dla A

T

 – jedno z kół jest odseparowane (C'

3

) i zdegenerowane 

    znajduje się w nim dokładnie jedna wartość własna (

=2)  - najlepsze oszacowanie

background image

8

Metody wyznaczania wartości 

i wektorów własnych

macierzy

Proste metody 

iteracyjne

Metoda potęgowa 

I jej modyfikacje

Macierze o elementach

I wartościach 

własnych rzeczywistych

Macierze symetryczne

(w tym 

hermitowskie

)

Macierze niesymetryczne

Redukcja macierzy 

do postaci 

trójdiagonalnej

Redukcja macierzy do

postaci Hessenberga

Metoda 

Hauseholdera

Redukcja do postaci 

diagonalnej

metodą Jacobiego

Metoda

Givensa

Rozkład QR

Metoda Lanczosa

(

macierze rzadkie

)

Wartości i wektory 

własne

background image

9

Metoda potęgowa wyznaczania 
pojedynczych wartości własnych 
i wektorów własnych.

Załóżmy że istnieje n liniowo niezależnych 
wektorów własnych macierzy A, stanowią 
bazę przestrzeni liniowej

Wówczas dla dowolnego wektora v

0

 

Jeśli



i

  stanowią wartości własne macierzy 

Zakładamy że wartości własne tworzą ciąg

Jeśli 

1

 jest dominującą wartością własną, 

oraz wektor v

0

 ma składową w kierunku x

1

 to 

wówczas zachodzi

vvv

0

=

n

X

i=1

a

i

x

x

x

i

x

x

x

1

; x

x

x

2

; x

x

x

3

; : : : ; x

x

x

n

Avvv

0

=

n

X

i=1

a

i

¸

i

x

x

x

i

vvv

m

= A

m

vvv

0

=

n

X

i=1

a

i

¸

m

i

x

x

x

i

1

j ¸ j¸

2

j ¸ j¸

3

j ¸ : : : ¸ j¸

n

j

Z czego można wysnuć wniosek że wartość własną 
można obliczyć następująco

Dla dowolnego wektora y nieortogonalnego do x

1

.

Zazwyczaj y ma 1 na pozycji elementu o największym 
module w v

m+1

 a na pozostałych 0.

Jaka jest zbieżność metody? 

Zależy od (

i

/

1

)

 ale również od współczynników a

i

 

czyli od wyboru v

0

. Jeśli wartość własna o największym 

module jest zespolona to ciąg nie jest zbieżny.

Jak wyznaczyć wektor własny x

1

?

Ponieważ 

więc unormowany wektor własny będzie miał postać

lim

m

!1

A

m

vvv

0

¸

m

1

= a

1

x

x

x

1

¸

1

= lim

m

!1

yyy

T

vvv

m+1

yyy

T

vvv

m

x

x

x

1

=

vvv

m

jvvv

m

j

vvv

m

= ¸

m

1

"

a

1

x

x

x

1

+

n

X

i=2

µ

¸

i

¸

1

m

a

i

x

x

x

i

#

vvv

m

¼ ¸

m

1

a

1

x

x

x

1

background image

10

Jeśli wartość własna jest pierwiastkiem 
wielokrotnym rówanania charakterystycznego to 
metoda jest zbieżna bo składnik z 

1

 dominuje 

Uwaga: problem pojawia się gdy 

1

=-

j

 

 tj. 

identyczne moduły generują oscylacje (wtedy 
wybieramy ciąg wektorów v

2k

)

Przyśpieszanie zbieżności

1. Proces 

2

 

Po podzieleniu licznika i mianownika przez b

1

1

m

 a 

następnie rozwinięciu M w szereg potęgowy 
dostajemy

Co dla dużych m można przybliżyć

vvv

m

= A

m

vvv

0

= ¸

m

1

k

X

i=1

a

i

x

x

x

i

+

n

X

i=k+1

¸

m

i

a

i

x

x

x

i

r =

¸

2

¸

1

¸

1

¼ R

m+2

¡

(¢R

m+1

)

2

¢

2

R

m

Wykorzystując różnice progresywne 1 i 2 rzędu można 
znaleźć lepsze przybliżenie

2. Iloraz Rayleigha

Jeśli macierz jest symetryczna to jej wektory własne 
są ortogonalne. Załóżmy że są również ortonormalne 

Wtedy

Oraz

co daje lepszą zbieżność niż wariant podstawowy 
metody

vvv

T

m

Avvv

m

= vvv

T

m

vvv

m+1

=

n

X

i=1

a

2

i

¸

2m+1

i

x

x

x

T

i

x

x

x

j

= ±

ij

vvv

T

m

vvv

m

=

n

X

i=1

a

2

i

¸

2m

i

vvv

T

m

vvv

m+1

vvv

T

m

vvv

m

= ¸

1

+ O

¸

i

¸

1

2m

#

¸

1

¼ R

m

¡ ¯

2

r

m

eee

T

j

vvv

m

=

n

X

i=1

a

i

¸

m

i

eee

T

j

x

x

x

i

=

n

X

i=1

b

i

¸

m

i

eee

T

j

vvv

m+1

eee

T

j

vvv

m

=

b

1

¸

m+1

1

+

P

n
i=2

b

i

¸

m+1
i

b

1

¸

m

1

+

P

n
i=2

b

i

¸

m

i

eee

T

j

vvv

m+1

eee

T

j

vvv

m

= ¸

1

+ ¯

2

µ

¸

2

¸

1

m

+ ¯

3

µ

¸

3

¸

1

m

+ : : :

R

m

=

eee

T

j

vvv

m+1

eee

T

j

vvv

m

background image

11

Wyznaczanie pozostałych wartości własnych

a) metoda redukcji wektora

Jeśli znamy wartość własną o największym 
module to możemy wykorzystać ten fakt przy 
wyznaczaniu kolejnej największej co do modułu 
wartości własnej tj. 

2

Ponieważ wektory v

m+1

 oraz 

1

v

m

 są bliskie więc 

metoda może być w pewnych przypadkach 
nieużyteczna.  

b) metoda zerowania składowej 

Znając 

1

 możemy zdefiniować wektor

Czyli z v

0

 usuwamy składową w kierunku x

1

Wektor w

0

 nie ma składowej w kierunku x

1

 więc 

ciąg w

m

  powinien być zbieżny do  

2

 oraz x

2

. Ze 

względu na błędy zaokrągleń jednak usuwa się co 
pewną ilość iteracji składową x

1

 tj.

w

w

w

m

= vvv

m+1

¡ ¸

1

vvv

m

=

n

X

i=2

a

i

i

¡ ¸

1

m

i

x

x

x

i

w

w

w

0

= (A

¡ ¸

1

I)vvv

0

zzz

m

= (A

¡ ¸

1

I)w

w

w

m

Aby wyznaczyć kolejne wartości własne korzystamy z 
wektora

Taki proces staje się mało wydajny dla kolejnych 
wartości własnych.

Uwaga:

Metoda redukcji składowej jest odpowiednikiem  

metody metody czasu urojonego z 
ortogonalizacją Schmidta

 stosowaną w mechanice 

kwantowej do rozwiązania równ. Schrodingera. 

Zaletami są: prostota oraz fakt że nie trzeba 
pamietać postaci macierzy (elementy mogą być 
wyznaczane na bieżąco).

w

w

w

0

= (A

¡ ¸

1

I)(A

¡ ¸

2

I)vvv

0

background image

12

Redukcja macierzy

Tw. Jeśli 

1

 jest wartością własną macierzy A i  x

1

 

odpowiadającym jej wektorem własnym oraz dla 
dowolnego wektora v o własności

Macierz zredukowana 

Ma te same wartości co macierz A oprócz 

1

 która 

jest zerem.

Pełny dowód – A. Ralston „Wstęp do analizy...”

Zakładamy 

Definiujemy wektor u

(e

i

 - i -ty wersor w n-wymiarowej przestrzeni 

kartezjańskiej, lub i-ta kolumna macierzy 
jednostkowej)
oraz macierz T

W

1

= A

¡ ¸

1

x

x

x

1

vvv

T

x

(1)
1

= 1

u

u

u = x

x

x

1

¡ eee

1

T = I + u

u

ueee

T

1

T T

¡1

= I

) T

¡1

= I

¡ uuueee

T

1

Przy pomocy T konstruujemy macierz podobną do A

pierwsza kolumna C jest równa
(C

10

 – pierwszy wiersz macierzy C)

Ponadto definiujemy macierz D

w której wiersze od 2 do n są równe 0, natomiast 
wiersz pierwszy ma postać

Wobec powyższego można rozpatrywać macierz 

która ma te same wartości własne co A za wyjątkiem 

1

 która jest równa 0.

C = T

¡1

AT

C

01

= ¸

1

eee

1

D = T

¡1

x

x

x

1

vvv

T

T

vvv

T

x

x

x

1

= 1

C

¡ ¸

1

D = T

¡1

(A

¡ ¸

1

x

x

x

1

vvv

T

)T = T

¡1

W

1

T

D

10

= vvv

T

+ (1

¡ v

(1)

)eee

T

1

background image

13

Związek wektorów własnych macierzy W i A 

Algorytm wykorzystujący redukcję macierzy

1. wyznaczamy metodą potęgową 

1

 oraz x

1

2. konstruujemy macierz W

1

 i znajdujemy 

2

 

oraz w

2

 

3. konstruujemy macierz

i szukamy  

3

 oraz w

3

Uwaga: 
jeśli wartosć własna jest k-krotna to 
zastosowanie powyższej procedury k-krotnie da 
tę samą wartość własną oraz k różnych 
wektorów.

x

x

x

i

= (¸

i

¡ ¸

1

)w

w

w

i

+ ¸

1

(vvv

T

w

w

w

i

)x

x

x

1

¸

i

6= ¸

1

x

x

x

i

= w

w

w

i

¸

i

= ¸

1

;

(A

¡ ¸

1

I)w

w

w

i

= 0

W

2

= W

1

¡ ¸

2

w

w

w

2

vvv

T

1

vvv

T

1

w

w

w

2

= 1

Metoda redukcji Hotellinga

Za wektor v przyjmujemy lewy wektor własny 
przynależny do wartości własnej 

1

. Ale na ogół nie 

znamy lewych wektorów.

Metoda jest więc skutecza tylko w  przypadku macierzy 
symetrycznych, wtedy lewe wektory są identyczne z 
prawymi 

lub rekurencyjnie

vvv = x

x

x

1

W

1

= A

¡ ¸

1

x

x

x

1

x

x

x

T

1

W

0

= A

W

i

= W

i

¡1

¡ ¸

1

x

x

x

i

x

x

x

T

i

i = 1; 2; : : : ; n

¡ 1

background image

14

Metoda redukcji Wielandta

Wektor v definiujemy następująco

Gdzie: x

1

(j)

 jest j-tą składową wektora x

1

A

j0

 jest j-tym wierszem macierzy A

Z powyższego wzoru wynika że j-ty wiersz macierzy 
W

1

 jest równy

vvv

T

=

A

j0

¸

1

x

(j)
1

vvv

T

x

x

x

1

=

A

j0

x

x

x

1

¸

1

x

(j)
1

=

¸

1

x

(j)
1

¸

1

x

(j)
1

= 1

W

1

= A

¡ ¸

1

x

x

x

1

vvv

T

= A

¡

x

x

x

1

A

j0

x

(j)
1

(W

1

)

j0

= A

j0

¡

x

(j)
1

A

j0

x

(j)
1

= 0

Wnioski:

a) wektory własne odpowiadające niezerowym 

wartościom własnym mają zerową j-tą 
składową

b) wartości własne W

1

 można obliczać (

metodą 

potęgową

) skreślając uprzednio j-ty wiersz i j-

tą kolumnę, powstała macierz jest stopnia n-1

c) wartość j dobieramy tak by odpowiadała 

pozycji elementu o największym module w x

1

d) w każdym kolejnym kroku metody mamy do 

czynienia z macierzami coraz niższego stopnia

background image

15

Redukcja macierzy gęstej do prostszej 
postaci - macierzy trójdiagonalnej lub 
Hessenberga

Macierz pierwotną przekształcamy 
iteracyjnie

tak aby końcowa macierz B była macierzą 
podobną do A

0

a następnie w łatwy sposób wyznaczamy 
wartości i wektory własne macierzy B

Uwagi:

a) nakład obliczeń  potrzebnych do 
wyznaczenia x

i

 oraz 

i

 macierzy B powinien 

być jak najmniejszy

A = A

0

! A

1

! A

2

! : : : ! A

m

b) algorytm numerycznego rozwiązania problemu 
własnego macierzy B nie może być gorzej 
uwarunkowany niż dla A

Dla

problem własny B będzie gorzej uwarunkowany niż 
dla A. Ale ponieważ

więc 

uwarunkowanie algorytmu zależy od postaci

 P

i

A

i

= P

¡1

i

A

i

¡1

P

i

;

i = 1; 2; : : : ; m

B

=

A

m

= P

¡1

AP

P

=

P

1

P

2

: : : P

m

B = P

¡1

AP

¢A = P ¢BP

¡1

cond(P ) >> 1

cond(P )

=

cond(P

1

: : : P

m

)

· cond(P

1

) : : : cond(P

m

)

jBj · jP

¡1

j ¢ jAj ¢ jP j = cond(P )jAj

j¢Aj · cond(P )j¢Bj

j¢Aj

jAj

· (cond(P ))

2

j¢Bj

jBj

B + ¢B = P

¡1

(A + ¢A)P

Byyy

=

¸yyy

x

x

x

=

Pyyy = P

1

: : : P

m

yyy

Ax

x

x

=

¸x

x

x

background image

16

Redukcja macierzy hermitowskiej do 
postaci trójdiagonalnej metodą 
Hauseholdera

Pierwotna macierz hermitowska

Transformujemy

Przy pomocy 

macierzy Hauseholdera

Dokonujemy transformacji A

i-1

 do A

i

. A

i-1

 ma 

postać

A

H

= A = A

0

A

H

= (A

T

)

¤

a

a

a

i

=

2
6

6

6

6

4

®

i+1;i

¢
¢
¢

®

ni

3
7

7

7

7

5

Zauważmy że w kolejnym kroku należy dokonać 
transformacji 
a) (n-i) elementowego wektora a

i

 

b) macierzy kwadratowej

stopnia (n-i)

~

A

i

¡1

2

4

J

i

¡1

ccc

ccc

H

±

i

3

5 =

2
6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

4

±

1

¹

°

2

¢

¢

¢

0

0

°

2

. .. ...

¢

¢

¢

. .. ... ...

¢

¢

¢

. .. ... ...

¢

¢

¢

. .. ... ¹°

i

¡1

0

0

¢

¢

¢

°

i

¡1

±

i

¡1

¹

°

i

0

¢

¢

¢

0

°

i

±

i

3
7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

5

A

i

= P

¡1

i

A

i

¡1

P

i

P

i

= P

H

i

= P

¡1

i

= I

¡ ¯

i

u

u

u

i

u

u

u

H

i

A

i

¡1

=

2
6

6

6

4

J

i

¡1

ccc

0

ccc

H

±

a

a

a

H

i

0

a

a

a

i

~

A

i

¡1

3
7

7

7

5

=

jk

]

background image

17

Trzeba utworzyć macierz Hauseholdera 
stopnia (n-i)

taką  aby transformowała wektor a

i

Macierz Hauseholdera konstruujemy 
następująco

¯ =

8

<

:

1

¾(¾+

i+1;i

j)

¾

6= 0

0

¾ = 0

¾ =

kaaa

i

k

2

=

v

u

u

t

n

X

j=i+1

ji

j

2

k =

¡¾e

i'

() ®

i+1;i

= e

i'

i+1;i

j

u

u

u =

2
6

6

6

6

6

6

4

e

i'

(¾ +

i+1;i

j)

®

i+2;i

..

.

®

n;i

3
7

7

7

7

7

7

5

Dokonujemy transformacji

~

P

i

~

P

i

a

a

a

i

= k

¢ eee

1

~

P

i

= I

¡ ¯uu

uu

uu

H

P

¡1

i

A

i

¡1

P

i

= P

i

A

i

¡1

P

i

=

=

2
6

6

6

4

J

i

¡1

ccc

0

ccc

H

±

i

a

a

a

H

i

~

P

i

0

~

P

i

a

a

a

i

~

P

i

~

A

i

¡1

~

P

i

3
7

7

7

5

=

2
6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

4

±

1

¹

°

2

0

0

°

2

. .. ...

..

.

0

. .. ... ¹°

i

¡1

0

0

°

i

¡1

±

i

¡1

¹

°

i

0

: : :

0

°

i

±

i

~

°

i+1

0 : : : 0

°

i+1

0

0

..

.

~

P

i

~

A

i

¡1

~

P

i

0

3
7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

5

P

i

=

2
6

6

6

4

I

0

0

0

1

0

0

0

~

P

i

3
7

7

7

5

gi ¡ 1

background image

18

Jak dokonać efektywnego mnożenia macierzy 
(zwykłe jest nieekonomiczne)?

Definiujemy dwa wektory pomocnicze

ponieważ 

oraz z faktu, że A

i-1

 jest hermitowska wynika

°

i+1

= k

p

pp = ¯ ~

A

i

¡1

u

u

u

qqq = pp

p

¡

¯

2

(ppp

H

u

u

u)u

u

u

¯

¸ 0

p

pp

H

u

u

u = ¯u

u

u

H

~

A

i

¡1

u

u

u = (ppp

H

u

u

u)

H

Przeprowadzając proces przekształcania macierzy 
do końca dostaniemy macierz B

Metoda Hauseholdera jest stabilna numerycznie.

B = A

n

¡2

=

2
6

6

6

6

4

±

1

~

°

2

¢ ¢ ¢

0

°

2

. .. ... ...

..

.

. .. ... ~°

n

0

¢ ¢ ¢

°

n

±

n

3
7

7

7

7

5

~

P

i

= I

¡ ¯uu

uu

uu

H

~

P

i

~

A

i

¡1

~

P

i

=

(I

¡ ¯uu

uu

uu

H

) ~

A

i

¡1

(I

¡ ¯uu

uu

uu

H

)

=

~

A

i

¡1

¡ ¯ ~

A

i

¡1

uu

uu

uu

H

¡ ¯uu

uu

uu

H

~

A

i

¡1

+

¯

2

uu

uu

uu

H

~

A

i

¡1

uu

uu

uu

H

~

P

i

~

A

i

¡1

~

P

i

= ~

A

i

¡1

¡ pu

pu

pu

H

¡ up

up

up

H

+ ¯up

up

up

H

uu

uu

uu

H

= ~

A

i

¡1

¡ uuu

µ

ppp

¡

¯

2

(ppp

H

u

u

u)u

u

u

H

¡

µ

ppp

¡

¯

2

(ppp

H

u

u

u)u

u

u

u

u

u

H

= ~

A

i

¡1

¡ uq

uq

uq

H

¡ qu

qu

qu

H

background image

19

Redukcja macierzy (rzadkiej) hermitowskiej 
do postaci trójdiagonalnej metodą Lanczosa

Naszym zadaniem jest znalezienie wartości i 
wektorów własnych macierzy stopnia n

Jeśli jednak n jest bardzo duże (np.rzędu ~10

5

 ) a 

nas interesuje jedynie mały wycinek widma 
wartości (wektorów) własnych (np. m=50) to 

wtedy należy zredukować macierz do postaci 
trójdiagonalnej stopnia m

.

W metodzie Lanczosa wykorzystujemy do tego 
celu 

podprzestrzeń Kryłowa

Zakładamy że wektory rozpinające podprzestrzeń

są 

ortonormalne 

¸

1

; ¸

2

; : : : ; ¸

n

2 R

x

x

x

1

; x

x

x

2

; : : : ; x

x

x

n

2 C

n

qqq

2 C

n

K

m

(qqq; A)

=

span[qqq; Aqqq; : : : ; A

m

¡1

qqq]

m

¸ 1

K

0

(qqq; A)

=

f000g

dimK

m

=

m

K

m

= span[qqq

1

; qqq

2

; : : : ; qqq

m

]

qqq

H

i

qqq

j

= ±

ij

=

½

1

i = j

0

i

6= j

Metoda Lanczosa jest metodą iteracyjną – w każdej 
iteracji poszukujemy nowego wektora bazowego 
(wymiar podprzestrzeni zwiększa się o 1)

Sposób postępowania.
Wybieramy dowolny wektor startowy

zakładamy taże

a następnie wykonujemy rekurencyjnie ciąg obliczeń

Procedura zatrzyma się gdy

Wówczas wymiar wygenerowanej podprzestrzeni

qqq

1

2 C

n

;

qqq

1

6= 000

°

1

qqq

0

= 0

Aqqq

i

= °

i

qqq

i

¡1

+ ±

i

qqq

i

+ °

i+1

qqq

i+1

i

¸ 1

±

i

= qqq

H

i

Aqqq

i

°

i+1

=

krrr

i

k

rrr

i

= Aqqq

i

¡ ±

i

qqq

i

¡ °

i

qqq

i

¡1

qqq

i+1

=

rrr

i

°

i+1

, °

i+1

6= 0

±

i

; °

i

2 R

°

i+1

= 0

m = max

i

dim K

i

(qqq; A)

background image

20

Schemat rekurencyjny można zapisać 
macierzowo

Gdy warunek 

jest spełniony wówczas zachodzi

Q

i

= [qqq

1

; : : : ; qqq

i

]

J

i

=

2
6

6

6

6

6

6

4

±

1

°

2

0

°

2

±

2

. ..

. .. ... °

i

0

°

i

±

i

3
7

7

7

7

7

7

5

AQ

i

=

Q

i

J

i

+ [0; : : : ; 0; °

i+1

q

i+1

]

=

Q

i

J

i

+ °

i+1

qqq

i+1

eee

T

i

i

=

1; 2; : : : ; m

eee

i

= [0; 0; : : : ; 1]

T

2 R

i

°

i+1

= 0

Q

H

i

Q

i

= I

AQ

m

= Q

m

J

m

Wartości własne J

m

 są wartościami własnymi A

W skrajnym przypadku metoda nie zatrzyma się sama
i wówczas

Macierz Q

n

 jest macierzą unitarną, a macierz J

n

 jest 

macierzą trójdiagonalną podobną do A.

Uwagi praktyczne:

1)można zredukować liczbę operacji wprowadzając  
    wektor pomocniczy

wtedy
 

ponieważ

J

m

zzz = ¸zzz;

zzz

6= 000

m = n

rrr

i

= u

u

u

i

¡ ±

i

qqq

i

±

i

= qqq

H

i

Aqqq

i

= qqq

H

i

u

u

u

i

qqq

H

i

qqq

i

¡1

= 0

u

u

u

i

= Aqqq

i

¡ °

i

qqq

i

¡1

x

x

x = Q

m

zzz

6= 0

Ax

x

x = AQ

m

zzz = Q

m

J

m

zzz = ¸Q

m

zzz = ¸x

x

x

background image

21

2. Jeśli nie interesują nas wektory własne A to 
można nie przechowywać wektorów q

i

 

W każdej iteracji (i->i+1) procedura wymaga obliczenia 5 
iloczynów skalarnych oraz 1 przemnożenia wektora przez 
macierz A.

3. Jeśli interesuje nas jedynie 

Wartości (wektorów) własnych macierzy A to wówczas nie 
czekamy aż 

i+1

=0 ale przerywamy procedurę wcześniej. 

Wykorzystujemy tu fakt, że zbieżność metody dla skrajnych 
wartości własnych jest bardzo szybka.

4. Ze względu na błędy zaokrągleń wektory rozpinające 
podprzestrzeń Kryłowa  przestają być ortogonalne (gdy 
rośnie wymiar podprzestrzeni). Konieczne staje się 
przeprowadzenie tzw. 

reortogonalizacji 

nowego wektora 

bazy

do już wyznaczonych
 

Przeprowadzenie reortogonalizacji jest kosztowne.

m << n

~

qqq

i+1

qqq

i+1

= ~

qqq

i+1

¡

i

X

j=1

(qqq

H

j

~

qqq

i+1

)qqq

j

background image

22

Redukcja macierzy do postaci macierzy 
Hessenberga (górnej) metodą 
eliminacji Gaussa

Naszym zadaniem jest przekształcenie w 
sposób rekurencyjny macierzy A

do postaci Hessenberga

T – macierz trójdiagonalna
U – macierz trójkątna górna

Uwagi: 

1) każda macierz kwadratowa jest 

ortogonalnie podobna do macierzy 
Hessenberga

2) ponieważ macierz A jest 

przekształcana przez podobieństwo 
więc przekształcenie macierzy 

hermitowskiej

 do postaci 

Hessenberga prowadzi do uzyskania 
macierzy trójdiagonalnej ze względu 
na symetrię A. Macierz trójdiagonalna 
jest hermitowska.

A = A

0

! A

1

! : : : ! A

n

¡2

= B

H

A

i

= P

¡1

i

A

i

¡1

P

i

W każdej iteracji macierz przekształcenia P

i

 

konstruujemy z dwóch macierzy:

a) macierzy permutacji

macierz do niej odwrotna

¼

rs

=

2
6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

4

1

0

. ..

1

0

1

1

. ..

1

1

0

1

. ..

1

3
7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

5

àr

às

¼

¡1

rs

= ¼

rs

B

H

=

2
6

6

6

6

4

¤ ¤ ¤ ¤ ¤
¤ ¤ ¤ ¤ ¤
0

¤ ¤ ¤ ¤

0

0

¤ ¤ ¤

0

0

0

¤ ¤

3
7

7

7

7

5

= T + U

background image

23

b) macierzy eliminacji

z  warunkiem

Macierz odwrotna G

j

-1

G

j

=

2
6

6

6

6

6

6

6

6

6

6

6

6

4

1

. ..

1

l

j+1;j

1

..

.

. ..

l

nj

1

3
7

7

7

7

7

7

7

7

7

7

7

7

5

l

ij

· 1

G

¡1

j

=

2
6

6

6

6

6

6

6

6

6

6

6

6

4

1

. ..

1

¡l

j+1;j

1

..

.

. ..

¡l

nj

1

3
7

7

7

7

7

7

7

7

7

7

7

7

5

Jakie są skutki mnożenia macierzy G

j

 oraz 

rs

  z A?

1.                              - zamienia wiersze r i s  w A

2.                              - zamienia kolumny r i s w A         
    

2.                              - odejmuje wiersz j-ty
                                   przemnożony przez l

rj

 od

                                   wiersza r w A

3.                              - przemnaża kolumnę r-tą przez
                                   l

rj

 a następnie dodaje do

                                   kolumny j-tej w A 

¼

¡1

rs

A

G

¡1

j

A

AG

j

rs

background image

24

Załóżmy, że w macierzy A

i-1

 (i-1) pierwszych 

kolumn posiada postać Hessenberga

A

i

¡1

=

2
6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

6

4

¤ ¢ ¢ ¢ ¢ ¤ ¤ ¢ ¢ ¢ ¤

¤ ¢

¢

¢

¢

¢ ¢

¢

¢

¢

¢ ¢

¢

¢

¢

¢ ¢

¢

¢

¢

0

¤ ¤ ¤ ¢ ¢ ¢ ¤

0

¢ ¢ ¢ 0 ¤ ¤ ¢ ¢ ¢ ¤

¢

¢

¢

¢

¢

¢

0

¤ ¢ ¢ ¢ ¤

3
7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

7

5

=

2
6

6

6

4

B

i

¡1

d

d

d

b

A

i

¡1

ccc

±

i

bbb

0

a

a

a

~

A

i

¡1

3
7

7

7

5

a

a

a =

2
6

6

6

4

®

i+1;i

..

.

®

ni

3
7

7

7

5

background image

25

Dokonujemy przekształcenia

macierz przekształcenia jest zdefiniowana w 
poniższy sposób

Macierz przekształcenia P

i

 nie zmienia postaci 

Hessenberga w (i-1) pierwszych kolumnach.
Zmienia natomiast kolumnę i-tą tj. zeruje 
elementy w tej kolumnie od i+2 do n.

A

i

= P

¡1

i

A

i

¡1

P

i

P

¡1

i

= G

¡1

i+1

¼

¡1

r;i+1

P

i

= ¼

r;i+1

G

i+1

r

¸ i + 1

¹

a

a

a =

2

6

6

6

4

¤
0

..

.

0

3

7

7

7

5

Aby określić które wiersze należy zamienić wierszami, 
najpierw trzeba określić element o największym module 
a

Po wyznaczeniu wartości r zamieniamy miejscami 
wiersze (r ↔ i+1)  oraz kolumny (r ↔ i+1) w A

i-1

Dzięki temu element o największym module w wektorze 
a znajduje się teraz na jego pierwszym miejscu.

Jak określić G

j+1

Warunek: trzeba wyzerować wszystkie elementy w a 
poza jego pierwszym elementem

Przeprowadzając drugi etap przekształcenia

Dostajemy macierz Hessenberga w i kolumnach 
macierzy A

i

ri

j = max

i+1

·j·n

ji

j;

r

¸ i + 1

P

¡1

i

¢

2

6

6

6

4

d

d

d

±

i

a

a

a

3

7

7

7

5

=

2

6

6

6

4

d

d

d

±

i

¹

a

a

a

3

7

7

7

5

A

0

= ¼

¡1

r;i+1

A

i

¡1

¼

r;i+1

= [®

0

jk

]

l

j;i+1

=

8

<

:

®

0
ji

®

0
i+1;i

®

i+1;i

6= 0

0

®

i+1;i

= 0

j

=

i + 2; i + 3; : : : ; n

A

i

= G

¡1

i+1

A

0

G

i+1

= P

¡1

i

A

¡1

P

i

background image

26

Uwagi:

1) Po n-2 iteracjach macierz A zostaje 

przekształcona do postaci Hessenberga 
(górnej).

2) Metoda eliminacji Gaussa jest stabilna 

numerycznie i wymaga wykonania

operacji mnożenia 

3) Redukcję macierzy do postaci Hessenberga 

można przeprowadzić także przy użyciu 
metody Hauseholdera lub Givensa
W metodzie Givensa trzeba wykonać 

a w metodzie Hausholdera

operacji mnożenia
Ze względu na błędy numeryczne, w obu 
metodach macierz Hessenberga jest 
macierzą podobną do A+E 

M =

2
3

n

3

+ O(n

2

)

M =

10

3

n

3

+ O(n

2

)

M =

5
3

n

3

+ O(n

2

)

kE

Givens

k · k

1

"n

3=2

kAk;

k

1

» 1

kE

Hauseholder

k · k

2

"n

2

kAk;

k

2

» 1

background image

27

Wyznaczanie wartości własnych 
macierzy

Wyznaczanie wartości własnych 
macierzy trójdiagonalnej hermitowskiej
metodą bisekcji

Po zredukowaniu macierzy A do postaci

chcemy znaleźć sposób na wyznaczenie 
wartości własnych J.
Jeśli warunek

to macierz jest 

nieredukowalna

.

W przeciwnym wypadku można ją zapisać w 
postaci

J =

2
6

6

6

6

6

6

4

±

1

¹

°

2

0

°

2

. .. ...
. .. ... ¹°

n

0

°

n

±

n

3
7

7

7

7

7

7

5

°

i

6= 0; i = 2; : : : ; n

Macierze

są nieredukowalne. Ponadto widmo wartości własnych J 
pokrywa się z widmem macierzy nieredukowalnych J

i

można więc badać je osobno.
Szukamy wielomianu charakterystycznego J

i

: W(

) 

rozwijając wyznacznik względem kolejnych kolumn 
macierzy 

Macierz jest hermitowska więc 

Do znalezienia wartości własnych można wykorzystać 
metodę bisekcji.

J

i

; k = 1; : : : ; k

J =

2
6

6

6

6

6

6

4

J

1

0

J

2

. ..

0

J

k

3
7

7

7

7

7

7

5

±

i

;

i

j

2

2 R

!

i

(¸) = det(J

i

¡ ¸I)

!

0

(¸)

=

1

!

1

(¸)

=

±

1

¡ ¸

: : :

: : :

: : :

!

i

(¸)

=

i

¡ ¸)!

i

¡1

(¸)

¡ j°

2

i

j!

i

¡2

(¸)

i

=

2; 3; : : : ; n

W (¸)

=

!

n

(¸)

background image

28

Sposób postępowania:

1) wybieramy dowolną liczbę 

 i obliczamy 

wartość wielomianu charakterystycznego 
rekurencyjnie

2) następnie korzystamy z poniższych twierdzeń:

 
Tw. Jeżeli elementy 

3

, ..., 

n

  (pozadiagonalne) są 

niezerowe, to wartości własne macierzy J są 
pojedyncze.

Tw. Jeżeli elementy 

2

3

, ..., 

n

  (pozadiagonalne) są 

niezerowe, to ciąg wartości 
               



 

n

 

spełnia warunki:

a) jeżeli 

i

  dla pewnego i<n, to 

            

i-1



i+1



b) jeżeli 

n

    jest różne od 0, to liczba zmian 

znaków sąsiednich liczb
                

0



1



n

  

jest równa liczbie wartości własnych macierzy J 
mniejszych od 

 

c) Jeżeli 

n

, to  jest wartością własną macierzy 

J, a ponadto jest tyle wartości   własnych 
mniejszych niż 

 ile nastąpiło zmian znaków w 

ciągu  

0



1



n-1



Metoda bisekcji jest bardzo dokładna. Wadą jest 
uzyskiwanie dużych wartości ciągu:  

0



1



n

  jeśli  znacznie różni się od 

wartości własnych J. Zaletą natomiast możliwość 
obliczenia wartości własnej o określonym indeksie 
k. Liczba iteracji potrzebna do wyznaczenia 

k

  

wynosi:



 

przedział poszukiwań wartości własnej



 

dokładność wyznaczenia wartości własnej

Wektory własne
Znając wartość własną macierzy J wektor własny 
wyznaczamy według wzorów:

gdzie: 

i

 – elementy diagonalne J

           

i

 – elementy pozadiagonalne J

IT = log

2

¯

0

¡ ®

0

½

x

1

= 1

x

2

=

¸

¡ ±

1

°

2

x

i+1

=

¡ ±

i

)x

i

¡ °

i

x

i

¡1

°

i+1

;

i = 2; 3; : : : ; n

¡ 1

background image

29

Algorytm poszukiwania 

k-tej

 wartości własnych metodą 

bisekcji

¯

0

=

kJk

®

0

=

¡kJk

background image

30

Wyznaczanie wartości własnych 
metodą LR

W metodzie tej iteracyjnie przekształcamy 
macierz A uzyskując ciąg

W którym ostatni element stanowi macierz 
trójkątna górna (Hessenberga). Elementy 
diagonalne macierzy A

m

 stanowią natomiast 

ciąg wartości własnych macierzy A czyli

W każdej iteracji wyznaczamy rozkład A

i

 na 

iloczyn macierzy trójkątnej dolnej L z 
jedynkami na diagonali oraz macierzy 
trójkątnej górnej R:

Przekształcamy macierz w następujący sposób

Macierze A

i

 oraz A

i+1

 są podobne 

A = A

1

! A

2

! A

3

! : : : ! A

m

A

i

= L

i

R

i

A

i+1

= L

i+1

R

i+1

= R

i

L

i

lim

i

!1

a

(i)
jj

= ¸

j

A

i

= L

i

R

i

;

L

¡1

i

A

i

= R

i

;

A

i

R

¡1

i

= L

i

A

i+1

= R

i

L

i

= L

¡1

i

AL

i

= R

i

AR

¡1

i

Wadą jest to że rozkład L

i

R

i

 może nie istnieć lub może 

istnieć ale jego znalezienie jest źle uwarunkowane. W 
obu przypadkach algorytm przerywa działanie.

Wyznaczanie wartości własnych 
metodą QR

Metoda wywodzi się z metody LR, przy czym macierz L 
zastąpiono 

macierzą ortogonalną

 Q przez co metoda 

jest stabilna numerycznie.

gdzie: R

i

 jest macierzą trójkątna górną, a Q

i

 jest 

macierzą ortogonalną
W metodzie QR otrzymujemy ciąg macierzy 

Macierze te są do siebie podobne więc mają te same 
wartości własne. Jeśli m jest duże wówczas 
spodziewamy się że na diagonali A

m

 będą znajdować się 

wartości własne A.

Tw. Jeżeli macierz A jest symetryczna i dodatnio 
określona, to algorytm QR generuje ciąg macierzy  A

(1

)

,A

(2)

,...  zbieżny do macierzy diagonalnej.

Wadą metody QR jest wolna zbieżność dla macierzy 
pełnych. Metoda jest szybkozbieżna dla 
macierzy trójdiagonalnych i macierzy Hessenberga.

Q

H

Q = I

A

i

= Q

i

R

i

A

1

= A

A = A

1

! A

2

! A

3

! : : : ! A

m

A

i+1

= R

i

Q

i

background image

31

Algorytm QR dla macierzy 
Hessenberga.

Macierzą (górną) Hessenberga jest macierz, 
którą można zapisać w postaci: 

Każda macierz kwadratowa jest 
ortogonalnie podobna do macierzy 
Hessenberga, więc ma ona to samo widmo 
wartości własnych co macierz H.

Wyznaczamy ciąg macierzy:

Wszystkie H

i

 są macierzami Hessenberga

wobec czego wszystkie macierze H

i

i=1,2,3,.. są podobne. Elementy na 
diagonali kolejnych H

i

 dążą do wartości 

własnych macierzy H

1

=H.

Jeżeli macierz H ma pojedyncze wartości 
własne takie, że są pomiędzy nimi pary
sprzężone o identycznych modułach, to nie 
wszystkie elementy poddiagonalne dążą 
do 0. W granicy otrzymuje się macierz:  

H = T + U

H

i

= Q

i

R

i

H

i+1

= R

i

Q

i

;

i = 1; 2; 3; : : :

H

i+1

= Q

T

i

H

i

Q

i

H

1

=

2
6

6

6

6

4

¤ ¤ ¢

¢

¢

¤ ¤ ¢

¢

¢

0

0

¤ ¢

¢

0

0

0

¤ ¤

0

0

0

¤ ¤

3
7

7

7

7

5

gwiazdki - elementy ustalone
kropki – elementy, których wartości mogą się 
zmieniać.

Wartości własne leżą wtedy bezpośrednio na 
diagonali lub są wartościami własnymi macierzy 2x2 
leżących na diagonali. 

Wadą metody może być powolna zbieżność. 
Zwiększenie wydajności uzyskuje się dokonując 
przesunięcia widma wartości własnych 

Wartość k

i

 wybiera się jako równe otrzymanym już 

wartościom własnym  

lub wartości k

i

 oraz k

i+1

 jako równe wartościom 

własnym macierzy 2x2 znajdujących się w prawym 
dolnym rogu macierzy H

i

H

i

¡ k

i

I = Q

i

R

i

H

i+1

= R

i

Q

i

+ k

i

I

k

i

= a

i

nn

background image

32

Metoda Hauseholdera rozkładu QR

Definiujemy macierz Hauseholdera w 
postaci

która ma własność

Macierz Hausholdera posłuży do znalezienia 
rozkładu QR.
Określamy ciąg macierzy P

(1)

,P

(2)

,P

(3)

,...P

(n-1)

 

przy pomocy których definiujemy macierz R
(górną trójkatną):

¿ =

kzk

2

(

kzk

2

¡ ®)

P

(n

¡1)

P

(n

¡2)

: : : P

(1)

A = R

Zakładamy

Macierz H

(1)

 jest macierzą Hauseholdera sprowadzającą 

pierwsza kolumnę macierzy A (a

(1)

1

)do postaci:

Przez P

(2)

 oznaczamy

Macierz H

(2)

 sprowadza pierwszą kolumnę macierzy o 

wymiarze (n-1)x(n-1) utworzonej z wierszy i kolumn o 
numerach 2,3,4,...,n macierzy P

(1)

A (a

(2)

1

) do postaci

Postępujemy w ten sposób otrzymując kolejno P

(3)

,P

(4)

itd.

P

(1)

= H

(1)

®

(1)

ka

(1)
1

k

2

[1; 0; 0; : : : ; 0]

T

1

£n

P

(2)

=

2

6

6

4

1

..

.

0

: : :

: : :

: : :

0

..

.

H

(2)

3

7

7

5

®

2

ka

(2)
1

k

2

[1; 0; 0; : : : ; 0]

T

1

£(n¡1)

H = I

¡

1

¿

u

u

uu

u

u

H

u

u

u = zzz

¡ ®kzzzk

2

eee

1

zzz = [z

1

; z

2

; : : : ; z

n

]

T

eee

1

= [1; 0; : : : ; 0]

T

Hzzz = ®

kzzzk

2

[1; 0; 0; : : : ; 0]

T

® =

§1 = ¡sign(z

1

)

background image

33

Macierz P

(n-1)

 ma postać

A macierz H

(n-1)

 ma wymiar 2x2.

Po wyznaczeniu wszystkich macierzy P

(i)

rozkład A=QR wyznaczamy według 
wzorów

Liczba operacji potrzebnych do uzyskania 
rozkładu Hauseholdera wynosi

a) mnożenie

b) dodawanie

c) pierwiastkowanie

P

(n

¡1)

=

2

6

6

4

1

1

1

H

(n

¡1)

3

7

7

5

Q = P

(1)

P

(2)

P

(3)

: : : P

(n

¡1)

R = Q

T

A

M =

2
3

n

3

+ 2n

2

¡

2
3

n

¡ 2

D =

2
3

n

3

+

10

3

n

¡ 4

p = n

¡ 1

Uwagi:

1) inne metody pozwalające uzyskać rozkład QR to 

metoda Grama-Schmidta  i Givensa

2) jeśli macierz A jest symetryczna to rozkład QR 

zachowuje jej symetrię natomiast rozkład LR nie

3) jeśli macierz A jest symetryczna i dodatniookreślona 

to można stosować rozkład LL

 wówczas algorytm LR 

zachowuje symetrię macierzy A 

4) metoda szybka dla macierzy Hessenberga i 

trójdiagonalnej

Wyznaczanie wektorów własnych dla rozkładu QR

Uwaga: 

nie zakładamy postaci macierzy pierwotnej A

W metodzie QR dostajemy ciąg przekształceń

co można uogólnić

A

i

= Q

i

R

i

! Q

¡1

i

A

i

= R

i

A

i+1

= R

i

Q

i

= Q

¡1

i

AQ

i

= Q

i+1

R

i+1

A

i+2

= R

i+1

Q

i+1

= Q

¡1

i+1

Q

¡1

i

AQ

i

Q

i+1

A

k+1

= Q

¡1

k

Q

¡1

k

¡1

: : : Q

¡1

1

AQ

1

Q

2

: : : Q

k

Q

¡1

i+1

Q

¡1

i

AQ

i

= R

i+1

background image

34

Oznaczmy 

Wówczas można bezpośrednio pokazać że 

Macierz A

k+1

 jest macierzą trójkątną górną z 

wartościami własnymi na diagonali.
Wektor własny x

(i)

 tej macierzy (

H

przynależny do wartości własnej 

wyznacza się według poniższych wzorów

P = P

k

= Q

1

Q

2

: : : Q

k

P

¡1

= P

¡1

k

= Q

¡1

k

Q

¡1

k

¡1

: : : Q

¡1

1

P

¡1

AP = A

k+1

= H

¸

i

= h

ii

x

(i)
j

=

0;

j = n; n

¡ 1; : : : ; i + 1

x

(i)
i

=

1

x

(i)
j

=

¡

P

i
k=j+1

h

jk

x

(i)
k

h

jj

¡ h

ii

;

j = i

¡ 1; i ¡ 2; : : : ; 1

Mając wyznaczone wektory własne macierzy H możemy 
obliczyć wektory własne macierzy pierwotnej A

P

¡1

AP = H

Hx

x

x = ¸x

x

x

P

¡1

APx

x

x = Hx

x

x = ¸x

x

x

A(Px

x

x) = ¸Px

x

x

yyy = Px

x

x

Ayyy = ¸yyy

background image

35

Uogólniony problem własny

Uogólniony problem własny definiujemy 
następująco

gdzie: A i B są macierzami kwadratowymi.

Najprościej byłoby przekształcić powyższe 
równanie tak aby przeprowadzić je do zwykłego 
problemu własnego tj.

Problemem pozostaje jednak jak znaleźć B

-1

?

W przypadku, gdy B oraz A są  macierzami 
symetrycznymi możemy posłużyć się rozkładem 
LL

(w ogólnym przypadku można skorzystać z 

rozkładu LU)

wówczas wykorzystując rozkład LL

T

 można 

znaleźć macierz podobną do B

-1

A

B = LL

T

BB

¡1

= I = LL

T

(L

T

)

¡1

L

¡1

B

¡1

= (L

T

)

¡1

L

¡1

Dzięki temu przekształceniu, macierz G jest 
symetryczna jak A i posiada identyczne widmo 
wartości własnych (ale inne wektory własne).

Jak znaleźć G?
Najpierw należy znaleźć macierz F

  

rozwiązując układ równań

a następnie wyznaczamy G

rozwiązując układ równań

Rozkład LL

T

 wymaga wykonania n

3

/6 mnożeń a 

wyznaczenie macierzy G (2/3)n

3

. Macierz G jest 

symetryczna więc w celu wyznaczenia jej wartości i 
wektorów własnych korzystamy z metod 
przeznaczonych dla tej klasy macierzy.

       Wektory własne macierzy A wyznaczamy    
       przekształcając wektory macierzy G lub 
       rozwiązując układ

G = L

¡1

F

LG = F

Ax

x

x = ¸Bx

x

x

B

¡1

Ax

x

x = Cx

x

x = ¸x

x

x

L

T

(B

¡1

A)(L

T

)

¡1

=

L

T

(L

T

)

¡1

L

¡1

A(L

T

)

¡1

=

L

¡1

A(L

T

)

¡1

=

G

F = A(L

T

)

¡1

F L

T

= A =

) LF

T

= A

T

= A

Gyyy = ¸yyy

L

T

x

x

x = yyy