METODO DI EULERO Il metodo di Eulero è il più
semplice metodo one-step esplicito per
risolvere una ODE. Questo metodo assume che, entro un
piccolo incremento della variabile indipendente dal
punto iniziale x0 , la derivata y' si
mantenga costante al valore f(x0,y0)
:
e,nota
la
y', : |
|
Il processo può
essere ripetuto usando (x1, y1)
come punto di partenza in modo da calcolare la y'(y1,
x1) e quindi una nuova y2 (è uso
scrivere h al posto di Δx): E così via. La formula generale del metodo di Eulero è: |
|
Si può osservare come il metodo di Eulero
in fondo corrisponde all'espansione di Taylor troncata
al primo ordine : l'errore di troncamento è quindi dell'ordine di O(h²). Per migliorare l'approssimazione
si può diminuire il valore di h e in questo caso
diminuisce anche l'errore di troncamento in modo
proporzionale ad h². Tuttavia, poiché, per il
dato insieme di x entro il quale l'integrazione è
eseguita, il numero di passi è proporzionale a 1/h ()
allora l'errore totale diventa del tipo O(h) ovvero
proporzionale ad h.
Si verifica spesso che il metodo di Eulero non fornisce un’approssimazione accurata della soluzione e l’errore cresce allontanandosi dal punto x1 in cui è assegnata la condizione iniziale. Bisogna scegliere un passo h molto piccolo per avere un errore globale accettabile. |
|
METODO DI HEUN Il difetto principale del metodo
di Eulero è che la pendenza all'inizio dell'intervallo
f(xi,yi) non è
rappresentativa del suo valore sull'intero intervallo.
Una soluzione è sostituirla con la media
aritmetica tra il valore all'inizio e quello alla fine
:
Tuttavia non conosciamo la yi+1 ( infatti le derivate f(x,y) in generale sono una funzione delle x e delle y). Una prima stima di yi+1 può essere fornita dal metodo di Eulero: |
METODO DEI POLIGONI O DEL PUNTO
MEDIO Nel metodo dei poligoni o del
punto medio si migliora la tecnica per calcolare
il predictor
.
Infatti si utilizza il valore della pendenza a metà intervallo utilizzando il metodo di Eulero: e con questo predictor si calcola la pendenza corretta. |
|
La formula finale del metodo
dei poligoni è: Anche in questo caso l'errore locale è di tipo O(h³) per cui l'errore totale dovuto all'applicazione su più intervallini è O(h²). |
Dal confronto si ricava che deve
essere: Si ottiene comunque un sistema di
3 equazioni a 4 incognite che non è determinato. Si
può scegliere il coefficiente a2, che
dipende dal punto dove viene calcolata la pendenza k2.
Alcune possibili scelte:
|
|
|
|
Runge-Kutta genera il metodo del punto medio o dei poligoni |
Runge-Kutta genera il metodo di Ralston |
Questi comunque sono tutti sviluppi di
Taylor fino al secondo ordine. RK al terzo ordine. Un metodo di RK al terzo ordine genera un errore locale è di tipo O(h4) e un l'errore totale dovuto all'applicazione su più intervalli di tipo O(h³). Nel RK al terzo ordine la pendenza media pesata è: dove : |
RK al quarto ordine. La versione del metodo di RK del quarto ordine è data da: Le quattro pendenze sono così
determinate:
In
questo caso l'errore locale è del tipo O(h5),
per cui l'errore totale dovuto all'applicazione su più
intervalli è O(h4)
|
PDE ELLITTICHE. Consideriamo
l'equazione di Laplace nel caso della distribuzione di
temperatura in una lastra conduttrice dovuta alla conduzione
del calore: Si vuole conoscere la temperatura in una serie di punti (nodi) equispaziati che formano un reticolo (la variabile indipendente in più genera un piano a cui punti è associato un valore di temperatura. Una prima approssimazione è rendere discreto l'insieme di punti). In questo schema la coordinata x è la stessa per tutti i nodi che hanno lo stesso indice i, anche se hanno diverso indice j. Allo stesso modo, tutti i nodi che hanno lo stesso indice j hanno la stessa coordinata y. Questo rende estremamente semplice la costruzione di schemi alle differenze finite per il calcolo PDE rispetto a x o y. Ci sono intervalli nelle ascisse positive e intervalli nelle ordinate positive. |
|
Una first central
approximation di
nel punto (i,j) è: in modo simile si ottiene la : Per semplificare ancora assumiamo che , sostituendo le precedenti differenze finite nell'equazione di Laplace otteniamo l'equazione numerica: Con questa tecnica possiamo scrivere un'equazione numerica per ogni nodo in modo da ottenere un sistema lineare la cui soluzione fornisce una soluzione approssimata della PDE (anche perché gli errori di troncamento O(h) sono stati trascurati). Se h è sufficientemente piccolo la soluzione approssimata sarà sempre più accurata (la limitazione è data dal tempo macchina disponibile e dalla memoria del computer). |
|
Tuttavia il sistema ancora non è determinato
se non vengono definite le costanti di integrazioni ovvero
le condizioni al contorno. Ci possono essere due
possibili scelte:
|
|
|
|
PDE PARABOLICHE. Consideriamo
ora il flusso di calore non stazionario e unidimensionale
attraverso una barra di un materiale con una conducibilità
termica k costante. La PDE che descrive la distribuzione di
temperatura T(x,t) nello spazio e nel tempo è:
Qui
occorre fissare le condizioni al contorno :
(valori di temperatura agli estremi della barra dove sono
posti i termostati) che non dipendono dal tempo, e le
condizioni iniziali
all'istante t=0 (distribuzione della temperatura nella barra
all'istante iniziale) Dividiamo la barra in N= L/Δx nodi e consideriamo un'evoluzione temporale discreta, ad intervalli Δt. Ad un certo tempo la derivata seconda spaziale può essere calcolata applicando il metodo esplicito (riprendendo la first central approximation delle derivate seconde qui chiamata CFDA, Central Finite Difference Approximation ) : mentre la derivata prima temporale tra l'istante e nel nodo è : Sostituendo nella PDE parabolica: Da cui esplicitando la temperatura nel nodo all'istante di tempo si trova la temperatura in ogni nodo nell'istante di tempo successivo: Si può dimostrare che il metodo esplicito fornisce soluzioni stabili solo se . Ora consideriamo la PDE parabolica con la derivata seconda spaziale al tempo : ponendo , si può scrivere: E si possono ottenere, in modo implicito, le temperature nei nodi al tempo note quelle al tempo risolvendo un sistema lineare algebrico triangolare senza che ci siano problemi di stabilità. Un miglioramento ulteriore al metodo implicito è quello di Crank-Nicolson che utilizza come derivata spaziale la media delle derivate seconde spaziali al tempo e : Da cui la PDE : e l'equazione generatrice del sistema triangolare algebrico: PDE paraboliche bidimensionali. Ora se si ha un flusso di calore non stazionario bidimensionale attraverso una lastra di un materiale con una conducibilità termica k costante, la PDE è : |
|
Anche in questo caso è possibile applicare la CFDA nel nodo : e risolvere l'equazione trovando in modo esplicito la temperatura nel nodo all'istante di tempo . Ma anche qui il metodo esplicito fornisce soluzioni stabili solo se . La formulazione implicita del metodo è: |
|
con :
E quella di Cranck-Nicolson è: con : |
Sia il metodo implicito che quello di
Cranck-Nicolson sono stabili ma possono portare alla
risoluzione di sistemi algebrici con un numero di equazioni
molto elevato. In figura è mostrata la matrice di
coefficienti per l'applicazione del metodo di Crank-Nicolson
ad una griglia di 4x4 nodi. Col metodo ADI (Alternating Direct Implicit) si effettuano due mezzi passi temporali uno in ogni direzione. Uno in direzione y: e, da quell'istante di tempo , uno in direzione x : |
|
Posto
e
si ottiene la relazione relativa al primo mezzo passo
temporale: In modo simile si ricava l’espressione relativa al secondo mezzo passo temporale: |
|
Queste espressioni generano un sistema lineare più semplice perché puramente tridiagonale che può essere risolto con l'algoritmo di Thomas. |