Metoda relaxarii succesive

Metoda relaxarii succesive

1 Breviar teoretic

Metoda relaxarii succesive este o metoda de rezolvare numerica a sistemelor de tip Cramer, prin aproximatii succesive. Aceasta metoda se deosebeste de metoda Gauss-Seidel prin aceea ca se introduc corectiile

Δx^(k) = x^(k+1) x^(k) , k = 0, 1, 2, (2.39)

Matricea A a sistemului se descompune in suma L + D + U, unde L este matrice triunghiulara subdiagonala, D matrice diagonala si U matrice triunghiulara supradiagonala.

Pentru un vector x⁽⁰⁾ Rⁿ, sirul de vectori x^(k) definit prin:

k = 0, 1, 2, . . . (2.40)

se numeste traiectoria vectorului x⁽⁰⁾ obtinuta prin relaxari succesive.

Traiectoria vectorului x⁽⁰⁾ obtinuta prin relaxari succesive converge daca si numai daca raza spectrala ρ a matricei

(2.41)

este strict subunitara.

In caz de convergenta, componentele ale vectorului situat pe traiectoria vectorului x⁽⁰⁾ obtinuta prin relaxari succesive sunt date de relatiile:

(2.42)

(2.43)

Observatie

Metoda Gauss-Seidel este un caz particular al metodei relaxarii succesive, pentru care ω=1.

2 Problema rezolvata

Sa se gaseasca primele 3 elemente ale traiectoriei vectorului (0, 0)^T folosind metoda relaxarii succesive cu ω = 0.5, pentru sistemul:

Rezolvare

Sistemul se poate scrie sub forma Ax = b, unde

.

Matricea A se descompune in suma L + D + U, cu

Algoritmul converge daca raza spectrala a matricei

este strict subunitara. Efectuand calculele, obtinem

(M) = 0.75 < 1

si deci algoritmul este convergent.

In continuare, aplicam formulele (2.42)-(2.43), plecand de la punctele x⁽⁰⁾ = 0, y ⁽⁰⁾ = 0, si obtinem:

x y⁽¹⁾ =

x y⁽²⁾ =

x y⁽³⁾ =

Pentru comparatie, determinam solutia exacta a sistemului considerat:

x y =

3 Implementare

A. Algoritm

Observatie. Deoarece metoda relaxarii succesive este o metoda iterativa, trebuie specificata o conditie de oprire a algoritmului. Aceasta conditie este similara cu cea folosita in paragrafele anterioare, si anume:

(2.44)

Date de intrare: un sistem de ecuatii (o multime de ecuatii), un punct initial, x⁽⁰⁾, o relaxare, ω, o eroare ε.

Date de iesire: solutia aproximativa a sistemului, obtinuta in urma aplicarii traiectoriei vectorului x⁽⁰⁾ obtinuta prin relaxari succesive pana cand este indeplinita conditia (2.44).

Algoritmul consta in urmatoarele etape:

1. generarea matricei A a sistemului si a vectorului b

n - numarul de necunoscute (numarul de linii ale matricei A)

2. generarea matricelor L, D, U si verificarea convergentei metodei:

si 0 < ω < 2

3. constructia traiectoriei relaxarilor succesive

repeta

(2.45)

pana cand

B. Programe MAPLE si rezultate

relsucc := proc(eq::set(equation), init::vector, rlx::numeric, eps::float)

local var, n, AA, A, b, l, d, u, i, j, m, lst, xo, test, k, x;

var:=[op(indets(eq))];

n:=nops(var);

if vectdim(init)<>n then

ERROR('numarul de necunoscute nu este egal cu dimensiunea vectorului initial')

fi;

if rlx<=0 or rlx>=2 then

ERROR('parametrul rlx trebuie sa fie in intervalul (0,2)')

fi;

AA:=genmatrix(eq, var, flag);

A:=delcols(AA, n+1..n+1);

b:=col(AA, n+1);

l:=matrix(n, n, 0):

u:=matrix(n, n, 0):

d:=matrix(n, n, 0):

for i from 1 to n do

for j from 1 to i-1 do

l[i, j]:=A[i, j];

od;

d[i, i]:=A[i, i];

for j from i+1 to n do

u[i, j]:=A[i, j];

od;

od;

# conditia de convergenta

m:=multiply( inverse(matadd(l, d, 1, 1/rlx)), matadd(d, u, -(1-1/rlx), -1));

lst:=[eigenvals(m)];

if evalf(max(seq(abs(lst[k]), k=1..nops(lst))))>=1 then

ERROR('Algoritmul nu converge');

fi;

# algoritmul propriu-zis

for i from 1 to n do

xo[i]:=init[i]

od;

test:=1;

while test>=evalf(eps*sqrt(n)) do

x[1]:=evalf((1-rlx)*xo[1] + rlx/A[1,1]*( b[1]-sum(A[1,k]*xo[k],k=2..n)));

for i from 2 to n do

x[i]:=evalf((1-rlx)*xo[i]+1/A[i,i]*( b[i]-sum(A[i,k]*x[k],k=1..i-1)-sum(A[i,k]*xo[k],

k=i+1..n)));

od;

test:=evalf(sqrt( sum( (x[k]-xo[k])^2, k=1..n ) ));

for i from 1 to n do

xo[i]:=x[i];

od;

od;

RETURN(seq(var[i]=x[i],i=1..n));

end:

debug(relsucc):

relsucc(, vector(2, [0,0]), 0.99, 0.01);

, array(1 ..2,[(1)=0,(2)=0]), .99, .1e-1

var := [x, y]

n

b := [5, 5]

d u_1,2 := 1 l_2,1 := 1 d_2,2 := 2

xo xo₂ := 0

test

x x₂ := 1.675000000

test

xo xo₂ := 1.675000000

x x₂ := 1.959875000

test

xo xo₂ := 1.959875000

x x₂ := 2.012409375

test

xo xo₂ := 2.012409375

x x₂ := 2.022099747

test

xo xo₂ := 2.022099747

x x₂ := 2.023887212

test

xo xo₂ := 2.023887212

<-- exit relsucc (now at top level) = x = .9926675704, y =2.023887212}

x = 0.9926675704, y = 2.023887212

Prezentam comparativ rezultatele aplicarii metodei relaxarii succesive pentru diferite valori ale lui ω:

> relsucc( , vector(2, [0,0]), 0.99, 0.0001);

x = 0.9919288541, y = 2.024277742

> relsucc( , vector(2,[0,0]), 0.9, 0.0001);

x = 0.9091141587, y = 2.272710330

> relsucc( , vector(2,[0,0]), 1.01, 0.0001);

x = 1.007912178, y = 1.976281288

> relsucc( , vector(2,[0,0]), 1.1, 0.0001);

x = 1.071425885, y = 1.785716899

> relsucc( , vector(2,[0,0]), 0.7, 0.0001);

x = 0.6250722364, y = 3.124921273

> relsucc( , vector(2,[0,0]), 1.3, 0.0001);

x = 1.176464330, y = 1.470600344

> relsucc( , vector(2,[0,0]), 1, 0.0001); #

> gauss-seidel

x = 1.000014289, y = 1.999992856