Metodo di Eulero

Disambiguazione – Se stai cercando l'omonimo metodo di fattorizzazione, vedi Metodo di fattorizzazione di Eulero.

In matematica e informatica, il metodo di Eulero è una procedura numerica del primo ordine per risolvere equazioni differenziali ordinarie (ODE) una volta fornito un valore iniziale. Si tratta del più basilare dei metodi espliciti per l'integrazione numerica di equazioni differenziali ordinarie, ed è il più semplice metodo Runge-Kutta. Prende il nome da Leonhard Euler, il quale lo espose nel suo libro Institutionum calculi integralis pubblicato nel 1768–70.

Introduzione

Illustrazione del metodo di Eulero. La curva sconosciuta è in blu, e la sua approssimazione polinomiale in rosso.

Si consideri il problema del calcolo della forma di una curva sconosciuta che inizia in un certo punto e soddisfa una data equazione differenziale. In questo caso, un'equazione differenziale ordinaria può essere immaginata come una formula tramite cui il coefficiente angolare della retta tangente alla curva può essere calcolato in ogni punto A 0 {\displaystyle A_{0}} della curva. Si prenda quindi un piccolo incremento sulla retta tangente, da A 0 {\displaystyle A_{0}} fino ad un punto A 1 {\displaystyle A_{1}} sufficientemente vicino. Si può supporre che in questo intervallo il coefficiente angolare non cambi significativamente, quindi se si assume che A 1 {\displaystyle A_{1}} è ancora sulla curva si può utilizzare nuovamente lo stesso ragionamento fatto per il punto A 0 {\displaystyle A_{0}} in modo che, dopo diversi passi, viene generata una curva poligonale A 0 A 1 A 2 A 3 {\displaystyle A_{0}A_{1}A_{2}A_{3}\dots } come mostrato nella figura a lato. Si tratta di una curva che non diverge troppo dalla curva originale sconosciuta; l'errore tra le due curve può essere diminuito se l'incremento è piccolo abbastanza e l'intervallo di calcolo è finito.

Formulazione

Si supponga di voler approssimare la soluzione del problema di Cauchy:

y ( t ) = f ( t , y ( t ) ) y ( t 0 ) = y 0 {\displaystyle y'(t)=f(t,y(t))\qquad y(t_{0})=y_{0}}

discretizzando la variabile t {\displaystyle t} , quindi definendo t n = t 0 + n h {\displaystyle t_{n}=t_{0}+nh} , con h {\displaystyle h} la dimensione di ogni intervallo. Tra t n {\displaystyle t_{n}} e t n + 1 = t n + h {\displaystyle t_{n+1}=t_{n}+h} il comportamento della soluzione può essere approssimato stimando:

y n + 1 = y n + h f ( t n , y n ) {\displaystyle y_{n+1}=y_{n}+hf(t_{n},y_{n})}

dove il valore di y n y ( t n ) {\displaystyle y_{n}\approx y(t_{n})} risulta essere un'approssimazione della soluzione della ODE al tempo t n {\displaystyle t_{n}} . Il metodo di Eulero è esplicito, ovvero la soluzione y n + 1 {\displaystyle y_{n+1}} è una funzione esplicita di y i {\displaystyle y_{i}} per i n {\displaystyle i\leq n} .

Il metodo integra una equazione differenziale ordinaria del primo ordine; tuttavia per un'equazione differenziale di ordine N:

y ( N ) ( t ) = f ( t , y ( t ) , y ( t ) , , y ( N 1 ) ( t ) ) {\displaystyle y^{(N)}(t)=f(t,y(t),y'(t),\ldots ,y^{(N-1)}(t))}

si introducono le variabili ausiliari z 1 ( t ) = y ( t ) , z 2 ( t ) = y ( t ) , , z N ( t ) = y ( N 1 ) ( t ) {\displaystyle z_{1}(t)=y(t),z_{2}(t)=y'(t),\ldots ,z_{N}(t)=y^{(N-1)}(t)} in modo da ottenere l'equazione equivalente:

z ( t ) = ( z 1 ( t ) z N 1 ( t ) z N ( t ) ) = ( y ( t ) y ( N 1 ) ( t ) y ( N ) ( t ) ) = ( z 2 ( t ) z N ( t ) f ( t , z 1 ( t ) , , z N ( t ) ) ) {\displaystyle \mathbf {z} '(t)={\begin{pmatrix}z_{1}'(t)\\\vdots \\z_{N-1}'(t)\\z_{N}'(t)\end{pmatrix}}={\begin{pmatrix}y'(t)\\\vdots \\y^{(N-1)}(t)\\y^{(N)}(t)\end{pmatrix}}={\begin{pmatrix}z_{2}(t)\\\vdots \\z_{N}(t)\\f(t,z_{1}(t),\ldots ,z_{N}(t))\end{pmatrix}}}

Si tratta di un sistema del primo ordine nella variabile z ( t ) {\displaystyle \mathbf {z} (t)} , e può essere approssimato col metodo di Eulero o con qualsiasi altro metodo per sistemi del primo ordine.

Esempio

Dato il valore iniziale del problema:

y = y y ( 0 ) = 1 {\displaystyle y'=y\qquad y(0)=1}

si vuole utilizzare il metodo di Eulero per approssimare y ( 4 ) {\displaystyle y(4)} .

Incremento h = 1

Illustrazione della funzione y = y , y ( 0 ) = 1 {\displaystyle y'=y,y(0)=1} : in blu il risultato del metodo di Eulero; in verde la regola del punto medio; in rosso la funzione esatta y = e t {\displaystyle y=e^{t}} . L'incremento è h = 1.0.

Il metodo di Eulero prevede:

y n + 1 = y n + h f ( t n , y n ) {\displaystyle y_{n+1}=y_{n}+hf(t_{n},y_{n})}

Per calcolare f ( t 0 , y 0 ) {\displaystyle f(t_{0},y_{0})} , essendo la funzione f {\displaystyle f} definita da f ( t , y ) = y {\displaystyle f(t,y)=y} si ha:

f ( t 0 , y 0 ) = f ( 0 , 1 ) = 1 h f ( y 0 ) = 1 1 = 1 {\displaystyle f(t_{0},y_{0})=f(0,1)=1\qquad h\cdot f(y_{0})=1\cdot 1=1}

Per ottenere il successivo valore da utilizzare nei calcoli:

y 0 + h f ( y 0 ) = y 1 = 1 + 1 1 = 2 {\displaystyle y_{0}+hf(y_{0})=y_{1}=1+1\cdot 1=2}

e allo stesso modo con y 2 {\displaystyle y_{2}} , y 3 {\displaystyle y_{3}} e y 4 {\displaystyle y_{4}} :

y 2 = y 1 + h f ( y 1 ) = 2 + 1 2 = 4 , y 3 = y 2 + h f ( y 2 ) = 4 + 1 4 = 8 , y 4 = y 3 + h f ( y 3 ) = 8 + 1 8 = 16. {\displaystyle {\begin{aligned}y_{2}&=y_{1}+hf(y_{1})=2+1\cdot 2=4,\\y_{3}&=y_{2}+hf(y_{2})=4+1\cdot 4=8,\\y_{4}&=y_{3}+hf(y_{3})=8+1\cdot 8=16.\end{aligned}}}

A causa della natura ripetitiva di questo algoritmo, può essere utile organizzare i calcoli sotto forma di grafico:

n {\displaystyle n} y n {\displaystyle y_{n}} t n {\displaystyle t_{n}} f ( t n , y n ) {\displaystyle f(t_{n},y_{n})} h {\displaystyle h} Δ y {\displaystyle \Delta y} y n + 1 {\displaystyle y_{n+1}}
0 1 0 1 1 1 2
1 2 1 2 1 2 4
2 4 2 4 1 4 8
3 8 3 8 1 8 16

La conclusione di questo calcolo è y 4 = 16 {\displaystyle y_{4}=16} . La soluzione esatta della equazione differenziale è y ( t ) = e t {\displaystyle y(t)=e^{t}} , quindi y ( 4 ) = e 4 54.598 {\displaystyle y(4)=e^{4}\approx 54.598} : l'approssimazione del metodo di Eulero con h = 1 {\displaystyle h=1} non è buona in questo caso, sebbene il suo comportamento è qualitativamente corretto. Una stima migliore si ottiene riducendo h {\displaystyle h} .

Diversi valori per l'incremento

Come cambia il risultato del metodo precedente con incremento h = 0,25.

Il metodo di Eulero è maggiormente accurato al diminuire del passo h {\displaystyle h} ; la tabella sottostante mostra i diversi risultati sfruttando passi di diverse dimensioni. La prima riga corrisponde all'esempio della sezione precedente mentre la seconda riga è illustrata in figura:

passo risultato con Eulero errore
1 16 38,598
0,25 35,53 19,07
0,1 45,26 9,34
0,05 49,56 5,04
0,025 51,98 2,62
0,0125 53,26 1,34

L'errore riportato nell'ultima colonna della tabella indica la differenza tra la soluzione esatta per t = 4 {\displaystyle t=4} e l'approssimazione con il metodo di Eulero. Si può constatare che al dimezzarsi del passo corrisponde approssimativamente un dimezzamento dell'errore, questo suggerisce che l'errore è all'incirca proporzionale alla dimensione del passo, almeno per piccoli valori di h {\displaystyle h} . Questo è vero in generale ed anche per altre equazioni.

Altri metodi, come la regola del punto medio, sono generalmente più precisi: sfruttando il punto medio l'errore è all'incirca proporzionale al quadrato del passo. Per questo motivo il metodo di Eulero viene definito del primo ordine, mentre il metodo del punto medio del secondo ordine.

Dalla tabella è possibile affermare che il passo necessario per approssimare un risultato fino alla terza cifra decimale è circa 0.00001: il che significa che saranno necessari 400000 passi. Questo alto numero di passi implica un alto costo computazionale. Per questo motivo generalmente si preferiscono metodi numerici più precisi (cioè di ordine maggiore) come ad esempio il metodo di Runge-Kutta oppure i metodi lineari multistep, soprattutto se è necessaria grande accuratezza.

Derivazione ed errore di troncamento locale

Ci sono diversi modi per ricavare il metodo di Eulero; ad esempio si può considerare l'espansione di Taylor di y {\displaystyle y} intorno a t 0 {\displaystyle t_{0}} :

y ( t 0 + h ) = y ( t 0 ) + h y ( t 0 ) + 1 2 h 2 y ( t 0 ) + O ( h 3 ) {\displaystyle y(t_{0}+h)=y(t_{0})+hy'(t_{0})+{\frac {1}{2}}h^{2}y''(t_{0})+O(h^{3})}

e quindi sostituire y = f ( t , y ) {\displaystyle y'=f(t,y)} ignorando i termini al quadrato e di ordine superiore. L'errore di troncamento locale L T E {\displaystyle \mathrm {LTE} } che si ha con tale approssimazione è dato dalla differenza tra la soluzione y 1 {\displaystyle y_{1}} :

y 1 = y 0 + h f ( t 0 , y 0 ) {\displaystyle y_{1}=y_{0}+hf(t_{0},y_{0})}

e la soluzione al tempo t 1 = t 0 + h {\displaystyle t_{1}=t_{0}+h} fornita con il precedente sviluppo, ovvero:

L T E = y ( t 0 + h ) y 1 = 1 2 h 2 y ( t 0 ) + O ( h 3 ) {\displaystyle \mathrm {LTE} =y(t_{0}+h)-y_{1}={\frac {1}{2}}h^{2}y''(t_{0})+O(h^{3})}

dove la derivata terza di y {\displaystyle y} è supposta limitata. L'errore è quindi proporzionale, per piccoli h {\displaystyle h} , ad h 2 {\displaystyle h^{2}} : l'estensione del metodo di Eulero, che si ha ad esempio con i metodi di Runge-Kutta, consente di avere stime più accurate, in cui l'errore è proporzionale a potenze maggiori di h {\displaystyle h} ,

Un'altra derivazione del metodo utilizza le differenze finite per scrivere la derivata:

y ( t 0 ) y ( t 0 + h ) y ( t 0 ) h {\displaystyle y'(t_{0})\approx {\frac {y(t_{0}+h)-y(t_{0})}{h}}}

in modo che sostituendo tale espressione in y = f ( t , y ) {\displaystyle y'=f(t,y)} si ottiene il metodo di Eulero.

Si può infine integrare direttamente l'equazione differenziale tra t 0 {\displaystyle t_{0}} e t 0 + h {\displaystyle t_{0}+h} e applicare il teorema fondamentale del calcolo integrale:

y ( t 0 + h ) y ( t 0 ) = t 0 t 0 + h f ( t , y ( t ) ) d t {\displaystyle y(t_{0}+h)-y(t_{0})=\int _{t_{0}}^{t_{0}+h}f(t,y(t))\,\mathrm {d} t}

Quindi si approssima l'integrale con la regola del rettangolo:

t 0 t 0 + h f ( t , y ( t ) ) d t h f ( t 0 , y ( t 0 ) ) {\displaystyle \int _{t_{0}}^{t_{0}+h}f(t,y(t))\,\mathrm {d} t\approx hf(t_{0},y(t_{0}))}

e combinando entrambe le equazioni si trova il metodo di Eulero.

Implementazione

Di seguito una possibile implementazione del metodo in linguaggio Matlab

function [th,uh]=eulero_avanti(f,t_0,t_max,y0,h)

th=t_0:h:t_max;

uh=zeros(length(th),1);
uh(1)=y0;

for i=1:length(th)-1
    u_new = uh(i) + h*f(th(i),uh(i));
    uh(i+1)=u_new;
end
end

Bibliografia

  • (EN) Kendall A. Atkinson, An Introduction to Numerical Analysis, 2nd, New York, John Wiley & Sons, 1989, ISBN 978-0-471-50023-0.
  • (EN) Uri M. Ascher e Linda R. Petzold, Computer Methods for Ordinary Differential Equations and Differential-Algebraic Equations, Philadelphia, Society for Industrial and Applied Mathematics, 1998, ISBN 978-0-89871-412-8.
  • (EN) John C. Butcher, Numerical Methods for Ordinary Differential Equations, New York, John Wiley & Sons, 2003, ISBN 978-0-471-96758-3.
  • (EN) Ernst Hairer, Syvert Paul Nørsett e Gerhard Wanner, Solving ordinary differential equations I: Nonstiff problems, Berlin, New York, Springer-Verlag, 1993, ISBN 978-3-540-56670-0.
  • (EN) Arieh Iserles, A First Course in the Numerical Analysis of Differential Equations, Cambridge University Press, 1996, ISBN 978-0-521-55655-2.
  • (EN) Josef Stoer e Roland Bulirsch, Introduction to Numerical Analysis, 3rd, Berlin, New York, Springer-Verlag, 2002, ISBN 978-0-387-95452-3.
  • (EN) Taras I. Lakoba, Simple Euler method and its modifications (PDF) (Lecture notes for MATH334, University of Vermont), 2012. URL consultato il 29 febbraio 2012.

Voci correlate

Altri progetti

Altri progetti

  • Wikimedia Commons
  • Collabora a Wikimedia Commons Wikimedia Commons contiene immagini o altri file su metodo di Eulero

Collegamenti esterni

  • (EN) Euler's Method for O.D.E.'sArchiviato il 16 luglio 2013 in Internet Archive., by John H. Matthews, California State University at Fullerton.
  Portale Matematica: accedi alle voci di Wikipedia che trattano di matematica