Curs 5 MN Sisteme Neliniare

8
Metode Numerice Curs 5 Sisteme neliniare Prof. dr. mat. Dumitru NICOARĂ 1 SISTEME DE ECUATII ALGEBRICE (continuare) Fie sistemul de n- ecuatii algebrice neliniare: 0 ) , ..... , , ( . . 0 ) , ..... , , ( 0 ) , ..... , , ( 2 1 2 1 2 2 1 1 n n n n x x x f x x x f x x x f (4. 25) Fie sistemul de doua ecuaţii: 0 1 75 . 0 25 . 0 ) , ( 0 2 ) , ( 2 2 2 1 2 1 2 2 2 2 1 2 1 1 x x x x f x x x x f (4.26) Vom detrmina soluţiile prin mai multe metode. Metoda 1 Metoda grafică Ecuaţiile din (4.26) reprezintă nişte suprafeţe. Soluţiile se vor găsi la intrsecţia linilor conturului celor doua suprateţe, Fig.4.26. Figura 4.26

Transcript of Curs 5 MN Sisteme Neliniare

Page 1: Curs 5 MN Sisteme Neliniare

Metode Numerice Curs 5 Sisteme neliniare Prof. dr. mat. Dumitru NICOARĂ

1

SISTEME DE ECUATII ALGEBRICE (continuare)

Fie sistemul de n- ecuatii algebrice neliniare:

0),.....,,(

.

.

0),.....,,(

0),.....,,(

21

212

211

nn

n

n

xxxf

xxxf

xxxf

(4. 25)

Fie sistemul de doua ecuaţii:

0175.025.0),(

02),(

2

2

2

1212

2

2

2

1211

xxxxf

xxxxf (4.26)

Vom detrmina soluţiile prin mai multe metode.

Metoda 1 Metoda grafică

Ecuaţiile din (4.26) reprezintă nişte suprafeţe. Soluţiile se vor găsi la intrsecţia linilor

conturului celor doua suprateţe, Fig.4.26.

Figura 4.26

Page 2: Curs 5 MN Sisteme Neliniare

Metode Numerice Curs 5 Sisteme neliniare Prof. dr. mat. Dumitru NICOARĂ

2

Pentru a desena conturul trebuie să gasim valorile punctelor de intersecţie x1 şi x2 .

Aceste valori vor fi stocate ȋn două matrice, una conţinând valorile x1 iar cealaltă valorile

x2. Vom folosi comanda contour (>>help contour).

Rezolvare o vom face ȋn MATLAB. Vom scrie un program ȋn fila scrtipt:

sys_nl1_grafic.m

% Solutia grafica

%------------------

clear

clc

format compact

% informatii pentru axe si generarea matricelor cu valorile

pe contur

% axele

x1 = -2:0.1:2;

x2 = -2:0.1:2;

%matricele de valori

[X1 X2] = meshgrid(x1,x2);

f1 = X1.^2 + X2.^2 -2;

f2 = 0.25*X1.^2 + 0.75*X2.^2 - 1;

%------------------------------

contour(x1,x2,f1,[2 2],'r')

axis square

% traseaza conturul pentru f1 la valoarea zero [0 0]

% pentru un singur contur in MATLAB trebuie sa specificam

un vector

% cu 2 elemnte

hold on

contour(x1,x2,f2,[2 2],'k')

xlabel('\bfx_1 ','FontSize',12,'FontName','Arial', ...

'Color','b','VerticalAlignment','top')

ylabel('\bfx_2 ','FontSize',12,'FontName','Arial', ...

'Color','b')

title('\bf\itSolutia grafica a sistemului x1^2+x2^2-2=0;

0.25x1^2+0.75x2^2-1=0','Color','k')

grid

legend('f_1 = 0','f_2 = 0',1)

Page 3: Curs 5 MN Sisteme Neliniare

Metode Numerice Curs 5 Sisteme neliniare Prof. dr. mat. Dumitru NICOARĂ

3

Soluţiile sunt:

x1 x2

------------------

1 1

-1 1

1 -1

-1 -1

Metoda 2 Soluţia obţinută folosind comanda fsolve >> help fsolve

Descriere : fsolve determină rădăcinile (zerourile) unui sistem de ecuaţii neliniare

Sintaxa:

a) x = fsolve(fun,x0) ; start din x0 şi ȋncearcă să rezolve ecuaţia descrisă ȋn fun.

b) [x, fval] = fsolve(fun,x0) ; start din x0 şi ȋncearcă să rezolve ecuaţia descrisă ȋn

fun.m; returnează valuarea funcţiei obiectiv fun

pentru soluţia x.

c) x = fsolve(fun,x0,options) ; [x,fval] = fsolve(fun,x0,options)

options=optimset('Display','iter');

Paşii de parcurs:

10. Scriem o M-filă tip funcţie, pentru introducerea sistemului

x1

x2

Solutia grafica a sistemului x12+x22-2=0; 0.25x12+0.75x22-1=0

-2 -1.5 -1 -0.5 0 0.5 1 1.5 2-2

-1.5

-1

-0.5

0

0.5

1

1.5

2

f1 = 0

f2 = 0

Page 4: Curs 5 MN Sisteme Neliniare

Metode Numerice Curs 5 Sisteme neliniare Prof. dr. mat. Dumitru NICOARĂ

4

function F = myfun(x)

ecuaţiile....

Salvăm această M-fila ca myfun.m. (myfun=numele dat de utilizator )

20. Setăm punctul iniţial de start

>> x0 = [x10, x20]

şi optiunile (options) pentru afişarea raspunsului (output):

>>options=optimset('Display','iter');

30. Aplelăm fsolve:

>>[x,fval] = fsolve(@myfun,x0,options)

Observaţie.

Recomandăm ca paşii 2 -3 să fie scrişi tot ȋntr-o M-fila (script) pe care să o apelăm din

fereastra de comenzi Matlab.

Aplicaţie. Cosiderăm sistemul rezolvat, mai sus, prin metoda grafică, (4.26):

0175.025.0),(

02),(

2

2

2

1212

2

2

2

1211

xxxxf

xxxxf

10. Scriem o M-filă tip funcţie, pentru introducerea sistemului de ecuaţii, salvată ca

sys_nl1.m:

function f=sys_nl1(x)

f=[(x(1)^2+x(2)^2-2),(0.25*x(1)^2+0.75*x(2)^2-1)];

Pentru paşii 2 şi trei vom scrie un program Matlab pe care-l salvăm ȋntr-o filă script

sol_fsolve_sys_nl1.m .

% Solutia numerica folosind MATLAB function fsolve

% f1 = x1*x1 + x2*x2 - 2;

% f2 = 0.25*x1*x1 +0.75*x2*x2 -1;

%------------------ %

clear

clc

format compact

% [0.5 0.5] punctul initial

% problem is set up in sys_nl1.m

x = fsolve(@sys_nl1,[0.5 0.5],optimset('Display','off'));

% Tipareste

fprintf('Solutia numerica folosind MATLAB function fsolve

\n')

fprintf('Sistemul de ecuatii este in sys_nl1.m\n\n')

fprintf(' x1 x2 \n');

fprintf('------------------\n')

disp([double(x(1)) double(x(2))]);

Page 5: Curs 5 MN Sisteme Neliniare

Metode Numerice Curs 5 Sisteme neliniare Prof. dr. mat. Dumitru NICOARĂ

5

Metoda 3 Soluţia obţinută folosind metoda Newton -Raphson

La metoda Newton pentru o ecuaţie de o singură variabilă, cursul 3, plecând de formula

Taylor, am obţinut formula de iteraţie

)(

)(1

i

iii

xf

xfxx

sau

)(

)(

i

i

xf

xfx

(4.27)

Pentru un sistem neliniar de forma (4.25) se poate proceda la fel.

Pentru simplitate, vom exemplifica metoda pentru un sistem de doua ecuaţii neliniare cu

două variabile

0),(

0),(

212

211

xxf

xxf (4.28)

Pentru fiecare ecuaţie scriem formula Taylor de ordinul unu şi obţinem sistemul

2

,2

,21,2

1

,2

,11,1,21,2

2

,1

,21,2

1

,1

,11,1,11,1

)()(

)()(

x

fxx

x

fxxff

x

fxx

x

fxxff

i

ii

i

iiii

i

ii

i

iiii

(4.29)

La acest pas, i+1, aproximarea rădăcinilor se obţine pentru 0,0 1,21,1 ii ff . Valorile

de la pasul anterior i fiind cunoscute, din (4.29) obţinem sistemul de două ecuaţii pentru

determinarea a două necunoscte 1,1 ix şi 1,2 ix :

2

,2

,2

1

,2

,1,21,2

2

,2

1,1

1

,2

2

,1

,2

1

,1

,1,11,2

2

,1

1,1

1

,1

x

fx

x

fxfx

x

fx

x

f

x

fx

x

fxfx

x

fx

x

f

i

i

i

iii

i

i

i

i

i

i

iii

i

i

i

(4.30)

Soluţia sistemului, obţinută de exemlu folosind regula lui Cramer, este:

Page 6: Curs 5 MN Sisteme Neliniare

Metode Numerice Curs 5 Sisteme neliniare Prof. dr. mat. Dumitru NICOARĂ

6

1

,2

2

,1

2

,2

1

,1

1

,2

,1

1

,1

2

,21,2

1

,2

2

,1

2

,2

1

,1

2

,1

,2

2

,2

,1

,11,1

x

f

x

f

x

f

x

f

x

ff

x

ff

xx

x

f

x

f

x

f

x

f

x

ff

x

ff

xx

iiii

i

i

i

i

ii

iiii

i

i

i

i

ii

(4.31)

Ecuaţiile (4.31) sunt ecuaţiile iterative ale metodei Newton-Raphson pentru un sistem

neliniar cu două ecuaţii şi două necunoscute .

Observăm că numitorul este valoarea determinantului matricei Jacobian a sistemului

)(detdetdet

2

2

1

2

2

1

1

1

Xgrad

x

f

x

f

x

f

x

f

J

(4.32)

Prin analogie cu relaţia (4.27), putem rescrie relaţiile (4.31) ȋn forma vectorială

)()()()( 1 xFXgradinvxFXgradX (4.33)

unde

)(

)()(

2

1

Xf

XfXF ;

2

2

1

2

2

1

1

1

)(

x

f

x

f

x

f

x

f

Xgrad (4.34)

Page 7: Curs 5 MN Sisteme Neliniare

Metode Numerice Curs 5 Sisteme neliniare Prof. dr. mat. Dumitru NICOARĂ

7

Algoritmul metodei Newton-Raphson pentru sisteme neliniare:

Step 1: Alegem: X (initial), N, nMax, eps,

i = 1

Step 2: Calculam F(X), Grad(X)

del X = -inv(Grad(X)) * F(X)

X = X + delX

if norm(f(X)) < eps

exit - convergent

if i == nMax

restart

i = i + 1

Implemetarea algoritmului s-a făcut în M- fila sol_newton_sys_nl1.m

% MATLAB Code metoda Newton - Raphson pentru sisteme neliniare

% sol_newton_sys_nl1.m

%-------------------------------------%

% fie sistemul

% f1 = x1*x1 + x2*x2 - 2;

% f2 = 0.25*x1*x1 +0.75*x2*x2 -1;

%-------------------------------------%

clear

clc

format compact

nmax = 10; eps1 = 1.0e-08;

X= [0.3 0.2];

fprintf('Solutia numerica cu metoda Newton\n')

fprintf('-------------------------------------\n')

fprintf(' x(1) x(2) f(1) f(2) dx(1) dx(2) \n');

fprintf('-------------------------------------------------\n')

for i = 1:nmax

F = sys_nl1(X); % calculeaza functia

if norm(F) < eps1

break

end

GradF = Grad_sys_nl1(X);%calculeaza matricea Jacobian (grad)

valX = -inv(GradF)*F';

if i == nmax

fprintf('Numarul maximum de iteratii depasit\n')

break

end

disp([X F valX']);

X = X + valX';

end

Page 8: Curs 5 MN Sisteme Neliniare

Metode Numerice Curs 5 Sisteme neliniare Prof. dr. mat. Dumitru NICOARĂ

8

Ecuaţiile sistemului se obţin apelând M-fila sys_nl1.m

function f=sys_nl1(x)

f=[(x(1)^2+x(2)^2-2),(0.25*x(1)^2+0.75*x(2)^2-1)];

iar matricea jacobian (gradientul) se calcululează apelând M-fila Grad_sys_nl1.m

function df = Grad_sys_nl1(x)

% x is a vector of x's [x(1) x(2) .. .. x(n)]

% f1 = x1*x1 + x2*x2 - 2;

% f2 = 0.25*x1*x1 +0.75*x2*x2 -1;

df = [2*x(1) 2.0*x(2) ; ...

0.5*x(1) 1.5*x(2)];

Obţinem soluţia prin executarea M-filei sol_newton_sys_nl1.m, prin comanda:

>> sol_newton_sys_nl1

Rezultă

Solutia numerica cu metoda Newton

-------------------------------------

x(1) x(2) f(1) f(2) dx(1) dx(2)

--------------------------------------------------------------

0.3000 0.2000 -1.8700 -0.9475 1.5167 2.4000

1.8167 2.6000 8.0603 4.8951 -0.6331 -1.1077

1.1836 1.4923 1.6278 1.0204 -0.1693 -0.4111

1.0142 1.0812 0.1977 0.1339 -0.0141 -0.0782

1.0001 1.0030 0.0063 0.0046 -0.0001 -0.0030

1.0000 1.0000 0.0000 0.0000 -0.0000 -0.0000