Oscillatore armonico con costante che varia bruscamente

[pier.paolo]

Ciao a tutti, vi scrivo perché sto cercando un buon metodo numerico per un problema di Cauchy lineare.

Siano $k_1 \ne k_2$ reali positivi, $t_0 > 0$ e consideriamo il seguente problema di Cauchy in $[0,+infty)$:
\begin{equation} \label{eq:cauchy}
\begin{cases}
y(t) + k(t)y''(t) = 0 \\
y(0) = 1/\sqrt{k_1} \\
y'(0) = 0,
\end{cases}
\end{equation}
dove
\[
k(t) =
\begin{cases}
k_1 \hspace{1 cm} 0 \leq t \leq t_0 \\
k_2 \hspace{1 cm} t > t_0
\end{cases}
\]

Poiché $k$ non è continua, il problema di Cauchy non soddisfa le ipotesi del teorema di esistenza e unicità. Tuttavia, è chiaro che esiste un'unica soluzione $\tilde y$ in $[0,t_0]$. Utilizzando i valori $\tilde y(t_0)$ e $\tilde y' (t_0)$ come dati iniziali, è possibile prolungare tale soluzione a un'unica soluzione $y$ di classe $C^1$ del problema (\ref{eq:cauchy}) nel senso seguente: $y$ soddisfa le condizioni iniziali e soddisfa l'equazione differenziale quasi ovunque, precisamente per $t \ne t_0$. Questa è la soluzione a cui farò riferimento da ora in poi.

Un calcolo esplicito dà la seguente forma per $y$:
\[
y(t) =
\begin{cases}
\dfrac{1}{\sqrt{k_1}} \mathrm{cos} \bigg(\dfrac{t}{\sqrt{k_1}} \bigg) \hspace{1 cm} 0 \leq t \leq t_0; \\
A \ \mathrm{cos}\bigg( \dfrac{t}{\sqrt{k_2}} + \phi \bigg) \hspace{1 cm} t > t_0,
\end{cases}
\]
dove $A$ e $\phi$ si ottengono imponendo che $y$ sia di classe $C^1$. Infine, facendo due conti si scopre che l'ampiezza resta costante, i.e. $A = 1/sqrt{k_1}$, se e solo se $t_0 = n \pi \sqrt{k_1}$ con $n \in \mathbb N$.

Ora, mettiamoci nel caso in cui l'ampiezza sia costante. Ho provato a risolvere numericamente il problema con vari metodi a passo fisso, in modo che $t_0 = n \pi \sqrt{k_1}$ sia uno dei nodi. Quello che ho osservato è che, quando $k_2$ è molto maggiore di $k_1$, la seconda ampiezza cresce a dismisura. Ho anche provato a far variare meno bruscamente la funzione $k$, cioè a farla variare linearmente da $k_1$ a $k_2$ in 2-3 passi, ma non ho avuto miglioramenti significativi. Mi rendo conto che la discontinuità di $k$ può creare serie difficoltà ai metodi numerici standard, per cui vi chiedo: mi consigliate un metodo che riesca a trattare efficacemente questa singolarità?

[feddy]

Quindi quello che hai fatto tu è: fare in modo di fermare l'integrazione a $t_0$ e poi farla ripartire con il valore nuovo di $k$, giusto? Una cosa che potresti provare è quella di raffinare il passo dove hai la discontinuità, ma non è una grande idea in realtà. Oltretutto, i metodi a passo costante sono un pessima idea: non risolvono la discontinuità, e ovviamente prendere un passo troppo piccolo è troppo oneroso. Andrebbero usati metodi adattivi, cioè a passo variabile, che però potrebbero ridurre di troppo il passo o cadere in altre situazioni sfortunate, specie se espliciti. E' un classico esempio di problema stiff, e per questo tipo di situazioni sono più indicati metodi impliciti.


In MatLab potresti usare cose di questo tipo: https://it.mathworks.com/help/matlab/ma ... ation.html
Ad esempio, ode15s potrebbe fare al caso tuo, è un solutore di sistemi stiff e potresti passare come options il tempo del salto.

[pier.paolo]

Quindi quello che hai fatto tu è: fare in modo di fermare l'integrazione a $t_0$ e poi farla ripartire con il valore nuovo di $k$, giusto?

Esatto.

Oltretutto, i metodi a passo costante sono un pessima idea: non risolvono la discontinuità, e ovviamente prendere un passo troppo piccolo è troppo oneroso. Andrebbero usati metodi adattivi, cioè a passo variabile, che però potrebbero ridurre di troppo il passo o cadere in altre situazioni sfortunate, specie se espliciti. E' un classico esempio di problema stiff, e per questo tipo di situazioni sono più indicati metodi impliciti.

Immaginavo. Purtroppo stavo modellando delle situazioni in cui le grandezze vengono in genere campionate a passo costante.

In MatLab potresti usare cose di questo tipo: https://it.mathworks.com/help/matlab/ma ... ation.html
Ad esempio, ode15s potrebbe fare al caso tuo, è un solutore di sistemi stiff e potresti passare come options il tempo del salto.

Darò un'occhiata, grazie mille!

[feddy]

Prego! Sei riuscito alla fine ?

[pier.paolo]

In realtà sto approfondendo un po' di teoria, usando i risolutori automatici rischio di non capire bene cosa c'è dietro i problemi. Ci tornerò più in là