Ellen®rz® kérdések
4. Sajátérték-feladatok numerikus megoldásamegoldása
4.3. A sajátértékeket egyszerre közelít® eljárások
Blokk háromszögmátrix deáció
Egy A∈Rn×n mátrixot blokk fels® háromszögmátrixnak nevezünk, ha vannak olyan Akl(k, l= 1, . . . , m)mátrixok (A11, . . . ,Amm kvadratikus mátrixok), melyekkel
A blokk alsó háromszögmátrix hasonlóan deniálható. Ha egy mátrix blokk háromszögmátrix alakú, akkor a sajátértékei a diagonálisban szerepl® A11, . . . ,Amm mátrixok sajátértékeinek uniójaként adódnak. (Ez könnyen igazolható a karakterisztikus polinom vizsgálatával.) Így a sajátérték-meghatározást visszavezethetjük kisebb méret¶ mátrixok sajátértékeinek meghatáro-zására. Különösen gyakran alkalmazzuk az eljárást akkor, ham= 2, és A11 (n−1)×(n−1)-es, A22 pedig egy 1×1-es mátrix. Gyakori az az eset is, hogy a mátrix egyszer¶ permutációkkal végzett hasonlósági transzformációkkal blokk háromszögmátrix alakra hozható (lásd pl. 4.5.2.
feladat).
4.3. A sajátértékeket egyszerre közelít® eljárások
Most rátérünk azokra a módszerekre, amelyek a mátrixok minden sajátértékére egyszerre mon-danak közelítést. Ezen módszerek közös alapgondolata az, hogy ha egy mátrixszal hasonlósági transzformációt végzünk, akkor annak sajátértékei nem változnak meg. Tehát ha pl. el® tudnánk állítani egy mátrix Schur-felbontását, azaz ha fel tudnánk írni az A∈Rn×nmátrixot A=STSH alakban, ahol T fels® háromszögmátrix és S unitér mátrix, akkor az SHAS=T egyenl®ség miatt A sajátértékei megegyeznének T diagonális elemeivel. A f® problémát az S unitér mátrix meg-keresése jelenti. Direkt módszerrel S általában nem határozható meg. Így iterációs módszereket kell találnunk, melyek az S mátrixot egy mátrixsorozat határértékeként állítják el®. Két módszert ismertetünk most ennek megvalósítására. A továbbiakban ismét valós szimmetrikus mátrixokkal foglalkozunk csak, így a Schur-felbontás tulajdonképpen az A=SΛST felbontást jelenti, ahol S az ortonormált sajátvektorok mátrixa,Λpedig a sajátértékek diagonálmátrixa.
4.3.1. A Jacobi-módszer
A most ismertetend® módszert el®ször Jacobi írta le 1845-ben.3Legyen A∈Rn×n egy valós szim-metrikus mátrix. Keresünk egy olyan S ortogonális mátrixot ésΛdiagonális mátrixot, melyekkel STAS=Λ.
Az S mátrixot 2×2-es szimmetrikus mátrixokra könnyen meghatározhatjuk. Legyen pl.
A=
3Elméleti eredményként jó darabig feledésbe merült a módszer. A számítógépek elterjedésével az 1950-es évek-t®l kezdték nagyméret¶ mátrixok sajárértékeinek meghatározására használni. Manapság párhuzamosíthatósága miatt népszer¶.
alakban, ahol s = sinθ és c = cosθ, valamilyen θ paraméterrel. (Az S mátrix elemeinek ezen megválasztása mellett S valóban ortogonális mátrix lesz.) A c éss értékeknek olyanoknak kell lenniük, hogy a mátrix diagonális legyen, azaz igaz legyen a
0 =sca−s2b+c2b−scd=sc(a−d) +b(c2−s2) = (a−d) sin(2θ)/2 +bcos(2θ) egyenl®ség. Továbbá nyilvánvalóans2+c2= 1.
Haa=d, akkor acos(2θ) = 0egyenl®ségnek kell teljesülnie. Ezt és azs2+c2= 1összefüggést felhasználvaséscmeghatározható.
Ha a6=d, akkor tg(2θ) =−2b/(a−d). Ebb®l és az s2+c2 = 1 összefüggésb®l sés c ismét meghatározható (lásd4.5.8.feladat).
Ezzel a 2×2-es szimmetrikus mátrixokhoz megadtunk olyan ortogonális mátrixot, amely a hasonlósági transzformáció során azt diagonális mátrixba viszi. Nagyobb méret¶ mátrixoknál azonban ez nehéz feladat. A Jacobi-módszer alapötlete az, hogy nagyobb méret¶ mátrixok esetén is a 2×2-es esetre levezetett diagonalizáló mátrixot használjuk, csak sokszor egymás után. Az algoritmus a következ®.
• Válasszuk ki a mátrix f®átlója felett a legnagyobb abszolút érték¶ elemet. Legyen ez azaij
elem (i < j).
• Határozzuk meg a diagonalizáló S mátrixot az aii aij
aij ajj
2×2-es szimmetrikus mátrixra, azaz adjuk meg a benne szerepl®sésc paramétereket.
• Deniáljuk az alábbin×n-es ortogonális mátrixot:
Sij = [e1, . . . ,ei−1, cei−sej,ei+1, . . . ,ej−1, sei+cej,ej+1, . . . ,en],
Azt mondjuk, hogy a fenti eljárással végrehajtottunk egy (i, j) elemhez tartozó Jacobi-transzformációt.
4.3. A sajátértékeket egyszerre közelít® eljárások 121
• Ezután visszalépünk az eljárás elejére újra kiválasztva a f®átló felett a legnagyobb abszolút érték¶ elemet.
Könny¶ látni, hogy minden lépésben az (i, j) és (j, i) helyeken álló f®átlón kívüli elemek lenullázódnak, de egy következ® lépés újra nullától különböz®vé teheti ®ket. Mutassuk meg, hogy az algoritmus valóban a várt eredményt adja! Ehhez igazolnunk kell, hogy az algoritmussal nyert Sij mátrixok szorzata határértékben éppen az S diagonalizáló mátrixhoz fog tartani.
4.3.1. tétel.
Tegyük fel, hogy a fent ismertetett módon végrehajtottunk egy(i, j)elemhez tartozó Jacobi-transzformációt. Ekkor az A(1) mátrix f®átlón kívüli elemeinek négyzetösszege2a2ij-tel keve-sebb lesz, mint az A mátrix esetén.
Bizonyítás: A Jacobi-transzformáció csak azi-edik és aj-edik sort és oszlopot változtathatja meg. Nyilvánvalóanl6=i, j esetén ezen oszlopok és sorok elemei az új mátrixban
a(1)il =cail−sajl, a(1)jl =cajl+sail,
módon látható, a másik pedig hasonlóan igazolható. Számítsuk ki most a négyzetösszeg-változást a mátrix f®átló feletti részén:
−a2ij+
Így a teljes f®átlón kívüli rész négyzetösszeg-változása valóban−2a2ij.
Legyen a f®átlón kívüli elemeket tartalmazó mátrix ak-adik lépésben Bk=A(k)−diag(diag(A(k))). Ekkor
kBkk2F =kBk−1k2F−2(a(k−1)ij )2≤ kBk−1k2F−2kBk−1k2F n(n−1)
=kBk−1k2F
1− 2
n(n−1)
, így
kBkk2F ≤ kB0k2F
1− 2 n(n−1)
k .
Ez a becslés mutatja, hogy a módszer konvergenciarendje legalább egy (valójában másodrend¶).
A Jacobi-transzformációkat végrehajtva tehát a mátrix f®átlón kívüli elemei nullához tartanak.
Így az Sij mátrixok szorzata valóban az S diagonalizáló mátrixhoz fog tartani. Elegend®en sok lépést végrehajtva a f®átlón kívüli elemek annyira lecsökennek, hogy a f®átlóban megkapjuk a sajátértékek közelítéseit a Gersgorin-tételt alkalmazva.
4.3.2. megjegyzés. A legnagyobb elem kiválasztásan(n−1)/2összehasonlítást jelent, ami sokáig tart. Ezért inkább az a szokás, hogy a f®átló feletti részen soronként haladunk végig a Jacobi-transzformációkkal (azaz ilyenkor nem a legnagyobb abszolút érték¶ elemet választjuk ki ez a módosított Jacobi-módszer). A tapasztalat szerint 5-szöri végighaladás elég szokott lenni egy megfelel® sajátérték-közelítéshez.
4.3.3. megjegyzés. A módosított Jacobi-módszer esetén egy Jacobi-transzformáció4nszorzást és2nösszeadást igényel, ésn(n−1)/2elemen 5-ször kell végighaladni. Ez kb.15n3op m¶veletet igényel.
4.3.4. megjegyzés. Ha egymás után az Si1j1, . . . ,Sikjk ortogonális mátrixokat használtuk a Jacobi-transzformációhoz, akkor a λj sajátértékhez tartozó sajátvektor az Si1j1 ·. . .·Sikjkej
módon közelíthet®.
4.3.5. megjegyzés. A Jacobi-módszer olyan mátrixokra is konvergál, melyek nem feltétlenül szimmetrikusak, de minden sajátértékük valós.
4.3.2. QR-iteráció
A módszer alapötlete az, hogy az ortogonális diagonalizáló mátrixot az A mátrix QR-felbontásának Q mátrixával közelítjük. Az algoritmus a következ®:
QR-iteráció, A adott mátrix, A(0):=A fork:= 1 :kmax do
Készítsük el A(k−1)QR-felbontását: A(k−1)=Q(k−1)R(k−1) A(k):=R(k−1)Q(k−1)
end for
El®ször is látnunk kell, hogy az iteráció A-hoz hasonló mátrixokat készít, hiszen A(k)=R(k−1)Q(k−1)= (Q(k−1))TQ(k−1)R(k−1)Q(k−1)= (Q(k−1))TA(k−1)Q(k−1). Tehát
A(k)= (Q(k−1))T. . .(Q(0))TA Q(0). . .Q(k−1)
| {z }
=:Qk
=QTkAQk,
4.3. A sajátértékeket egyszerre közelít® eljárások 123
és így a hasonlósági transzformációk miatt A(k) sajátértékei megegyeznek A sajátértékeivel.
Minden lépésben szükség van egy QR-felbontásra. Ez elég költséges eljárás, hiszen Householder-tükrözésekkel4n3/3op a m¶veletigénye. A gyakorlatban általában úgy érdemes eljárni, hogy el®-ször a mátrixot Hessenberg-alakra transzformáljuk hasonlósági transzformáció segítségével (hogy a sajátértékek ne változzanak), majd pedig ezzel az új mátrixszal hajtjuk végre a QR-iterációt. Az iteráció során a Hessenberg-struktúra megmarad. A QR-felbontást így minden lépésben Givens-forgatásokkal hajtjuk végre (3n2 op iterációs lépésenként).
A Hessenberg alakra való transzformációt Householder-tükrözésekkel végezzük el úgy, hogy el®ször H˜1 legyen az A(2 : n,1) oszlophoz tartozó Householder-tükrözési mátrix, majd H1 = diag(1,H˜1), H˜2. Legyen A(2) =H1AH1. Ezután legyen H˜2 az A(2)(3 :n,2) oszlophoz tartozó Householder-tükrözési mátrix, majd H2 = diag(1,1,H˜2), és deniáljuk az A(3) = H2A(1)H2
mátrixot. Hasonlóan eljárva kapjukn−1lépés után az A(n)=Hn−1. . .H2H1AH1H2. . .Hn−1 Hessenberg alakú, A-val hasonló mátrixot. Az eljárást, melynek m¶veletigénye4n3/3op, sema-tikusan az alábbi módon foglalhatjuk össze:
Igazoljuk, hogy ha A egy fels® Hessenberg-mátrix, akkor a QR-iteráció megtartja ezt a struk-túrát. Legyen A=QRfels® Hessenberg-mátrix. Ekkor az
A(1)=QTAQ=QTQRQ=RQ=RQRR−1=RAR−1
mátrix is fels® mátrix, hiszen fels® mátrixok szorzata is fels® Hessenberg-mátrix lesz. Itt feltettük, hogy R, azaz A is invertálható mátrix. Különben a nulla sajátértéke lenne és deációs eljárással kisebb méret¶ mátrix vizsgálatára térhetnénk át. A továbbiakban tegyük fel, hogy az eredeti A mátrix már Hessenberg alakra van transzformálva.
Az általunk vizsgált szimmetrikus mátrixokra a Hessenberg alakra transzformált mátrix egy szimmetrikus, tridiagonális mátrix lesz. Egy ilyen mátrixszal kell ezután végrehajtani a QR-iterációt. A továbbiakban feltesszük, hogy a QR-iterációt már mindig egy Hessenberg alakra hozott mátrixszal indítjuk.
Az, hogy a Qk-val jelölt Q(0). . .Q(k−1) mátrix valóban az A-t diagonalizáló S mátrixhoz tart, egyáltalán nem nyilvánvaló az algoritmusból. Az alábbi konvergenciatételt bizonyítás nélkül közöljük.
4.3.6. tétel. (pl. QuarteroniSaccoSaleri [29], 202. oldal)
Ha az A∈Rn×n mátrixnak minden sajátértéke valós és abszolút értékben különböz® (|λ1|>
valamilyen megfelel® ˜aij konstansokkal, azaz a határértékmátrix egy fels® háromszögmátrix.
Ha A szimmetrikus, akkor{A(k)} diagonális mátrixhoz tart.
4.3.7. megjegyzés. Bár mi csak valós, szimmetrikus mátrixok sajátértékeivel foglalkozunk, a fenti tétel olyan mátrixokra is vonatkozik, melyeknek minden sajátértéke abszolút értékben kü-lönböz®.
4.3.8. megjegyzés. Igazolható, hogy hak elég nagy, akkor az A(k) mátrixok szubdiagonálbeli elemeire (a Hessenberg alak miatt csak ezek különbözhetnek nullától a f®átló alatt)
|a(k)i,i−1| ≈
Ez mutatja, hogy a QR-iteráció konvergenciarendje lineáris, és hogy a konvergenciasebességet az abszolút értékben egymást követ® sajátértékek hányadosa határozza meg. Ha vannak abszolút értékben egymáshoz közeli sajátértékek, akkor lassan konvergál csak a módszer.
4.3.9. példa. Határozzuk meg QR-iterációval aH6 Hilbert mátrix sajátértékeit! Állítsuk le akkor a QR-iterációt, ha a f®átlón kívüli elemeket tartalmazó mátrix maximumnormája már csak legfeljebb10−6-szorosa a f®átló legkisebb elem abszolút értékének! A 16. iteráció után az alábbi mátrixot nyerjük
Figyeljük meg, hogy az iteráció során abszolút értékben csökken® módon kerülnek a sajátérték-közelítések a f®átlóba. A legnagyobb abszolút érték¶ sajátértékr®l pl. a Gersgorin-tételt alkalmazva azt mondhatjuk, hogy értéke az1.61889985892434±7.2289×10−14 inter-vallumba esik (itt több tizedesjegyre pontosan írtuk ki a számokat, mint a mátrixban).
Láttuk, hogy a QR-iteráció els®rend¶ módszer sajátérték meghatározására. A módszer fel-gyorsítható az ún. eltolt QR-iteráció alkalmazásával. Az algoritmusban szerepl®µk értékek meg-választásának módjáról kés®bb lesz szó.
Eltolt QR-iteráció, A Hessenberg alakú mátrix, A(0):=A,µk adott paraméterek
fork:= 1 :kmax do
Készítsük el A(k−1)−µk−1E QR-felbontását: A(k−1)−µk−1E=Q(k−1)R(k−1) A(k):=R(k−1)Q(k−1)+µk−1E
end for
Az eljárás valóban mindig A-hoz hasonló mátrixot ad, ugyanis
A(k)=R(k−1)Q(k−1)+µk−1E= (Q(k−1))TQ(k−1)R(k−1)Q(k−1)+µk−1E
= (Q(k−1))T(A(k−1)−µk−1E)Q(k−1)+µk−1E= (Q(k−1))TA(k−1)Q(k−1)