Metody rozwiązywania układów równań
liniowych
n
n
nn
n
n
n
n
n
n
b
x
a
x
a
x
a
b
x
a
x
a
x
a
b
x
a
x
a
x
a
2
2
1
1
2
2
2
22
1
21
1
1
2
12
1
11
b
Ax
n
nn
n
n
n
n
b
b
b
a
a
a
a
a
a
a
a
a
2
1
2
1
2
22
21
1
12
11
b
A
0
det
A
Parę słów o macierzach
Macierz m
n: tablica m na n (m wierszy n kolumn) liczb (np.
tabliczka mnożenia).
Macierz kwadradowa: m=n
Macierz symetryczna (zawsze kwadratowa): a
ij
=a
ji
Macierz transponowana A
T
: (A
T
)
ij
=a
ji
Macierz nieosobliwa: macierz o niezerowym wyznaczniku.
Macierz dodatnio określona:
x
T
Ax>0 dla każdego niezerowego wektora x.
Norma euklidesowa macierzy:
Norma spektralna macierzy
Wskaźnik uwarunkowania macierzy
n
i
n
j
ij
a
1
1
2
A
i
i
max
A
1
cond
A
A
A
Metody skończone:
• Metoda Gaussa
• Metoda Gaussa-Jordana
• Metody Choleskiego
• Metoda Householdera
• Metoda sprzężonych gradientów
Metody iteracyjne dla dużych układów równań:
• Metoda Jacobiego
• Metoda Gaussa-Seidla
n
n
nn
n
n
n
n
c
x
r
c
x
r
x
r
c
x
r
x
r
x
r
2
2
2
22
1
1
2
12
1
11
Metoda eliminacji Gaussa z wyborem elementu głównego w
kolumnie
Układ równań sprowadzamy do postaci trójkątnej
Układ z macierzą trójkątną można następnie łatwo rozwiązać
zaczynając od obliczenia wartości x
n
z n-tego równania,
następnie wstawić x
n
do równania n-1 i wyliczyć z niego x
n-1
,
następnie wstawić x
n
oraz x
n-1
do równania n-2 i wyliczyć x
n-2
aż
do dotarcia do równania pierwszego i wyznaczenia x
1
.
1. Wybieramy równanie i takie, że |a
i1
| jest największym elementem w
pierwszej kolumnie po czym przestawiamy i-te równanie na początek i
eliminujemy x
1
z równań od 2 do n.
)
1
(
)
1
(
2
2
)
1
(
2
)
1
(
2
)
1
(
2
2
)
1
(
22
)
0
(
1
)
0
(
1
2
)
0
(
12
1
)
0
(
11
n
n
n
n
n
n
n
n
b
x
a
x
a
b
x
a
x
a
b
x
a
x
a
x
a
1
1
,
1
)
0
(
11
)
0
(
1
)
0
(
1
)
0
(
)
1
(
)
0
(
11
)
0
(
1
)
0
(
1
)
0
(
)
1
(
i
a
a
b
b
b
k
i
a
a
a
a
a
i
i
i
i
k
ik
ik
2. Procedurę powtarzamy z macierzą A
(1)
o rozmiarach (n-1)x(n-1) i
wektorem b
(1)
o rozmiarze n-1, eliminując z nich drugą zmienną i
otrzymując macierz A
(2)
o rozmiarach (n-2)x(n-2) i wektor b
(2)
o
rozmiarze n-2. W ten sam sposób postępujemy z kolejnymi macierzami
A
(2)
, A
(3)
,..., A
(n-1)
oraz wektorami b
(2)
, b
(3)
,..., b
(n-1)
.
j
i
a
a
b
b
b
j
k
j
i
a
a
a
a
a
j
jj
j
ij
j
j
j
i
j
i
j
jj
j
ij
j
jk
j
ik
j
ik
)
1
(
)
1
(
)
1
(
)
1
(
)
(
)
1
(
)
1
(
)
1
(
)
1
(
)
(
,
Dla j-tego kroku
Po zakończeniu operacji otrzymujemy układ równań z macierzą
trójkątną
)
1
(
)
1
(
22
)
0
(
11
)
1
(
)
1
(
)
2
(
3
)
0
(
3
3
)
2
(
33
)
1
(
2
)
0
(
2
3
)
1
(
23
2
)
1
(
22
)
0
(
1
)
0
(
1
3
)
0
(
13
2
)
0
(
12
1
)
0
(
11
)
1
(
det
n
nn
p
n
n
n
n
nn
n
n
n
n
n
n
a
a
a
b
x
a
b
x
a
x
a
b
x
a
x
a
x
a
b
x
a
x
a
x
a
x
a
A
p jest liczbą przestawień wierszy macierzy A podczas
sprowadzania układu równań do postaci trójkątnej.
1
,...,
2
,
1
1
)
1
(
)
1
(
)
1
(
)
1
(
)
1
(
)
1
(
n
n
j
x
a
a
a
b
x
a
b
x
k
n
j
k
j
jj
j
jk
j
jj
j
j
j
n
nn
n
n
n
3. Z otrzymanego układu równań z macierzą trójkątną
wyznaczamy po kolei x
n
, x
n-1
,..., x
1
.
Wysiłek obliczeniowy (liczba mnożeń i dzieleń) w metodzie
eliminacji Gaussa:
Faktoryzacja macierzy A: n(n
2
-1)/3 operacji
Przekształcenie wektora b: n(n-1)/2 operacji
Obliczenie x: n(n+1)/2 operacji.
Razem: n
3
/3+n
2
-n/3≈n
3
/3 operacji.
Zastosowania metody eliminacji Gaussa i innych metod
rozwiązywania układów równań liniowych.
1. Odwracanie macierzy
0
0
1
0
0
1
,
,
1
,
1
in
i
i
i
i
i
i
i
1
1
1
1
1
1
A
A
A
A
A
A
I
AA
Należy więc jednocześnie rozwiązać n układów równań z wektorami
jednostkowymi poszczególnych osi w przestrzeni n-wymiarowej jako
wyrazami wolnymi. Wszystkie wektory wyrazów wolnych
przekształca się jednocześnie w czasie faktoryzacji macierzy A.
Koszt obliczeniowy operacji wynosi (4/3)n
3
-n/3 operacji.
2. Mnożenie macierzy przez macierz odwrotną
B
AC
B
A
C
1
Jeżeli macierz A ma rozmiar nxn to macierz B musi mieć n
wierszy (liczba kolumn może być dowolna, np. m).
Zagadnienie sprowadza się do jednoczesnego rozwiązania m
układów równań liniowych, których wektory wyrazów wolnych
są równe kolejnym kolumnom macierzy B. Podobie jak przy
odwracaniu macierzy podczas faktoryzacji A jednocześnie
przekształcamy wszystkie kolumny macierzy B.
Wysiłek obliczeniowy wynosi n
3
/3-n/3+mn
2
operacji.
Metoda Gaussa-Jordana
W odróżnieniu od metody eliminacji Gaussa eliminujemy
poszczególne zmienne nie tylko z równań następujących po danym
równaniu ale również z równań je poprzedzających. W wyniku tej
operacji końcowy układ równań ma postać diagonalną a w kroku j-
1 ma on następującą postać:
)
(
)
(
1
)
(
1
,
)
(
2
)
(
1
)
(
1
,
2
)
(
)
(
1
)
(
1
1
)
(
1
,
1
1
)
0
(
11
j
n
n
j
nn
j
j
j
n
j
n
j
jn
j
j
j
j
j
jj
j
n
j
n
j
j
j
b
x
a
x
a
b
x
a
x
a
x
a
b
x
a
x
a
x
a
)
1
(
)
1
(
)
1
(
)
1
(
)
1
(
)
1
(
)
(
)
1
(
)
1
(
)
1
(
)
1
(
)
(
,...,
1
,
1
,...,
2
,
1
,...,
1
,
1
,...,
2
,
1
,...,
2
,
1
0
i
ii
n
i
i
j
jj
j
ij
j
j
j
i
j
i
j
jj
j
ij
j
jk
j
ik
j
ik
a
b
x
n
j
j
i
a
a
b
b
b
j
i
j
k
a
a
a
a
n
j
j
i
j
k
a
Dla dużych n metoda Gaussa-Jordana wymaga ok. n
3
/2
operacji.
Metody typu Choleskiego dla macierzy symetrycznych silnie
nieosobliwych
T
LDL
A
L
T
L
D
L
T
LL
A
klasyczna metoda Choleskiego
tylko dla macierzy dodatnio
określonych.
Postępowanie przy rozwiązywaniu układów równań metodą
faktoryzacji Choleskiego.
1. Wyznaczenie faktorów L i D. Układ przyjmuje postać
LDL
T
x=b
2. Obliczenie pomocniczego wektora w.
w=L
-1
b przez rozwiązanie układu równań Lw=b.
Ponieważ L jest macierzą trójkątną dolną układ ten rozwiązuje się
wyliczając kolejno w
1
, w
2
,…, w
n
podobnie jak w koncowym etapie
eliminacji Gaussa.
3. Obliczenie z=D
-1
w (D jest macierzą diagonalną więc po prostu
dzielimy w
i
przez d
ii
. Ten etap nie występuje w klasycznej metodzie
Choleskiego.
4. Obliczenie x poprzez rozwiązanie układu równań z macierzą
trójkątną górną
L
T
x=z
Ten etap jest identyczny z ostatnim etapem metody eliminacji
Gaussa.
Metoda wymaga ok. n
3
/6 operacji (2 razy mniej niż metoda eliminacji
Gaussa). Uwaga: klasyczna metoda Choleskiego wymaga ponadto
n pierwiastkowań.
Klasyczna faktoryzacja Choleskiego (A=LL
T
)
1
1
1
1
2
11
1
1
11
11
1
i
k
jk
ik
ij
ii
ji
i
k
jk
ii
ii
j
j
l
l
a
l
l
l
a
l
l
a
l
a
l
Faktoryzacja “bezpierwiastkowa”
1
0
1
0
0
1
1
,
1
21
2
1
n
n
n
n
l
l
l
d
d
d
L
D
i
n
k
ik
jk
k
ji
ji
j
j
n
k
nk
k
nn
n
i
k
ki
k
ii
i
d
l
l
d
a
l
d
a
l
l
d
a
d
l
d
a
d
a
d
1
1
1
1
1
1
1
2
1
1
2
11
1
/
Metoda Householdera
)
2
(
)
2
(
1
)
2
(
12
)
2
(
11
)
1
(
)
2
(
)
1
(
1
)
1
(
21
)
1
(
1
)
1
(
11
)
1
(
11
1
1
1
2
1
)
1
(
)
1
(
)
2
(
)
1
(
~
sgn
2
,
A
A
H
A
a
v
v
v
v
I
H
c
b
Q
Ax
Q
cRz
b
Q
Ax
Q
Rz
Q
Q
H
H
H
Q
QR
A
T
T
T
T
T
1
n
n
T
n
a
a
a
a
a
a
a
0
~
0
sgn
2
)
1
(
)
1
(
)
1
(
)
2
(
2
)
2
(
22
)
2
(
1
)
2
(
12
)
2
(
11
)
(
)
1
(
)
1
(
)
(
)
(
,
1
)
(
)
(
)
(
2
1
)
(
i
n
in
i
ii
n
n
i
i
i
i
ni
i
i
i
i
i
i
ii
i
ii
i
T
i
i
i
i
n
i
a
a
a
a
a
a
a
a
a
a
a
A
A
H
A
a
v
v
v
v
I
H
Iteracyjna poprawa rozwiązania układu równań liniowych
)
1
(
)
(
)
1
(
)
(
)
(
)
1
(
)
1
(
)
0
(
)
1
(
)
0
(
)
0
(
)
1
(
p
p
p
p
p
p
z
x
x
r
Ax
b
Az
z
x
x
r
Ax
b
Az
Proces kończymy gdy
)
(
)
1
(
max
max
p
i
i
p
i
i
x
z
Metoda sprzężonych gradientów
Rozwiązanie układu równań przedstawiamy jako rozwiązanie
zagadnienia minimalizacji funkcji E(x), której gradient jest różnicą
pomiędzy lewą a prawą stroną układu równań.
b
Ax
x
b
Ax
A
b
Ax
x
1
T
)
(
grad
2
1
)
(
E
E
Metoda polega na wykonaniu n kroków w kierunkach d
0
, d
1
,…, d
n-1
stopniowo zbliżających się do kierunku jednej z osi
hiperparaboloidy E(x); długość kroku w każdym kierunku jest tak
dobrana aby zlokalizować minimum E na tym właśnie kierunku.
Poszczególne kierunki są A-ortogonalne, tj.
d
k+1
Ad
k
=0
1. Wybieramy punkt startowy x
0
. Obliczamy d
0
=-g
0
=b-Ax
0
. Jeżeli
g
0
=0 procedura jest zakończona.
2. W kroku k=0,1,…,n-1 obliczamy
k
T
k
k
T
k
Ad
d
g
d
k
Wstawiamy x
k+1
=x
k
+
k
d
k
Wstawiamy g
k+1
=g
k
+
k
Ad
k
. Jeżeli ||g
k+|
||< kończymy proces,
w przeciwnym wypadku wstawiamy
k
T
k
k
T
1
k
Ad
d
Ad
g
k
a następnie wstawiamy d
k+1
=-g
k+1
+
k
d
k
Procedura jest teoretycznie zawsze zbieżna w n krokach (0,1,
…,n-1). Dla dobrze uwarunkowanego zadania zbieżność jest
szybsza a dla źle uwarunkowanego zadania często występują
oscylacje.
Metody iteracji Jacobiego oraz iteracji Gaussa-Seidla
Metody te stosuje się dla dużych układów równań, pojawiających
się przykładowo przy rozwiązywaniu równań różniczkowych
cząstkowych lub równań całkowych.
Iteracja Jacobiego
1. Wyznaczamy przybliżenie początkowe x
(0)
2. Kolejne przybliżenia x
(1)
, x
(2)
,…, obliczamy wyznaczając x
i
z i-tego
równania dla i=1,2,…,n a następnie podstawiając po prawej
stronie poprzednie przybliżenie x.
ii
i
n
i
k
k
p
k
ii
ik
p
i
a
b
x
a
a
x
1
)
(
)
1
(
3. Proces kończymy jeżeli ||x
(p+1)
-x
(p)
||<.
Metoda iteracji Gaussa-Seidla.
Metoda ta różni się od metody iteracji Jacobiego tym, że kolejne
przybliżenie i-tej współrzędnej wektora x obliczamy wykorzystując
przybliżenia x
1
,…,x
i-1
wyliczone w aktualnej iteracji oraz przybliżenia
x
i+1
,…,x
n
z poprzedniej iteracji.
ii
i
i
ii
ik
r
ik
ii
ik
l
ik
i
n
i
k
p
k
r
ik
i
k
p
k
l
ik
ii
i
n
i
k
p
k
ii
ik
i
k
p
k
ii
ik
p
i
a
b
c
i
k
a
a
i
k
b
i
k
i
k
a
a
b
c
x
b
x
b
a
b
x
a
a
x
a
a
x
0
0
)
(
)
(
1
)
(
)
(
1
1
)
1
(
)
(
1
)
(
1
1
)
1
(
)
1
(
Metoda iteracji Gaussa-Seidla jest zawsze zbieżna dla macierzy
dodatnio określonych.
)
1
(
)
(
)
(
)
1
(
2
)
(
)
(
)
1
(
2
1
)
(
1
1
)
(
1
1
)
1
(
)
(
)
1
(
max
1
1
1
1
2
1
1
2
p
k
p
k
p
x
p
k
k
p
p
p
opt
ii
i
p
i
n
k
p
i
ii
ik
i
k
p
i
ii
ik
p
i
p
i
x
x
x
x
q
q
q
a
b
x
x
a
a
x
a
a
x
x
Kod źródłowy metody Gaussa-Seidla.
Optymalna wartość czynnika
relaksacji;
1
jest największą
wartością własną macierzy
B=B
(l)
+B
(r)
Wartość
czynnika
relaksacji łatwa
do obliczenia.