Risoluzione numerica equazione dinamica del moto
Buonasera a tutti,
Mi sono imbattuto ultimamente in un problema che non so bene come affrontare. Si tratta della soluzione numerica di una equazione differenziale alle derivate parziali deputata a descrivere la vibrazione di un cavo teso soggetto all’azione del vento e della pioggia (rain-wind induced vibration).
L’equazione in oggetto presenta la forma sotto riportata dove m, T e c sono rispettivamente la massa per unità di lunghezza del cavo, la trazione in esso agente e la costante di smorzamento strutturale intrinseco:
$m(del^2w)/(delt^2)+T(del^2w)/(delx^2)+c(delw)/(delt)=F(x,t)$
La forzante F(x,t) presenta un’espressione molto complessa legata al problema specifico dell’interazione pioggia-vento, dipendente da numerosi parametri comunque noti o determinabili sperimentalmente.
In bibliografia, la suddetta equazione viene risolta caso per caso, ricorrendo al metodo numerico Runge-Kutta-Fehlberg (RKF45). Ricercando il suddetto metodo ne ho sostanzialmente compresa la filosofia e le modalità di applicazione. Ho tuttavia sempre riscontrato che il metodo è applicato ad una EDO nella forma seguente:
$dy/dt=f(y(t),t)$
Nel caso di una EDO di ordine superiore al primo, questo non è un problema in quanto è possibile ricondursi ad un sistema composto da un numero finito di EDO nella forma richiesta.
Mi trovo tuttavia in difficoltà nel capire come sia possibile ridurre l’equazione del moto in oggetto ad una forma tale da potervi applicare il metodo RKF45.
Non ho trovato alcun esempio di soluzione del problema, bensì solamente i risultati. Ragionando, forse potrei suddividere il cavo in un certo numero di elementi ed applicare l’equazione a ciascun nodo. In tal modo questa dovrebbe potersi ridurre ad una EDO del secondo ordine dipendente solamente dal tempo e, di conseguenza, ad un sistema di due EDO del primo ordine. Ad ogni modo, se così fosse, non mi è molto chiaro come dovrei procedere.
Nella speranza che qualcuno possa fornirmi un aiuto, vi ringrazio anticipatamente
Mi sono imbattuto ultimamente in un problema che non so bene come affrontare. Si tratta della soluzione numerica di una equazione differenziale alle derivate parziali deputata a descrivere la vibrazione di un cavo teso soggetto all’azione del vento e della pioggia (rain-wind induced vibration).
L’equazione in oggetto presenta la forma sotto riportata dove m, T e c sono rispettivamente la massa per unità di lunghezza del cavo, la trazione in esso agente e la costante di smorzamento strutturale intrinseco:
$m(del^2w)/(delt^2)+T(del^2w)/(delx^2)+c(delw)/(delt)=F(x,t)$
La forzante F(x,t) presenta un’espressione molto complessa legata al problema specifico dell’interazione pioggia-vento, dipendente da numerosi parametri comunque noti o determinabili sperimentalmente.
In bibliografia, la suddetta equazione viene risolta caso per caso, ricorrendo al metodo numerico Runge-Kutta-Fehlberg (RKF45). Ricercando il suddetto metodo ne ho sostanzialmente compresa la filosofia e le modalità di applicazione. Ho tuttavia sempre riscontrato che il metodo è applicato ad una EDO nella forma seguente:
$dy/dt=f(y(t),t)$
Nel caso di una EDO di ordine superiore al primo, questo non è un problema in quanto è possibile ricondursi ad un sistema composto da un numero finito di EDO nella forma richiesta.
Mi trovo tuttavia in difficoltà nel capire come sia possibile ridurre l’equazione del moto in oggetto ad una forma tale da potervi applicare il metodo RKF45.
Non ho trovato alcun esempio di soluzione del problema, bensì solamente i risultati. Ragionando, forse potrei suddividere il cavo in un certo numero di elementi ed applicare l’equazione a ciascun nodo. In tal modo questa dovrebbe potersi ridurre ad una EDO del secondo ordine dipendente solamente dal tempo e, di conseguenza, ad un sistema di due EDO del primo ordine. Ad ogni modo, se così fosse, non mi è molto chiaro come dovrei procedere.
Nella speranza che qualcuno possa fornirmi un aiuto, vi ringrazio anticipatamente
Risposte
Sei sulla buona strada, quello che si fa in questi casi è fare prima una discretizzazione spaziale, trasformando il problema continuo in un problema discreto con un numero finito di variabili, e poi per ciascuno di essi si scrive un'ODE normale a cui si applica un qualunque metodo di avanzamento in tempo.
Nella versione più ignorante di tutte, sostituisci alle derivate in tempo e in spazio delle differenze finite e poi parti da \(u_h^0 = u_h(t_0)\) e calcoli \(u_h^1\), poi con quello \(u_h^2\) e così via. Prova e vediamo che succede.
P.s. non serve passare per il sistema, basta che per la derivata del secondo ordine metti una differenza finita adeguata.
Nella versione più ignorante di tutte, sostituisci alle derivate in tempo e in spazio delle differenze finite e poi parti da \(u_h^0 = u_h(t_0)\) e calcoli \(u_h^1\), poi con quello \(u_h^2\) e così via. Prova e vediamo che succede.
P.s. non serve passare per il sistema, basta che per la derivata del secondo ordine metti una differenza finita adeguata.
Grazie per la riposta.
Da quel che ho capito leggendo in giro, l'integrazione del tempo è effettuata mediante il metodo RKF45 mentre quella nello spazio discretizzando il cavo in n elementi finiti. Tuttavia non mi è chiaro nella pratica come procedere dopo aver discretizzato il cavo, ossia come ridurmi ad un problema risolvibile mediante il metodo RKF45.
Da quel che ho capito leggendo in giro, l'integrazione del tempo è effettuata mediante il metodo RKF45 mentre quella nello spazio discretizzando il cavo in n elementi finiti. Tuttavia non mi è chiaro nella pratica come procedere dopo aver discretizzato il cavo, ossia come ridurmi ad un problema risolvibile mediante il metodo RKF45.
Lascia perdere Runge-Kutta, puoi farlo col metodo che vuoi. Ti ho già detto "nella pratica come procedere", leggi ed esegui 
Proverò a risolverlo come da te indicato. Comunque presumo che se in bibliografia viene comunemente impiegato il metodo Runge-Kutta ci sia un motivo.
Non sono un matematico e, sfortunatamente, sono un po' arrugginito in analisi numerica. Mi trovo a dover risolvere il suddetto problema senza stare a perderci troppo tempo (il fatto che a mezzanotte ci stia ancora ragionando non è il massimo
). Speravo in un spiegazione un po' più completa. A livello intuitivo mi è chiaro il ragionamento ma come una volta ho sentito dire "tra teorica e pratica, in teoria non c'é differenza, in pratica sì." ...
Non sono un matematico e, sfortunatamente, sono un po' arrugginito in analisi numerica. Mi trovo a dover risolvere il suddetto problema senza stare a perderci troppo tempo (il fatto che a mezzanotte ci stia ancora ragionando non è il massimo
). Speravo in un spiegazione un po' più completa. A livello intuitivo mi è chiaro il ragionamento ma come una volta ho sentito dire "tra teorica e pratica, in teoria non c'é differenza, in pratica sì." ...
Il tuo libro usa Runge-Kutta perché *tutti* usano RK4, tant'è che lo chiamano The Runge-Kutta. La spiegazione è già completa com'è, manca solo da scrivere le equazioni, che corrisponde a scegliere come fare le discretizzazioni in tempo e in spazio. Per un problema 1D "tutte le discretizzazioni in spazio sono uguali", quindi metti pure differenze finite.
Sto provando a risolvere il problema come suggerito ma non mi è molto chiaro come procedere. L'equazione alle differenze finite dovrebbe essere scritta come segue:
$(w(t_(i+1),x_j)-2w(t_i,x_j)+w(t_(i-1),x_j))/(tau^2)+T/m (w(t_i,x_(j+1))-2w(t_i,x_j)+w(t_i,x_(j-1)))/(h^2)+c/m (w(t_(i+1),x_j)-w(t_(i-1),x_j))/(2*tau)=(F(t_i,x_j))/m$
Le condizioni al contorno sono $w(t_i,x_0)=w(t_i,x_N)=0$ per ogni $t_i$, ossia che gli estremi del cavo sono fissi.
Le condizioni iniziali sono invece $w(t_0,x_j)=0$ per ogni $x_j$, ossia che il cavo è inizialmente rettilineo (segmento congiungente i due punti fissi).
Se adesso esplicito l'equazione nel nodo (1) ho dentro le incognite $w(t_2,x_1), w(t_1,x_2), w(t_1,x_1)$. La prima è quanto sto ricercando, ossia la soluzione all'istante successivo, ma le altre due sono incognite. Non mi è chiaro per quale motivo non è necessario scrivere il sistema.
Diciamo che una linea guida da seguire sarebbe utile.
$(w(t_(i+1),x_j)-2w(t_i,x_j)+w(t_(i-1),x_j))/(tau^2)+T/m (w(t_i,x_(j+1))-2w(t_i,x_j)+w(t_i,x_(j-1)))/(h^2)+c/m (w(t_(i+1),x_j)-w(t_(i-1),x_j))/(2*tau)=(F(t_i,x_j))/m$
Le condizioni al contorno sono $w(t_i,x_0)=w(t_i,x_N)=0$ per ogni $t_i$, ossia che gli estremi del cavo sono fissi.
Le condizioni iniziali sono invece $w(t_0,x_j)=0$ per ogni $x_j$, ossia che il cavo è inizialmente rettilineo (segmento congiungente i due punti fissi).
Se adesso esplicito l'equazione nel nodo (1) ho dentro le incognite $w(t_2,x_1), w(t_1,x_2), w(t_1,x_1)$. La prima è quanto sto ricercando, ossia la soluzione all'istante successivo, ma le altre due sono incognite. Non mi è chiaro per quale motivo non è necessario scrivere il sistema.
Diciamo che una linea guida da seguire sarebbe utile.
Stai facendo bene! Questa è una equazione del secondo ordine in tempo, quindi dovresti avere due condizioni iniziali, una sulla posizione al tempo iniziale, che hai già scritto, e una sulla sua derivata. Quest'ultima ti serve per ricavare la posizione al primo istante temporale, così puoi finalmente calcolare la soluzione dall'istante 2 fino alla fine.
Ci stavo provando. Considerando velocità iniziale nulla, ossia in termini di differenze finite:
$w(t_1,x_j)-w(t_-1,x_j)=0$
Scrivo l'equazione alle differenze finite nel punto $(t_0,x_j)$ e, sfruttando le suddette condizioni iniziali, ricavo la soluzione $w(t_1,x_j)$, ossia lo spostamento di tutti i nodi del cavo all'istante t1.
Poi scrivo l'equazione alle differenze finite nel punto $(t_1,x_j)$ e, come prima, ricavo la soluzione $w(t_2,x_j)$.
Procedo così per i successivi istanti temporali. In questo senso non devo risolvere il sistema perché risolvo tutto in cascata. Corretto?
In secondo luogo, la valutazione degli intervalli $tau$ e h come è opportuno effettuarla al fine di garantire convergenza?
$w(t_1,x_j)-w(t_-1,x_j)=0$
Scrivo l'equazione alle differenze finite nel punto $(t_0,x_j)$ e, sfruttando le suddette condizioni iniziali, ricavo la soluzione $w(t_1,x_j)$, ossia lo spostamento di tutti i nodi del cavo all'istante t1.
Poi scrivo l'equazione alle differenze finite nel punto $(t_1,x_j)$ e, come prima, ricavo la soluzione $w(t_2,x_j)$.
Procedo così per i successivi istanti temporali. In questo senso non devo risolvere il sistema perché risolvo tutto in cascata. Corretto?
In secondo luogo, la valutazione degli intervalli $tau$ e h come è opportuno effettuarla al fine di garantire convergenza?
"Casti91":
In questo senso non devo risolvere il sistema perché risolvo tutto in cascata. Corretto?
Corretto! Questo che tu chiami "a cascata" è la proprietà del metodo che hai scelto di essere esplicito. Ci sono molte altre scelte possibili di approssimazione della derivata temporale, e alcune di esse non permettono di ricavare i singoli valori individualmente; per quelle è necessario scrivere un sistema lineare e risolverlo.
"Casti91":
In secondo luogo, la valutazione degli intervalli $tau$ e h come è opportuno effettuarla al fine di garantire convergenza?
Questa è un'equazione iperbolica, quindi ci sarà da rispettare una condizione CFL. Questo dipende dal problema specifico e dal metodo numerico scelto, non ti posso dare una risposta universale. Ti posso dire però che il metodo che hai usato adesso è un metodo basico e non ha grandi proprietà di stabilità, per questo di solito si usano metodi più sofisticati, come il già citato Runge-Kutta di ordine 4. La cosa importante è che adesso che hai capito come si scrivono la discretizzazione spaziale e temporale ti sarà semplice sostituire le derivate in tempo con le corrispondenti espressioni del metodo Runge-Kutta.
Grazie dell'aiuto. Non mi è del tutto chiaro come vado ad applicare ora Runge-Kutta, è la prima volta che faccio questa cosa è la faccenda mi è un po' nebulosa. Non mi sembra di essere nella forma per la quale è prevista l'applicazione del metodo RKF45. Provo a spiegarmi meglio a costo di dire qualche cavolata: da quel che ho capito, il metodo si applica ad una EDO nella forma $dy/dt=f(t,y)$. La mia equazione alle differenze finite non mi sembra in questa forma e non riesco a riconoscere la funzione f.
Inoltre, le derivate nello spazio le tratto alle differenze finite o passo anche qui a Runge-Kutta?
Inoltre, le derivate nello spazio le tratto alle differenze finite o passo anche qui a Runge-Kutta?
"Casti91":
Grazie dell'aiuto. Non mi è del tutto chiaro come vado ad applicare ora Runge-Kutta, è la prima volta che faccio questa cosa è la faccenda mi è un po' nebulosa. Non mi sembra di essere nella forma per la quale è prevista l'applicazione del metodo RKF45. Provo a spiegarmi meglio a costo di dire qualche cavolata: da quel che ho capito, il metodo si applica ad una EDO nella forma $dy/dt=f(t,y)$. La mia equazione alle differenze finite non mi sembra in questa forma
Adesso devi usare un piccolo trucco con cui converti un'equazione del secondo ordine scalare in una del primo ordine vettoriale.
Scrivi
\[
v = w', \quad
\vec u = \begin{pmatrix}
w \\
v
\end{pmatrix}
\]
e riscrivi l'equazione differenziale nella nuova incognita vettoriale \(\vec u\) come
\[
\vec u' = F(\vec u)
\]
"Casti91":
e non riesco a riconoscere la funzione f.
\(f\) o \(F\) è tutto ciò che non è la derivata.
"Casti91":
Inoltre, le derivate nello spazio le tratto alle differenze finite o passo anche qui a Runge-Kutta?
In spazio hai un problema ai limiti, non un problema ai valori iniziali, quindi non si può usare Runge-Kutta. Rimani con differenze finite che va bene.
Ok, appena ho un po' di tempo mi metto a fare i conti e vedo cosa ottengo.
Ragionando sulla valutazione dell'intervallo h con cui discretizzare il dominio geometrico, stavo pensando se possa avere senso assumere il numero di nodi minimo per descrivere l'ultimo modo di vibrare a cui posso essere interessato. In altre parole, se sono interessato ai modi 1 e 2 mi sarebbero sufficienti 2 nodi interni (ossia 3 elementi). Assumendo ad esempio 10 elementi, ossia 9 nodi interni sono in grado di cogliere fino al nono modo che è fin troppo.
Può essere un scelta sensata oppure posso comunque incorrere in errori grossolani?
Ragionando sulla valutazione dell'intervallo h con cui discretizzare il dominio geometrico, stavo pensando se possa avere senso assumere il numero di nodi minimo per descrivere l'ultimo modo di vibrare a cui posso essere interessato. In altre parole, se sono interessato ai modi 1 e 2 mi sarebbero sufficienti 2 nodi interni (ossia 3 elementi). Assumendo ad esempio 10 elementi, ossia 9 nodi interni sono in grado di cogliere fino al nono modo che è fin troppo.
Può essere un scelta sensata oppure posso comunque incorrere in errori grossolani?
Non so come dimostrarlo, ma ho la sensazione che sia una scelta sbagliata. L'obiezione che mi viene è che tu non dovresti assumere niente sulla soluzione prima di averla calcolata, quindi io ti consiglio di mettere tanti nodi interni quanti il tuo pc può permettersi di calcolare in tempi ragionevoli. Siccome è un problema 1D, se lo scrivi in matlab male puoi facilmente gestire centinaia di nodi; se lo scrivi in matlab bene puoi arrivare a varie migliaia di nodi interni, sempre dipendendo dall'intervallo temporale.
Sfortunatamente non ho modo di usare Matlab. Posso ricorrere solamente ad Excel. Un approccio che utilizzo di solito per modellare problemi di statica o dinamica dei cavi mediante software ad elementi finiti e che mi ha condotto a buoni risultati è la suddivisione del cavo in circa una ventina di elementi. Non ho idea se possa ripetere qualcosa di simile in questo contesto.
Deciditi: o puoi usare solo excel o puoi usare "software a elementi finiti". In ogni caso, dipende da qual è il tuo fine ultimo. Secondo me con le differenze finite te la puoi cavare senza troppi sforzi, mentre se vuoi discretizzare in spazio con elementi finiti c'è da fare qualche conto in più.
Non mi sono spiegato: il problema in oggetto lo posso risolvere solo con Excel. Per problemi simili mi è capitato di usare software ad elementi finiti in cui una discretizzazione spaziale in circa 20 elementi mi ha condotto a buoni risultati (in media un cavo di 100m discretizzato in elementi da 5m). Mi chiedevo se un approccio simile possa condurre a qualcosa di sensato. Chiaramente non ho la possibilità di usare migliaia di elementi con Excel. Volevo capire se ci sono delle linee guida che mi possono aiutare a definire un numero ottimale di nodi, e procedere quindi più speditamente, oppure devo andare per tentativi e vedere cosa succede.
La linea guida è quella che ti ho dato: più nodi => minor errore. Mi rendo conto che non ti sto dando esattamente le risposte più utili possibile, però è che davvero a priori non c'è una scelta ottimale di nodi, di disposizione o altro. Ci sono tecniche avanzate che fanno cose simili ma direi che in questo caso è come usare un cannone per uccidere una zanzara, e comunque non sono fattibili con excel.
Accedi a tutti gli appunti
Tutor AI: studia meglio e in meno tempo