A-stabilità e passo di discretizzazione per un problema differenziale

Messaggioda feddy » 07/12/2017, 21:26

Buonasera a tutti,

dato questo sistema differenziale

$ { ( u''(x)=-9000u'(x)-5000(sin(u(x))+u(x)) ),( u(0)=pi ),( u'(0)=0 ):} , x \in (0,10) $

Lo si risolva con un opportuno metodo esplicito, giustificando la scelta del passo di discretizzazione.



Ho provato, visto che il testo chiedeva un metodo esplicito, a considerare per semplicità Eulero: $y_{n+1}=y_{n}+k*f(t_n,y_n)$, $k$ passo di discretizzazione dato da $k=t_{f}/(ts)$, dove $ts$ è il numero di passi temporali (timesteps) e $t_{f}$ è il tempo finale.

L'implementazione è molto semplice e risulta
Codice:
tstar=10 %tempo finale
ts=200; %passi temporali
k=tstar/ts;
t=linspace(0,tstar,ts);
y0=[pi,0]';
y=zeros(2,ts);
y(:,1)=y0;

f=@(t,y) [y(2);-9000*y(2)-5000*(sin(y(1))+y(1))]; %definizione funzione

for n=1:ts-1
    y(:,n+1)=y(:,n)+k*f(t(n),y(:,n));
end


Plottando la soluzione, risulta evidentemente sbagliata, come si vede

Immagine


Mi pare di intuire che questo sia dovuto al fatto che il passo temporale è troppo grande, in quanto ho visto stamane a lezione che per eulero esplicito e altri metodi NON tutti i passi temporali sono accettabili, e in particolare c'è una restrizione sul passo di discretizzazione $k$.
Pertanto, con lo stesso codice di prima, ma prendendo $k=1e-8$, ottengo

Immagine

che mi sembra accettabile come soluzione numerica, visto l'andamento a tangente orizzontale per un intorno di $x=0$.

La vera domanda è: come scegliere in modo accurato questo $k$? Esiste un modo "operativo".


L'ho visto fare per il caso banale $y'(t)=lambda y(t), y(t_0)=y_0$ di cui però conosco la soluzione analitica e diciamo che come esempio didattico funziona e i conti sono molto facili.
Da quel che mi pare di capire la strategia è (assumendo sempre di aver implementato correttamente il codice) notare che la soluzione compie oscillazioni troppo elevate, o comunque assume comportamenti errati, e perciò si prende un passo temporale più piccolo, in modo da entrare nella regione di assoluta stabilità.

Ringrazio chiunque abbia la cortesia di cimentarsi :smt026
Avatar utente
feddy
Moderatore
Moderatore
 
Messaggio: 1627 di 5934
Iscritto il: 26/06/2016, 00:25
Località: SISSA

Torna a Analisi Numerica e Ricerca Operativa

Chi c’è in linea

Visitano il forum: Nessuno e 1 ospite