Iterativno reˇsevanje sistemov linearnih enaˇcb

Iterativno reˇsevanje sistemov linearnih enaˇ cb

Numeriˇ cne metode, sistemi linearnih enaˇ cb

B. Jurˇ ciˇ c Zlobec

Numeriˇ cne metode FE, 2. december 2013

1 Iterativno reˇsevanje sistemov linearnih enaˇ cb Sistem linearnih enaˇ cb

Razcep matrike sistema Konvergenca

Jacobieva in Gauss-Sidlova iteracija

Jacobijeva in Gauss-Seidlova iteracija v akciji

Tridiagonalni sistemi

Vsebina

1 Iterativno reˇsevanje sistemov linearnih enaˇ cb Sistem linearnih enaˇ cb

Razcep matrike sistema Konvergenca

Jacobieva in Gauss-Sidlova iteracija

Jacobijeva in Gauss-Seidlova iteracija v akciji

Tridiagonalni sistemi

Sistem linearnih enaˇ cb z n neznankami.

a _i1 x ₁ + a _i2 x ₂ + · · · + a _in = b _i i = 1, . . . , n

V matriˇ cni obliki zapiˇsemo:

A x = b

Matrika A = (a _ij ) n×n je matrika sistema, (koeficieti pred nezanakami).

Vektor x = (x _i ) _n je vektor neznank in

Vsebina

1 Iterativno reˇsevanje sistemov linearnih enaˇ cb Sistem linearnih enaˇ cb

Razcep matrike sistema Konvergenca

Jacobieva in Gauss-Sidlova iteracija

Jacobijeva in Gauss-Seidlova iteracija v akciji

Tridiagonalni sistemi

Sploˇsno

Razcepimo matriko A na razliko matrik A = M − N.

Razcep je smiselen, ˇ ce lahko na preprost naˇ cin izraˇ cunamo inverzno matrtiko M ⁻¹ .

Zapiˇsemo sistem

(M − N) x = b M x = N x + b

x = M ⁻¹ N x + M ⁻¹ b

Iteracija

Oznaˇ cimo matriko S = M ⁻¹ N,

ki ji reˇ cemo iteracijska matrika in vektor c = M ⁻¹ b.

Sistem zapiˇsemo z novimi oznakami:

x = S x + c.

Izberemo zaˇ cetni pribliˇ zek x ⁽⁰⁾ in sproˇ zimo iteracijo

x ^(k+1) = S x ^(k) + c, k = 0, 1, . . .

Vsebina

1 Iterativno reˇsevanje sistemov linearnih enaˇ cb Sistem linearnih enaˇ cb

Razcep matrike sistema Konvergenca

Jacobieva in Gauss-Sidlova iteracija

Jacobijeva in Gauss-Seidlova iteracija v akciji

Tridiagonalni sistemi

Spektralni radij

Spektralni radij matrike je enak absolutno najveˇ cji lastni vrednosti matrike.

Spektralni radij matrike S oznaˇ cimo z ρ(S).

Spektralni radij matrike je manjˇsi ali enak normi matrike, ρ(S) ≤ kSk.

Zaporedje x ^(k) je konvergentno, ko je ρ(S) < 1.

Vsebina

1 Iterativno reˇsevanje sistemov linearnih enaˇ cb Sistem linearnih enaˇ cb

Razcep matrike sistema Konvergenca

Jacobieva in Gauss-Sidlova iteracija

Jacobijeva in Gauss-Seidlova iteracija v akciji

Tridiagonalni sistemi

Razcepimo matriko A na vsoto A = L + D + U.

Sumandi so strogo spodnjetrikotna, diagonalna in strogo zgornjetrikotna matrika matrike A.

Jacobijeva iteracija

M = D in N = −L − U Gauss-Seidlova iteracija

M = D + L in N = −U

Razcepimo matriko A na vsoto A = L + D + U.

Sumandi so strogo spodnjetrikotna, diagonalna in strogo zgornjetrikotna matrika matrike A.

Jacobijeva iteracija

M = D in N = −L − U

Gauss-Seidlova iteracija

M = D + L in N = −U

Razcepimo matriko A na vsoto A = L + D + U.

Sumandi so strogo spodnjetrikotna, diagonalna in strogo zgornjetrikotna matrika matrike A.

Jacobijeva iteracija

M = D in N = −L − U

Gauss-Seidlova iteracija

M = D + L in N = −U

Diagonalno dominantne matrike

Matrika je diagonalno dominantna, ˇ ce je absolutna vsota izvediagonalnih ˇ clenov matrike majˇsa od absolutne vrednosti diagonalnega ˇ clena.

|a _ii | > ^X

j 6=i

|a _ij |, i = 1, . . . , n

Norma matrike

kSk = sup

x

kSxk kxk .

Neskonˇ cna norma je maksimalna absolutna vrstiˇ cna vsota.

kSk _∞ = max

i



 X

j

|s _ij |



 .

Ce je matrika ˇ A diagonalno dominantna, potem je iteracijska matrika S Jacobijeve iteracije konvergentna. Njena neskonˇ cna norma je pod 1, kSk _∞ < 1.

Elementi iteracijske matrike so

s ij =

( 0, i = j

− ^a _a

ii

, i 6= j , i, j = 1, . . . , n.

Velja, (zaradi diagonalne dominantnosti matrike A):

X |s | < 1, i = 1, . . . , n.

Vsebina

1 Iterativno reˇsevanje sistemov linearnih enaˇ cb Sistem linearnih enaˇ cb

Razcep matrike sistema Konvergenca

Jacobieva in Gauss-Sidlova iteracija

Jacobijeva in Gauss-Seidlova iteracija v akciji

Tridiagonalni sistemi

Jacobijeva iteracija

a _ii x _i = − ^X

j 6=i

a _ij x _j + b _i , i = 1, . . . , n

x _i ^(k+1) = − ^X

j6=i

a _ij

a _ii x _j ^(k) + b _i a _ii

x _i ^(k+1) = ^X

j

s _ij x _j ^(k) + c _i

Gauss-Seidlova iteracija a _ii x _i + ^X

j <i

a _ij x _j = − ^X

j>i

a _ij x _j + b _i , i = 1, . . . , n

a _ii x _i ^(k+1) + ^X

j <i

a _ij x _j ^(k+1) = − ^X

j >i

a _ij x _j ^(k) + b _i , i = 1, . . . , n

x _i ^(k+1) = − ^X

j <i

a _ij a ii

x _j ^(k+1) − ^X

j >i

a _ij a ii

x _j ^(k) + b _i

a ii

Program

Dana je diagonalno dominantna matrika A in vektor b.

Reˇsi sistem Ax = b z Jacobijevo in Gauss-Seidlovo iteracijo.

Priprava

d=diag(A);S=-A;c=b./d;n=length(d);

eps=1e-5;m=1000;

for i=1:n

S(i,i)=0; S(i,:)=S(i,:)/d(i);

Jacobijeva iteracija

x=zeros(size(d));

for i=1:m xx=x;

x=S*x+c;

if abs(x-xx)<eps, break, end;

end;

Gauss-Seidlova iteracija

x=zeros(size(b));

for i=1:m xx=x;

for j=1:n

x(j)=S(j,:)*x+c(j);

end;

if abs(xx-x)<eps, break, end;

end;

Vsebina

1 Iterativno reˇsevanje sistemov linearnih enaˇ cb Sistem linearnih enaˇ cb

Razcep matrike sistema Konvergenca

Jacobieva in Gauss-Sidlova iteracija

Jacobijeva in Gauss-Seidlova iteracija v akciji

Tridiagonalni sistemi

Razcep po Jacobiju

d ₁ x ₁ = −u ₁ x ₂ + b ₁

d _i x _i = −l _i−1 x _i −1 − u _i x _i+1 + b _i i = 2, . . . n − 1 d n x n = −l _n−1 x n−1 + b n

Jacobijeva iteracija

x ₁ ^(k+1) = −u ₁ /d ₁ x ₂ ^(k) + b ₁ /d ₁

x _i ^(k+1) = −l _i−1 /d _i x _i−1 ^(k) − u _i /d _i x _i+1 ^(k) + b _i /d _i i = 2, . . . n − 1

(k+1) (k)

Razcep po Gauss-Seidlu

d ₁ x ₁ = −u ₁ x ₂ + b ₁

d _i x _i + l i−1 x i−1 = −u _i x _i+1 + b _i i = 2, . . . n − 1 d n x n + l n−1 x n−1 = b n

Gaus-Seidlova iteracija

x ₁ ^{(k +1)} = −u ₁ /d ₁ x ₂ ^(k) + b ₁ /d ₁

x _i ^{(k +1)} = −l _i−1 /d _i x _i−1 ^(k+1) − u _i /d _i x _i+1 ^(k) + b _i /d _i i = 2, . . . n − 1

x _n ^{(k +1)} = −l n−1 /d n x _n−1 ^{(k +1)} + b n /d n

Program

Podatki

n=5;

d=1.91*ones(1,n);

u=ones(1,n-1);

l=ones(1,n-1);

b=2+sin(1:n);

eps=1e-5; m=1000;

Priprava

Jacobiejeva iteracija

xx=zeros(size(c));

for k=1:m, x=xx;

xx(1)=-u(1)*x(2)+c(1);

for i=2:n-1,

xx(i)=-l(i-1)*x(i-1)-u(i)*x(i+1)+c(i);

end;

xx(n)=-l(n-1)*x(n-1)+c(n);

if abs(x-xx)<eps, break, x=xx; end;

end;

Gaus-Seidlova iteracija

x=zeros(size(c));

for k=1:m, xx=x;

x(1)=-u(1)*x(2)+c(1);

for i=2:n-1,

x(i)=-l(i-1)*x(i-1)-u(i)*x(i+1)+c(i);

end;

x(n)=-l(n-1)*x(n-1)+c(n);

if abs(x-xx)<eps, break, end;

end;