METODE NUMERICE DE SOLU|IONARE A ECUA|IILOR DIFEREN

Transcript METODE NUMERICE DE SOLU|IONARE A ECUA|IILOR DIFEREN

CURS 6 METODE NUMERICE DE REZOLVARE A ECUATIILOR DIFERENTIALE

Introducere O ecuatie diferentiala poate sa nu posede solutie sau, chiar daca are solutie, nu totdeauna aceasta se poate explicita. In multe situatii, mai ales in cazul ecuatiilor diferentiale neliniare, trebuie sa ne consideram multumiti daca obtinem o aproximatie a solutiei ecuatiei diferentiale. In cele ce urmeaza, utilizand metode numerice, se vor obtine seturi de puncte care, atunci cand se utilizeaza discretizari suficient de fine se poate aproxima solutia ecuatiei diferentiale considerate, asa cum se prezinta in fig. 1.

ECUATII DIFERENTIALE DE ORDINUL 1 Forma generala a acestei ecuatii este: dy dx = f(x,y) dat fiind faptul ca ecuatiile diferentiale de ordin superior pot fi reduse la sisteme de ecuatii diferentiale de ordinul 1.

EXEMPLE 1) Ecuatia de miscare a unei particule de masa m într-un câmp de forte F(x): m d 2 dt x 2 = F(x)



dp dt = F(x) 2) Ecuatia diferentiala a fibrei medii deformate: d 2 w dx 2



= dw dx Intr-adevar, considerand ecuatia diferentiala: d n y d t n

 

a 0 a n y



a 1 a n y

   

a n



1 a n y





 

 

si introducand variabilele: y



x 1 , y

 

x 2 , y

 

x 3 ,



, y





 

x n

rezulta ca: y

 

x 1

 

x 2 , y

 

x 2

 

x 3 ,



, y





 

x n



 

x n In acest fel ecuatia diferentiala de ordinul n scrisa sub forma (1) se poate echivala cu urmatorul sistem de n ecuatii diferentiale de ordinul 1:

         

x

n   

x x x x

1 2 3       n  1   

x x x x

2 3 4 n

a

x

1 

a

x

2   

a

n ,

a

n  1

x a

n n 

f

 

.

a) METODA EULER (METODA LINIILOR POLIGONALE)

Aceasta este una dintre cele mai simple tehnici de solu\ionare numerica a ecuatiilor diferentiale, fiind cunoscuta si sub denumirea de metoda tangentelor.

Se considera ecuatia diferentiala (inclusiv conditia initiala aferenta): y

 



x , y



, y

  

y 0 Daca se considera pe axa Ox o diviziune echidistanta de pas h se poate gasi un punct (x 1 , y 1 ) = (x 0 +h, y 1 ) pe tangenta la curba ce reprezinta solutia ecuatiei diferentiale in punctul (x 0 , y 0 ), asa cum se prezinta in figura alaturata. Se poate scrie: unde: y 0

 



x 0 , y 0

 

x y 1 0

 





0 x 0



y 0

 

y 1



y 0



h y 0



Daca se noteaza x 0 +h = x 1 , atunci punctul de coordonate (x 1 , y 1 ) situat pe tangenta considerata reprezinta a aproximatie a punctului de coordonate (x 1 , y(x 1 )) situat pe curba ce reprezinta solutia ecuatiei diferentiale. Este evident faptul ca eroarea metodei este cu atat mai redusa cu cat valoarea pasului diviziunii considerate - h este mai mic.

Se poate construi astfel un proces iterativ: y n





y n



h y n

 

y n



h f



x n , y n



, x n



x 0



n h a carui reprezentare grafica se prezinta in figura alaturata.

Aproximarea solutiei ecuatiei diferentale

APLICATIA 1

Considerand ecuatia diferentiala: y

 

2 x y si conditia initiala y(1)=1, sa se obtina o aproximatie pentru a gasi valoarea lui y in punctul de abscisa x = 1.5, utilizand un pas h cu valorile 0.1, 0.05 si 0.01, apoi sa se calculeze si erorile relative raportate la valoarea exacta, cu patru zecimale.

Solutie Solutia analitica a ecuatiei diferentiale considerate este: y

  

e x 2



x

1.00

1.10

1.20

1.30

1.40

1.50

Tab. 1 Valorile obtinute pentru h = 0.10

y

1.0000

1.2000

1.4640

1.8154

2.2874

2.9278

Exact

1.0000

1.2337

1.5527

1.9937

2.6117

3.4904

Eroare

0.0000

0.0337

0.0887

0.1784

0.3244

0.5625

Eroare %

0.00

2.73

5.71

8.95

12.42

16.12

x

1.00

1.10

1.20

1.30

1.40

1.50

Tab. 2 Valorile obtinute pentru h = 0.05

y

1.0000

1.2155

1.5044

1.8955

2.4311

3.1733

Exact

1.0000

1.2337

1.5527

1.9937

2.6117

3.4904

Eroare

0.0000

0.0182

0.0483

0.0982

0.1806

0.3171

Eroare %

0.00

1.47

3.11

4.93

7.98

9.08

x

1.00

1.10

1.20

1.30

1.40

1.50

Tab. 3 Valorile obtinute pentru h = 0.01

y

1.0000

1.2298

1.5423

1.9723

2.5719

3.4197

Exact

1.0000

1.2337

1.5527

1.9937

2.6117

3.4904

Eroare

0.0000

0.0039

0.0104

0.0214

0.0398

0.0707

Eroare %

0.00

0.31

0.67

1.07

1.52

2.03

b) METODA EULER-HEUN

Aceasta metoda este cunoscuta si sub denumirea de metoda Euler imbunatatita, situatia fiind prezentata in figura de mai jos pentru primul pas de integrare. Procesul iterativ pentru solutionarea numerica este descris de relatia urmatoare: y n



1 y * n



 



x n , y n y n



h y n



h f



x n , y n .

  

2 f



x n



1 , y * n





Analizand figura se poate constata ca valoarea lui y 1 este mai apropiata de solutie decat valoarea obtinuta pentru y 1 * , ceea ce poate incadra aceasta metoda in clasa metodelor predictor-corector.

Intr-adevar, aceasta situatie se poate interpreta in felul urmator: Pasul predictor:

y 0



f , y



Pasul corector: y

 * 1 

y 0

 

h f

h f



0 ,

, y

 





x 1 * 1 In continuare, cu referire la ecuatia diferentiala considerata anterior, se prezinta rezultatele obtinute prin aplicarea acestei metode pentru diferiti pasi de integrare si programul folosit pentru solutionarea numerica.

Tab. 4 Valorile obtinute pentru h = 0.1

x

n 1.00

1.10

1.20

1.30

1.40

1.50

y

n 1.0000

1.2320

1.5479

1.9832

2.5908

3.4509

Exact

1.0000

1.2337

1.5527

1.9937

2.6117

3.4904

Eroare

0.0000

0.0017

0.0048

0.0106

0.0209

0.0394

Eroare %

0.00

0.14

0.31

0.53

0.80

1.13

x

1.00

1.10

1.20

1.30

1.40

1.50

Tab.5 Valorile obtinute pentru h=0.05

y

1.0000

1.2332

1.5514

1.9909

2.6060

3.4795

Exact

1.0000

1.2337

1.5527

1.9937

2.6117

3.4904

Eroare

0.0000

0.0004

0.0013

0.0029

0.0057

0.0108

Eroare %

0.00

0.04

0.08

0.14

0.22

0.31

x

1.00

1.10

1.20

1.30

1.40

1.50

Tab.6 Valorile obtinute pentru h = 0.01

y

1.0000

1.2337

1.5527

1.9936

2.6115

3.4899

Exact

1.0000

1.2337

1.5527

1.9937

2.6117

3.4904

Eroare

0.0000

0.0001

0.0002

0.0005

Eroare %

0.00

0.005

0.010

0.025

c) METODA DEZVOLTARII IN SERIE TAYLOR

Aceasta metoda nu prezinta o larga aplicabilitate deoarece, in principal, rezultatele sunt comparabile cu cele obtinute prin aplicarea metodei Euler-Heun in condi\iile unor calcule mai laborioase.

Asa cum se stie, dezvoltarea unei functii y(x) in jurul unui punct x = a are expresia: y





  

x 1



a !







 





2 2 !

 

Pentru a apropia expresia de mai sus de problema curenta, daca se considera cazul in care a = x n si x = x n +h, se obtine: y



x n



h n n h





 

h 2 2

 

Daca se presupune ca functia y(x) este solutia ecuatiei diferentiale: y

 

Oprind numai doi termeni din se poate scrie: y



x n



    





x n









x n , y

  

h Se observa faptul ca daca in ecuatia de mai sus se inlocuieste y(x n +h) cu y n+1 si y(x n ) cu y n se obtine formula ce caracterizeaza metoda Euler: y n





y n



h f



x n , y n



Oprind trei termeni si facand inlocuirile antementionate se obtine urmatoarea formula: y n





y n



y n





y n



h 2 Metoda poate fi considerata drept cheia de bolta a majoritatii metodelor de rezolvare numerica a ecuatiilor diferentiale.

d) METODE TIP RUNGE-KUTTA

Aceasta clasa de metode reprezinta una dintre cele mai folosite in abordarea numerica a ecuatiilor diferentiale, imbinand numarul relativ redus de operatii elementare cu acuratetea rezultatelor.

Metoda Runge-Kutta de ordinul II consta in ga[sirea constantelor a, b astfel incat expresia: y n





y n



a k 1



b k 2 cu: k 1



h k 2



h f f

 

x x n , y n n

 



, , y n

 

k 1 .



sa se apropie de dezvoltarea in serie Taylor pentru cat mai multi termeni posibili. Meritul principal al acestei clase de metode rezida deci in aceea ca se apropie de acuratetea unei dezvoltari in serie Taylor fara insa a fi nevoie sa se calculeze si derivatele de ordin superior.

Se fac urmatoarele observatii: a) Daca: a





1 , b

 

1 2 , b

 

1 2 atunci relatia anterioara reprezinta dezvoltarea in serie Taylor, oprind primii trei termeni, a functiei considerate.

b) Daca: a



1 2 , b



1 2 ,

 

1 , b



1 atunci relatia anterioara se reduce la expresia metodei Euler-Heun.

c) Se poate constata ca metoda Euler este de fapt o procedura Runge-Kutta de ordinul I.

Folosind aceleasi consideratii se determina astfel algoritmul de calcul si in cazul metodei Runge-Kutta de ordinul IV (se opresc primii patru termeni din dezvoltarea in serie Taylor):

y n





k 1



h f y n





n 1 6 , y n





, 1



2 k 2



2 k 3



k 4



, k 2



h f

 

x n



1 2 h , y n



1 2 k 1

 

, k 3



k 4



h h



1 f f

  

x x n n 2



h , y h , n



y k n 3

 

2 1 k 2

 

, In continuare se prezinta rezultatele numerice obtinute in cazul aplicarii acestei metode pentru rezolvarea ecuatiei diferentiale considerate anterior si programul folosit in acest scop.

x

1.00

1.10

1.20

1.30

1.40

1.50

Tab. 7 Valorile obtinute pentru h = 0.1

y

1.0000

1.2337

1.5527

1.9937

2.6116

3.4902

Exact

1.0000

1.2337

1.5527

1.9937

2.6117

3.4904

Eroare

0.0000

0.0001

0.0002

Eroare %

0.00

0.0038

0.0076

Tab. 8 Valorile obtinute pentru h = 0.05

x

1.00

1.10

1.20

1.30

1.40

1.50

y

1.0000

1.2337

1.5527

1.9937

2.6117

3.4903

Exact

1.0000

1.2337

1.5527

1.9937

2.6117

3.4904

Eroare

0.0000

0.0001

Eroare %

0.00

0.0028

Tab. 9 Valorile obtinute pentru h = 0.01

x

1.00

1.10

1.20

1.30

1.40

1.50

y

1.0000

1.2337

1.5527

1.9937

2.6117

3.4903429

Exact

1.0000

1.2337

1.5527

1.9937

2.6117

3.4904

Eroare

0.0000

0.0000571

Eroare %

0.00

0.0016

e) METODA MILNE

Aceasta este o metoda[ care lucreaza in doi pasi

•

Pasul predictor: y * n





y n





4 h 3 2 y n

 

y n



 

2 y n





y n

 

f y n



 

y n



 

f f



x n , y n



1 n



, y ,



x n



2 , y n



1 n



 

, .

, n



•

Pasul corector: y n





y n





y h 3 n















1 , 4 y n



y * n



 

y n





Se remarca faptul ca valorile pentru y 0 , y 1 , y 2 si y 3 trebuie cunoscute pentru demararea procesului iterativ. Acestea se pot calcula cu oricare dintre metodele anterioare, recomandandu-se insa ca acest calcul initial si se faca cu un nivel inalt de acuratete.

CONCLUZII

Aplicarea uneia sau alteia dintre aceste metode in vederea solutionarii unei probleme concrete implica o analiza atenta a erorilor induse, in vederea realizarii unui compromis acceptabil intre precizie si efortul de calcul cerut.