Metodi numerici
per equazioni differenziali ordinarie - laboratorio
Docente:
Beatrice Meini
Lezione 8: Problemi al contorno risolti con il metodo alle differenze finite e BVP4C di Matlab
-
La function [t,y]=bvp_df(A,q,Ba,Bb,r,a,b,l,h)
applica il metodo alle differenze finite, basato sul metodo del punto di mezzo, al BVP lineare x'(t)=A(t)x(t)+q(t), a< t < b,
Ba*x(a)+Bb*x(b)=r, dove le funzioni A(t), q(t) dipendono dal parametro l, e h e' il passo di discretizzazione.
-
Consideriamo il BVP lineare x'(t)=A(t)x(t)+q(t), 0 < t < 1,
B0*x(0)+B1*x(1)=c, dove le funzioni A(t), q(t) dipendono da un parametro l,
visto alla lezione precedente (Example 7.2 del libro Ascher-Petzold).
Le functions [B0, B1,
b]=bc(l) , z=A(t,l),
w=q(t,l) definiscono le matrici B0, B1, il vettore b, e le
funzioni A(t) e q(t).
- Risolvere il problema con il metodo alle differenze finite, e confrontare il risuttato con la soluzione esatta, che e' la funzione u=(exp(l*(t-1))+ exp(2*l*(t-1))+exp(-l*t)) / (2+exp(-l))+cos(pi*t);
Ad esempio, dare i comandi l=10; h=0.01; [t,x]=bvp_df(@A, @q, B0, B1, b, 0, 1, l, h); u=(exp(l*(t-1))+ exp(2*l*(t-1))+exp(-l*t)) / (2+exp(-l))+cos(pi*t); max(abs(u-x(1,:)))
- Osservare sperimentalmente che la convergenza a zero dell'errore e' quadratica in h.
- Osservare che si riescono a risolvere problemi in cui il metodo di shooting falliva.
- La function bvp4c di Matlab risolve problemi al contorno utilizzando un metodo con ordine di convergenza 4, basato su tecniche di collocazione.
Il file BVP_TUTORIAL.zip contiene un tutorial dove si spiega il metodo utilizzato, il funzionamento della function e si fanno alcuni esempi.
L'uso di bvp4c e' diverso da quello delle function ode** per risolvere IVP.
L'uso base richiede in input, oltre alla funzione f(t,x) che definisce l'equazione differenziale, la funzione g(u,v) che definisce le condizioni al contorno, una discretizzazione iniziale dell'intervallo, e una approssimazione iniziale della soluzione. La discretizzazione e l'approssimazione iniziale vengono date utilizzando la function bvpinit.
- ESEMPIO: La function z=es_bvp(t,x) definisce l'equazione differenziale x'(t)=A(t)x(t)+q(t) del punto precedente, dove il parametreo l e' fissato a l=10. La function y=es_bvp_bc(u,v) definisce le condizioni al contorno B0*x0+B1*x1-b=0, relative al parametro l=10.
Come discretizzazione iniziale scelgo 5 punti equispaziati dell'intervallo [0 5], e come approssimazione iniziale scelgo la funzione costante uguale a 1, nelle 3 componenti. Questa inizializzazione si realizza con il comando solinit = bvpinit(linspace(0,1,5),ones(3,1)).
Il BVP viene risolto dando il comando
sol = bvp4c(@es_bvp,@es_bvp_bc,solinit)
L'output sol e' uno "struct array", tale che sol.x e' la discretizzazione dell'intervallo, sol.y e' l'approssimazione della soluzione, sol.yp e' la derivata.
Qual e' l'errore dell'approssimazione?
- Provare a usare gli esempi del tutorial (i codici sono gia' pronti)
- Provare a risolvere i BVP della lezione scorsa.