Iterační metody řešení soustav lineárních rovnic Matice je: n
● ●
diagonálně dominantní právě tehdy, když ∣a ii∣
∑
j=1, j≠i
∣aij∣
pozitivně definitní (symetrická matice) právě tehdy, když pro ∀ x ≠0 platí x , A x 0 .
Tyto vlastnosti budou důležité pro zaručení konvergence iteračních metod. Postup: ●
zvolí odhad řešení x 0 (pokud odhad neznáme tak libovolné číslo)
●
k1 =B k x k ck předpokládáme x
●
●
pro správné řešení musí platit x =B k x ck , tedy k1 −x =B k x k −x =B k B k−1 x k−1−x =B k B0 x 0−x požadujeme, aby pro x rostoucí k šlo k 0 B k B0 =0 pro konvergenci je tedy nutné a stačí, aby lim k ∞
Metody: ●
nestacionární – matice B v každém kroku jiná
●
stacionární – matice B pořád stejná
Nestacionární metoda – Řízená relaxace Řešíme soustavu A x =b . Na počátku odhadneme řešení libovolně x 0 . Spočítáme reziduum 0 A x −b a určíme jeho v absolutní hodnotě maximální složku. Potom volíme matici B tak, abychom tuto složku vynulovali. Například pro soustavu 3 ×3 , kde maximální je druhá složka, 0 0 je to a 21 x 0 1 a 22 x 2 a 23 x 3 −b 2 . 1 1 Požadujeme tedy, aby v následujícím kroku platilo a 21 x 1 1 a 22 x 2 a 23 x 3 −b 2 =0 . Z této 1 1 1 1 rovnice vyjádříme itou, tedy v tomto případě druhou složku x 2 = b 2 −a 21 x 1 −a 23 x 3 . a 22 Matici B a vektor c volíme tak, aby všechny složky až na itou zůstaly a pro itou složku aby platil předcházející vztah, tedy v našem případě
x 1 2 =
1 a 21
0
1 0 b 2 −a 21 x 0 1 −a 23 x 3 . To odpovídá volbě B= − a 22 a 22 0
1 0 0 1 ⋮ ⋮ ai2 obecném případě tedy B= ai1 − − aii a ii ⋮ ⋮ 0 0
⋮ ⋮ ⋮ 0 − ⋮ ⋮ ⋮
0 0 ⋮ ai n aii ⋮ 1
0 − 0
0 a 23 a 22 1
0 0 ⋮ a c = bi aii ⋮ 0
a c =
0 b2
a 22 0
. V
.
Tato metoda je výpočetně náročná, neboť je v každém kroku nutné hledat maximum rezidua. Proto se v praxi nepoužívá. Stacionární metody (Prostá iterace, Jakobiho, Gauss-Seidlova a Superrelaxační metoda) ●
Konvergence stacionárních metod
Vlastní číslo matice B je číslo takové, že k němu existuje nenulový vlastní vektor, pro který platí B x = x . Nutná a postačující podmínka konvergence stacionárních iteračních metod je, že spektrální poloměr matice B musí být menší než 1 . Jinými slovy pro všechna vlastní čísla i matice B musí platit B=max i=1,... , n∣i∣1 . Postačující podmínka konvergence metody pak je, že v některé maticové normě platí ∥B∥1 . (Z této podmínky samozřejmě plyne podmínka předcházející, neboť ∣∣∥x∥=∥ x∥=∥B x∥≤∥B∥⋅∥x∥ , tedy ∣∣≤∥B∥1 .) ¿Proč požadujeme max i=1,... , n∣i∣1 ? Nějaká věta říká, že pro matici B platí následující k ekvivalence lim B = 0 ⇔ B1 . Pro konvergenci jsme již dříve ukázaly podmínku k ∞
lim B k B0 = 0 , tedy ve stacionárním případě lim B k B0 =lim B k = 0 . k ∞ k ∞
k ∞
Přesnost řešení můžeme odhadnout jako k k −1 −1 k −1 k ∥x −x∥=∥x − A b∥=∥A A x −b∥≤∥A ∥⋅∥A x −b∥ . ●
Prostá iterace
Rovnici A x =b můžeme přepsat jako x = I− Ax b .
Prostá iterace je založena právě na tomto vzorci, tedy x k1= I− Ax kb . Jako matice B je tedy volena matice I− A a jako vektor c je volen vektor b . Přesnost můžeme odhadnout následovně: k k−1 b= x =B x
=B B x k−2b b= =B 2 x k−2B bb== =B k x 0B k−1 bB k−2 bb 2 m −1 Platí dále implikace : jestliže B1 , pak lim IBB B = I−B . m∞
(v podstatě jde o geometrickou posloupnost matic) Pro řešení x pak z předchozího platí x = A −1 b= I−B−1 b= IBB 2 b . ∥x k −x∥=∥B k x 0 IBB k−1 b− IBB 2 b∥= =∥B k x 0− B k B k1 b∥= =∥B k x 0− IBB 2 b∥≤ ≤∥B k∥⋅∥x 0− IBB 2 b∥≤ ≤∥B∥k⋅∥x 0− IBB 2 b∥≤
[
]
≤∥B∥k⋅∥x 0∥∥I∥∥B∥∥B 2∥∥b∥ ≤ k
[
0
≤∥B∥ ⋅ ∥x ∥
∥b∥ 1 −∥B∥
]
Metoda prosté iterace se pro systémy lineárních rovnic prakticky nepoužívá.
●
Jacobiho metoda
V čem metoda spočívá?
n
∑ aij x j −bi=0
Hledáme řešení rovnice A x =b . Napíšeme itou rovnici z této soustavy 1 této rovnice vyjádříme x i jako x i = aii
x i k1=
1 a ii
∑ n
j=1, j≠i
∑ n
aij x j −b i
aij x kj −b i
j=1, j≠i
j=1
a takto provádíme iterace
.
Jinak se to dá pro matici soustavy A (řádu 3) formulovat takto. Matice A se rozloží na součet 3 matic podle následujícího vzoru.
a 11 a 12 a13 d 11 0 0 0 r 12 r 13 0 0 0 A= a 21 a 22 a 23 = DLR= 0 d 22 0 l 21 0 0 0 0 r 23 l 31 l 32 0 a 31 a 32 a 33 0 0 d 33 0 0 0
Z těchto tří matic se pak vytvoří matice B
1 d 11
−
B= D −1 LR=
0 0
0
−
1 d 22
0
0
0
−
a12
−
a13
a11 a11 0 r 12 r 13 a a 23 0 × l 21 0 r 23 = − 21 0 − a 22 a 22 l 31 l 32 0 1 a 31 a 32 − − − 0 d 33 a 33 a 33
a vektor c se definuje jako
b1
−1 c = D b=
a11 b2
a 22 b3
.
a 33
Iterace se provádí podle obecného vzorce x k1=B x kc . Už z formulace úlohy je patrné, že matice A nesmí mít nulové diagonální prvky.
. Z
Pro Jakobiho metodu je konvergence zaručená pro matice soustavy A diagonálně dominantní: n
n
Diagonální dominance je definována ∣aii∣
∑
j=1, j≠i
∣a ij∣ . Z toho plyne
j=1, j≠i
maximovou normu matice B= D −1 L R platí ∥B∥ 1 .
●
∑
∣∣ a ij
aii
1 a tedy pro
Metoda GaussSeidlova
Téměř stejná, jako Jacobiho, jen vždy používáme už všechny nově spočítané složky vektoru x . Konvergence pro poměrně širokou třídu matic A , ale může být velmi pomalá. Porovnání s předchozími metodami: GaussSeidl
x k1 =x ki i
1 k1 −ai1 x k1 −−a i i−1 x i−1 −a i i x ki −−a i n x kn b i 1 aii
k1 =x k − DL −1 A x k −b x
B=− DL
−1
R
−1 c = DL b
konvergence zaručena A diagonálně dominantní nebo pozitivně definitní přímá iterace
k k x k1 =x i k−ai1 x k i 1 −a i2 x 2 −−a i n x n b i k1 =x k − A x k −b x
B= I− A
Jacobiho metoda
x k1 =x ki i
c =b 1 −ai1 xk1 −ai2 x k2 −−ai n xnkbi aii
k1 =x k − D −1 A x k −b x
B=−D
−1
L R
c = D
−1
b
konvergence zaručena A diagonálně dominantní
●
Superrelaxační metoda
Vlastně urychlení konvergence GaussSeidlovy metody. V GaussSeidlově metodě označíme x k =x k1−x k . V GaussSeidlově metodě tedy platí x k1=x k x k . Superrelaxační metoda pak tento vztah přeformuluje takto x k1=x k x k , kde je z intervalu 0,2 ( 1 subrelaxační, 1 superrelaxační, = 1 GaussSeidl). Relaxační faktor slouží tedy k urychlení konvergence a jeho optimální hodnotu lze vypočítat 2 −1 podle vztahu opt = , kde BGS =− DL R je matice B z Gauss 2 1 1 − BGS Seidlovy metody.
(Upozornění: V některých materiálech se matice L a R násobí −1 při odvozování Gauss Seidlovy a Superrelaxační metody.) Příklad na Superrelaxační metodu v PASCALU DEMRELAX.PAS . Matice MATREL1.DAT , MATREL.DAT .
Problém vlastních čísel Vlastní číslo : , ∃x ≠0 , že A x = x , vlastní číslo, x je vlastní vektor A . Řeší se problémy: úplný nebo částečný Vlastní čísle lze počítat z charakteristického polynomu matice A det A− I =0 . Pro matici n×n se jedná o polynom ntého stupně, který má tedy n kořenů (mohou být i násobné, komplexně sdružené). Ke každému vlastnímu číslu existuje alespoň jeden vlastní vektor. Matice je defektní – má méně než n LN vlastních vektorů Normální matice (splňuje vztah A AT = AT A ) má n LN vlastních vektorů. Symetrická matice ( A= AT ) má všechna vlastní čísla reálná. Trojúhelníková matice má vlastní čísla na diagonále. Matice podobné (matice A a P −1 A P ) mají stejná vlastní čísla, neboť −1 −1 −1 det P A P− I =det P A− I P =det P det A− I det P =det A− I . Ke každé matici existuje matice jí podobná v Jordanově normálním tvaru. K normálním maticím existuje dokonce podobná diagonální matice. Numerické metody řešení úplného problému– travsformace na přibližně diagonální nebo trojúhelníkový tvar. ●
Jacobiho transformace
Najde všechny vlastní čísla symetrické matice. Cílem metody je podobnostními transformacemi postupné zmenšování mimodiagonálních prvků. Postup: Nalezení v absolutní hodnotě maximálního mimodiagonálního prvku a pq .
Podobnostní transformace – pootočení os, aby se čtvercová matice
a pp a pq a qp a qq
prvků
1 . . . . . 0 . ⋱ . . . . . . . c −s . . matice A stala diagonální. Matice transformace je T pq = . . ⋮ 1 ⋮ . . . . s c . . . . . . . ⋱ . 0 . . . . . 1
.
Pro koeficienty c , s platí c=cos , s=sin . Po kté iteraci je prvek matice A akpq=c 2 −s 2 apqk−1−cs appk−1−aqqk−1 a abychom ho nulovaly, musíme volit 2 a pq úhel jako tak, že tan 2 = . a pp −a qq Příklad v PASCALU ZKJACOBI.PAS , matice JAC1.DAT , JAC2.DAT .
●
LU rozklad Pomalá konvergence, nutný velký počet operací. Provádí se rozklad A k =L k U k a dále se volí A k1=U k L k . Matice U k se volé jako T ortogonální, tedy U k U k =I . Z toho plyne, že matice A k a A k1 jsou si podobné. A k1 =U k L k =U k L k U k U Tk =L k U k = A k
●
Částečný problém vlastních čísel Hledá se nejčastěji extrémní vlastní číslo (v absolutní hodnotě největší, nejmenší). Například metoda: k1 = x
1 A x k , kde k =e T1 A x k ( k je první složka A x k ) k
k = a lim x k= x1 platí lim k ∞ k ∞
Nejmenší vlastní čísla matice A jsou největší matice A −1 . Druhé nejmenší vlastní číslo – redukce matice na řád n−1 . Vlastní číslo v nějaké oblasti – posun A I x =x . Příklad v PASCALU CASTP.PAS , matice např. MATREL.DAT , MATREL1.DAT ...