Metoda lui Euler#

Formulele care nu pornesc de la un concept intuitiv sunt imposibil de reținut fără un efort considerabil.

Deoarece vreau să pun accent pe înțelegere, vom dezvolta întâi o intuiție geometrică, iar abia apoi vom trece prin rigoarea matematică, vom analiza eroarea introdusă, vom exemplificia un cod în Matlab și vom trage câteva concluzii.

Intuiția geometrică#

Să considerăm următoarea problemă - vrem să aflăm forma unei curbe necunoscute care începe dintr-un punct \(P_0\) și care satisface o ecuație diferențială ordinară.

În această situație, ecuația diferențială poate fi văzută drept panta tangentei la graficul funcției în oricare punct al acesteia. Deși forma curbei rămâne momentan un mister, este cunoscut un punct \(P_0(t_0,y_0)\) de început din care ne putem deplasa pe tangenta la grafic cu un pas mic, \(h\in\mathbb{R}\), către punctul următor, să zicem \(P_1(t_0+h,y_1)\).

Deoarece funcția \(y\) este continuă, făcând acest pas, tangenta în noul punct nu se va modifica foarte mult, deci \(P_1\) este foarte apropiat de curba inițială.

Dacă îl acceptăm pe \(P_1\) ca fiind pe curbă, atunci putem să reluăm algoritmul, pornind din \(P_1(t_1,y_1)\) către un alt punct, să zicem \(P_2(t_1+h,y_2)\), apoi din \(P_2(t_2,y_2)\) în alt punct, \(P_3(t_2+h,y_3)\), ș.a.m.d.

Exemplu#

Să aplicăm metoda lui Euler pentru funcția \(y=\sin t\), pornind din \(t_0=0.25\) cu deplasarea \(h=0.5\):

Trecerea de la \(P_0\) la \(P_1\)

Trecerea \(P_0\to P_1\to P_2\to P_3\)

Intuim cum, în funcție de pasul \(h\) ales, vom avea o eroare mai mică sau mai mare la final. Detalii suplimentare vor fi discutate ulterior, însă ceea ce trebuie reținut este conceptul fundamental al acestei metode.

Forma analitică#

Să formalizăm ceea ce am discutat.

Fie o ecuație diferențială ordinară, de forma \(\frac{dy}{dt}=f(t, y)\).

Considerăm acum progresia aritmetică \((t_n)_{n\in\mathbb{N}}\) aleasă astfel încât \(t_k=t_0+kh\), \(\forall k\in\overline{1,n}\), unde \(t_0,h\in\mathbb{R}\) sunt constante cunoscute și presupunem cunoscută și valoarea \(y(t_0)\).

Metoda lui Euler presupune utilizarea următoarei recurențe de aproximare:

\[ \boxed{y\left(t_{k+1}\right)\approx y(t_k)+h\cdot f\left(t_k, y(t_k)\right)}\,,\,\forall k\in\overline{0,n-1} \]

Alternativ, dacă notăm cu \(y_k\) aproximarea lui \(y(t_k)\) (\(y_k\approx y(t_k)\)), \(\forall k\in\overline{0,n}\), formula de recurență devine:

\[ \boxed{ y_{k+1}=y_k+h\cdot f(t_k,y_k) }\,,\,\,\forall k\in\overline{0,n-1} \]

Demonstrație

Vom scrie mai întâi \(y(t_{k+1})\) sub formă de serie Taylor în jurul punctului \(t_k\), \(\forall\,k\in\overline{0,n-1}\), și vom obține:

\[\begin{split} y(t_{k+1})&= y(t_k)+y'(t_k)(t_{k+1}-t_k) + \frac{y''(t_k)}{2!}(t_{k+1}-t_k)^2 + \frac{y^{(3)}(t_k)}{3!}(t_{k+1}-t_k)^3 + \dots\\ &= y(t_k)+y'(t_k)h + \frac{y''(t_k)}{2!}h^2 + \frac{y^{(3)}(t_k)}{3!}h^3 + \dots\\ \end{split}\]

Cunoaștem \(y'=f(t,y)\), deci facem înlocuirea:

\[\begin{split} y(t_{k+1})&= y(t_k)+h\cdot f(t,y(t_k)) + \frac{y''(t_k)}{2!}h^2 + \frac{y^{(3)}(t_k)}{3!}h^3 + \dots\\ &= y(t_k)+h\cdot f(t,y(t_k)) + O(h^2) \end{split}\]

Dacă acceptăm aproximarea \(y_k\approx y(t_k)\), se va obține ecuația pe care trebuia să o demonstrăm.

Vizualizarea algoritmului#

Să considerăm următoarea ecuație diferențială:

\[ \frac{dy}{dt}=y\cdot sin^2t,\,y(0)=0.5 \]

Aceasta va avea ca soluție funcția:

\[ y(t)=\frac{1}{2}\exp\left\{\frac{1}{2}\left(t-\sin t\cdot\cos t\right)\right\} \]

Cu toate aceasta, funcția \(y(t)\) va fi considerată necunoscută pentru noi și va fi de datoria noastră să o "aflăm", cel puțin numeric. Urmăm recurența \(y_{k+1}=y_k+h\cdot f(t_k,y_k)\), considerând \(h=0.5\):

\[\begin{split} y_{k+1}&=y_k+h\cdot\left(y_k\cdot sin^2t_k\right)\\ &=y_k+\frac{1}{2}\cdot y_k\cdot \sin^2\left(t_0+\frac{k}{2}\right) \end{split}\]

Nu vom plictisi cititorii cu rezultatele numerice, însă, din punct de vedere grafic, am ilustrat curba calculată prin metoda lui Euler, respectiv eroarea sa:

Curba calculată prin metoda lui Euler pentru \(h=0.5\)

Eroarea metodei lui Euler pentru \(h=0.5\)

O posibilă îmbunătățire a rezultatului ar fi să dublăm numărul de puncte, adică să înjumătățim distanța dintre ele, deci \(h\to0.25\). Recurența devine:

\[ y_{k+1}= y_k+\frac{1}{4}\cdot y_k\cdot \sin^2\left(t_0+\frac{k}{4}\right) \]

Vom ilustra acum cele două aproximări simultan:

Curba calculată prin metoda lui Euler, atât pentru \(h=0.5\), cât și pentru \(h=0.25\)

Eroarea metodei lui Euler când \(h=0.5\), respectiv pentru \(h=0.25\)

Analiza erorii#

Acuratețea acestei metode este foarte slabă - eroarea la fiecare pas poate fi determinată matematic prin formula:

\[ \boxed{\varepsilon_k(h)=\big|y(t_k)-y_k\big|\leq\frac{hM}{2K}\left[e^{K(t_k-t_0)}-1\right]}\,,\,\,\forall k\in\overline{1,n} \]

unde \(K\) este o constantă astfel încât funcția \(f\) să respecte condiția Lipschitz, iar \(M\) este o constantă aleasă astfel încât \( |y''(t)|\leq M \), \(\forall t\in[t_0,t_n]\).

Nu vom demonstra această relație, întrucât matematica devine destul de complicată, însă o demonstrație pertinentă poate fi găsită în bibliografie[7]. Detalii suplimentare în subcapitolul dedicat.

Cod ilustriv#

Algoritmul lui Euler poate fi scris în Matlab sub forma funcției euler. Aceasta va primi drept parametri:

a - capătul de la care se pornește, echivalentul lui \(t_0\);
b - capătul la care se dorește să se ajungă, echivalentul lui \(t_n\);
n - numărul de pași pe care îi va face algoritmul (tehnic, impune pasul \(h=\frac{b-a}{n}\));
y0 - evaluarea funcției pe care o căutăm în punctul \(t_0\), echivalentul lui \(y_0=y(t_0)\);
f - partea dreaptă a ecuației diferențiale \(\frac{dy}{dt}=f(t,y)\).

Funcția va returna y, soluția ecuației diferențiale.

function y = euler(a, b, n, y0, f)	
    h = (b - a) / n;
    y = zeros(n + 1, 1);
    y(1) = y0;

    for k = [1 : n]
        t = a + (k - 1) * h;
        y(k+1) = y(k) + h * f(t, y(k));
    end
end

Concluzii#

În practică, această metodă nu este folosită niciodată, deoarece este foarte imprecisă. În cazul exemplului vizualizat, eroarea devine mai mare decât \(0.25\) după \(3\), respectiv după \(8\) iterații.

Licență#

The book "Metode Numerice", written by Valentin-Ioan Vintilă, is licensed under CC BY-NC-SA 4.0