3 - redukcja macierzy symetrycznej do postaci trojgonalnej

1

Redukcja macierzy symetrycznej do postaci trojdiagonalnej

1.1

Algorytm Householdera

Algorytm Householdera redukuje macierz symetryczna n  n A do postaci trojdiagonalnej poprzez n 2 operacji. Kazda operacja anihiluje rzadana czesc kolumny oraz odpowiadajacego jej rzedu. Podstawa algorytm jest macierz P postaci u
u P =1 2 (1) juj gdzie u jest dowolnym rzeczywistym wektorem spelniajacym a u
u to iloczyn zewnetrzny wektorow czyli macierz. Macierz P jest ortogonalna poniewaz u
u u
u P = (1 2 )(1 2 juj juj ) = 1 4 uju
uj + 4 u
(u ju
ju)
u (2) = 1 Czyli P = P . Poniewaz rowniez P = P , wiec P = P co dowodzi ortogonalnosci. Przypuscmy, ze x jest wektorem zbudowanym z pierwszej kolumny A. Wybieramy u = x  jxje (3) gdzie e to wektor jednostkowy [1; 0;


; 0] . Wtedy u
u u P
x = (1 2 juj )x = x 2 juj (x  jxje ) x = x 22ju(jxxj j 2jjxxjjxx ) = x u = jxje (4) To pokazuje ze macierz Hauseholdera P dzialajac na dany wektor x zeruje wszystkie jego elementy poza pierwszym. Aby zredukowac symetryczna macierz A do postaci trojdiagonalnej wybieramy jako wektor x dla pierwszej macierzy Householdera wektor zbudowany z n 1 dolnych elementow pierwszej kolumny A. Wtedy jej dolne n 2 elementow zostanie wyzerowanych 2 1 j 0


0 32 a j a


a 3 T

2

T

T

2

T

2

2

T

T

T

2

1

4

T

1

T

1

1

T

T

2

2

1

2

6 6 P1
A = 6 6 6 4

1

2

1

1

11

0 j .. . j 0 j 2

6 6 6 k = 666 0 6 .. 4 .

j

12

76 76 7 6 a21 76 7 6 .. 54 .

P10 a11

T

j j j

an1






a12

a1n

1n

7 7 7 7 7 5

3 7 7 7 7 7 7 7 5

j j nieistotne j 0 j gdzie P 0 to macierz Hauseholdera o rozmiarze (n 1)  (n 1) a powstale k to 1

0

k = P10 B @

a21

.. .

an1

Pelna ortogonalna transformacja ma postac 2 a j a 11

6 6 6 k 0 A =P
A
P =6 6 0 6 6 .. 4 .

0

12

j j j j

1

0

C B A = j @




nieistotne

1

a21

.. .

an1 a1n

1

C Aj 3

2 7 76 76 76 76 76 74 5

1 j 0


0 0 j .. . j 0 j

P10

3 7 7 7 7 7 5

2

j

a11

6 6 6 k = 666 0 6 .. 4 .

0

k






3

0

7 7

7 j 7 7 j nieistotne 7 7 5 j 0 j Widac, ze P AP zaczyna byc trojdiagonalne. Wykonujac w analogiczny sposob dalsze kroki iteracji po n 2 krokach dostaniemy macierz trojdiagonalna z macierza podobienstwa Q = P


P . Jezeli znajdziemy wektory wlasne macierzy A0 przemnozenie ich przez Q da nam wektory wlasne A. Zamiast mnozenia P
A
P w rzeczywistosci wpierw znajduje sie wektor A
u (5) p=2 juj a nastepnie u
u A
P = A
(1 2 (6) juj ) A0 = P
A
P = A p
u u
p + 2Ku
u (7) gdzie skalar K jest zde niowany jako u
p K= (8) juj Jezeli sie zapisze 1

1

1

n

2

2

T

2

T

T

T

T

2

q=p

otrzymuje sie

A0 = A

q
uT

co stanowi uzyteczna praktycznie formule.

2

(9)

Ku u
qT

(10)

Wektory i wartosci wlasne macierzy trojdiagonalnej

2.1

Algorytm QR i QL

Jest to metoda bardzo prosta od strony algorytmicznej lecz trudno wykazuje sie jej wlasnosci teoretyczne i numeryczne. Mozna ja stosowac rowniez dla macierzy pelnych. Wtedy jednak wymaga n operacji na iteracje zas w przypadku macierzy trojdiagonalnych liczba ta maleje do n operacji na iteracje. Metoda opiera sie na: 3

Twierdzenie Schura

Dla dowolnej macierzy A istnieje taka unitarna macierz Q, ze R = Q AQ (11) jest macierza gorna trojkatna. Wartosciami wlasnymi A sa zatem elementy diagonalne R. Macierzy Q nie da sie wyznaczyc przy pomocy skonczonej liczby dzialan arytmetycznych. Twierdzenie Schura sugeruje jednak mozliwosc budowy ciagu macierzy Q~ takich, ze Q~ AQ~ daza do macierzy gornej trojkatnej, ktorej diagonalnymi elementami sa wartosci wlasne  macierzy A. Mozna to uzyskac poprzez zde niowanie ciagu: A = A A = Q R k = 1; 2;


(12) A = RA Mozna zuwazyc, ze (13) A = Q Q R Q = Q =


= Q~ AQ~ 1

1

k

k

k

i

1

k+1

k

k

k

k+1

k

k

T k

k

k

1

k

k

2

1

k

k

gdzie

Q~ k

= Q Q


Q (14) Wszystkie macierze A sa wiec podobne do A. Jesli A ma n wartosci wlasnych o roznych modulach j j > j j >


> j j > 0 (15) to mozna pokazac, ze dla i = 1; 2;


n elementy a daza przy k ! inf do  , zas dla i > j elementy a sa zbiezne do zera. Dla dowolnych macierzy A konstruowany ciag A nie ma takich wlasnosci. Trudnosc stanowi rowniez wyznaczanie macierzy Q . 1

2

k

k

1

2

n

(k)

i

ii

(k)

k

ij

k

W przypadku macierzy trojdiagonalnej algorytm sprowadza sie do wymnazania A 2 3 a b


7 A=6 4 b c


5 .. . przez macierze obrotu   cos  sin  sin  cos  

=



a b b c   cos  sin  sin  cos 



b(cos2 

b(cos2 

cos

sin ) + (a 2

c) sin  cos 



sin ) + (a c) sin  cos  cos i dobieraniu  tak aby elementy niediagonalne znikaly cos  sin  = c a (16) sin  cos  b Ale cos  sin  = 1 tan  =  (17) sin  cos  2 tan  gdzie oznaczylismy  = . Oznaczajac dalej tan  = t mamy t + 2t 1 = 0 (18) Z wzorow Vieta widac, ze t t = 1, a zatem jeden z pierwiastkow spelnia jt j < 1. Ze wzgledow numerycznych bierzemy ten pierwiastek. Mamy p (19) t =  1+ i poniewaz t t = 1 mamy rowniez p1 t = (20) + 1+ p1 t =  1+ Aby uniknac odejmowania liczb zblizonych co do modulu w mianowniku bierzemy zawsze ten pierwiastek, w ktorym w mianowniku bedzie dodawanie. Daje to zawsze obrot o =4 sgn() p t= = tan  (21) jj + 1 +  i ostatecznie 2

2

2

2

2

c

2

a

b

2

1 2

i

2

1;2

1 2

1

2

2

2

2

cos  = p 1 1+t sin  = t
cos  = p

(22)

2

3

t

1+t

2