Appunti Fluidodinamica computazionale

Politecnico di Bari - Dipartimento di meccanica, matematica e management

Corso di laurea in ingegneria meccanica LM33

Fluidodinamica computazionale appunti

Primo esonero

Antonio Paolo Fontanella

Sbobinato per esame scritto in formula scritto+tema d’anno, comprendenti i seguenti macro-argomenti:

• Equazioni della fluidodinamica
• Curve caratteristiche
• Classificazione delle EDDP del secondo ordine
• Analisi di stabilità
• Metodi per risoluzione delle equazioni d’onda, calore, Laplace
• Metodi Multigrid

Appunti esclusivamente per uso personale, NON cedere ad altri.

Equazioni da ricordare

Le equazioni a cui siamo più interessati sono quelle della fluidodinamica. Vediamole scritte in forma differenziale:

• Equazione di conservazione della massa:

Abbiamo utilizzato la notazione vettoriale, in modo che sappiamo come scrivere questa equazione nel caso di riferimento (cartesiano, cilindrico, polare, …).

Il procedimento generale per scrivere queste equazioni parte dall’enunciato del trasporto di Reynolds: il volume di controllo scelto può essere fisso o mobile. Se scegliamo un volume mobile dobbiamo fare attenzione ad inserire la velocità v relativa nelle nostre equazioni. Se il volume di controllo è fisso, la velocità relativa v coincide con la velocità assoluta u.

Utilizzando il teorema di Reynolds siamo quindi in grado di scrivere le equazioni di conservazione, prima in forma integrale e poi, sfruttando il teorema di Gauss (teorema della divergenza) in forma differenziale.

Nel caso della massa si parte dal fatto che questa non cambia nel tempo: Possiamo ancora manipolare questa equazione.

Possiamo fare l’ipotesi di flusso incomprimibile:

• Equazione di conservazione della quantità di moto:

Scrivo la forma differenziale:

• Equazione di conservazione dell’energia:

Per chiudere il sistema abbiamo bisogno di varie equazioni: Facendo l’ipotesi di flusso incomprimibile possiamo scrivere una forma semplificata dell’equazione di conservazione della quantità di moto: Posso moltiplicare scalarmente per u ottenendo l’equazione di conservazione dell’energia cinetica:

Per ottenere l’equazione dell’energia interna devo sottrarre membro a membro l’equazione dell’energia totale meno l’equazione dell’energia cinetica.

Riscrivo l’equazione dell’energia totale: Distribuisco l’operazione di divergenza: Possiamo scrivere il termine di reversibilità come prodotto della viscosità per un termine di dissipazione: Posso quindi ottenere l’equazione dell’energia interna:

Voglio ora scrivere l’equazione dell’entropia. Introduco una delle relazioni notevoli della termodinamica: Riscrivo utilizzando un linguaggio più appropriato per i nostri scopi. Abbiamo trovato prima la forma di De/Dt, troviamo ora una forma per Dv/Dt tramite l’equazione di conservazione della massa. Posso quindi sostituire nell’equazione dell’entropia. Ricordiamoci di moltiplicare prima per ρ.

Abbiamo visto cosa accade per il flusso incomprimibile, vediamo cosa accade per un flusso di Euler o fluido ideale. L’equazione della massa rimane la stessa:

Nell’equazione di conservazione della quantità di moto devo trascurare la divergenza degli sforzi viscosi:

Nell’equazione dell’energia dobbiamo trascurare tutti i termini diffusivi. In generale dovremmo trascurare anche il termine relativo alle forze di volume f poiché abbiamo a che fare con un gas. Non è inusuale nella fluidodinamica scrivere le equazioni per il flusso comprimibile utilizzando invece dell’energia interna l’entalpia. Vediamo come fare:

Dall’equazione dell’energia nel caso ideale posso scrivere:

Equazione onda acustica

Possiamo scrivere la generica equazione d’onda come:

Posso scrivere la velocità come somma di due contributi: una velocità indisturbata più il disturbo. Nel caso di onda acustica il flusso si suppone in quiete, quindi la velocità indisturbata è nulla. Posso scrivere anche la densità come somma di due contributi. I disturbi u’ e ρ’ sono delle perturbazioni: sono piccole rispetto alle quantità indisturbate. Se così non fosse si avrebbe un urto.

Riscriviamo l’equazione della massa sostituendo la rappresentazione trovata sopra. Le grandezze indisturbate sono ovviamente costante. Vediamo ora l’equazione della quantità di moto. La scriviamo trascurando i termini viscosi perché le piccole perturbazioni causate dall’onda sonora sono fenomeni isentropici (Infatti abbiamo ricavato l’espressione della velocità del suono in funzione della relazione che c’è tra la pressione e la densità).

Anche la velocità del suono a può essere scritta in funzione della velocità indisturbata e del disturbo. Per scriverla devo usare l’espansione in serie di Taylor. Sostituiamo nell’equazione 1. Ho quindi già ottenuto il risultato, ma non siamo ancora arrivati all’equazione dell’onda. Per trovare questa equazione derivo rispetto al tempo l’equazione (A).

Posso facilmente ottenere la divergenza della variazione locale di u’ applicando l’operatore divergenza all’equazione (B). Sostituisco nell’equazione precedente e ottengo l’equazione dell’onda. Nell’equazione ottenuta il disturbo è un disturbo di densità. Utilizzando le relazioni fisiche possiamo ottenere equazioni con disturbi di pressione o velocità.

Flusso potenziale (o irrotazionale)

In un flusso potenziale la vorticità è nulla. Questa è una condizione più forte rispetto a quella di fluido ideale: le equazioni di Euler che governano il moto di un fluido ideale sono rotazionali, ammettono una soluzione con vorticità diversa da zero. Se la vorticità è nulla si può verificare in qualsiasi sistema di riferimento che la velocità u può essere scritta come gradiente di una funzione scalare φ.

Possiamo anche verificare l’opposto: se la velocità può essere espressa come gradiente di φ il rotore della velocità sarà nullo, quindi il flusso sarà potenziale. Se ipotizziamo che il flusso sia anche incomprimibile si verifica che il laplaciano di φ è nullo.

Curve caratteristiche

Le curve caratteristiche servono a fare una classificazione delle equazioni differenziali. Scriviamo un equazione differenziale del primo ordine nell’incognita u. a, b, c sono i coefficienti dell’equazione differenziale x e y sono le variabili indipendenti.

Il caso considerato è di una equazione quasi lineare in cui i coefficienti dipendono dalle variabili indipendenti e al più dall’incognita. Si dice “quasi lineare” perché i coefficienti contengono derivate di un ordine inferiore rispetto l’ordine massimo. Nel nostro caso l’ordine massimo è il primo, quindi i coefficienti non contengono derivate del primo ordine, contengono al più l’incognita. L’equazione si direbbe lineare se i coefficienti fossero costanti o al più dipenderebbero dalle variabili indipendenti (non dipendono dall’incognita u).

Tornando alla nostra forma, le proprietà della soluzione di questa equazione dipendono dalle proprietà di quella che viene chiamata “superficie integrale”. Una volta che conosco la soluzione di questa equazione, questa dipende da x e y.

Se consideriamo un punto sulla superficie integrale di coordinate (x, y, u) su questa superficie posso individuare un vettore infinitesimo tangente alla superficie stessa che avrà 3 componenti (dx, dy, du). Calcolo ora il differenziale di u. Posso pensare che quest’ultima relazione sia il prodotto scalare tra due vettori, il vettore tangente alla superficie (dx, dy, du) è un vettore di componenti (u_x, u_y, -1). Essendo il prodotto scalare nullo vuol dire che i due vettori sono perpendicolari.

Possiamo fare la stessa considerazione nell’equazione differenziale: questa può essere scritta come prodotto scalare di due vettori. Essendo il vettore (u_x, u_y, -1) perpendicolare ad entrambi i vettori, vuol dire che i vettori (dx, dy, du) e (a, b, c) sono complanari. Essendo complanari posso scrivere:

A ciascun punto sulla superficie integrale posso associare una famiglia di curve caratteristiche tangente alla direzione caratteristica (a, b, c). Queste equazioni della famiglia delle curve caratteristiche costituiscono un sistema di equazioni differenziali ordinarie che individua contemporaneamente le curve caratteristiche e come varia la soluzione lungo la curva caratteristica.

Per come è stata costruita le curve caratteristiche spazzano tutta la superficie integrale, quindi mettendo insieme tutte le curve ottengo per costruzione la soluzione. Per integrare queste equazioni ho bisogno di due condizioni (due costanti di integrazione). Prima di questo, per risolvere l’equazione differenziale di partenza (3) ho bisogno di una condizione al contorno per risolvere il problema ai valori iniziali (o problema di Cauchy). Questo ci dice che le due condizioni per integrare l’equazione delle curve caratteristiche sono in realtà vincolate al fatto che devo soddisfare la condizione al contorno dell’equazione differenziale di partenza, ovvero devo definire una curva dei valori iniziale che rappresenta la condizione al contorno. Devo quindi avere una curva in (x, y) di cui conosco la soluzione.

Ipotizzo che la curva dei valori iniziali sia parametrizzata con la variabile “t”. Devo avere assegnato (devo conoscere) x=x(t) e y=y(t) nonché la soluzione u=u(t). Questa è una delle condizioni che mi servono per integrare le curve caratteristiche t è l’ascissa curvilinea che si muove lungo la curva dei valori iniziali e la mia soluzione non dipenderà più da due parametri ma dipenderà soltanto dall’ascissa curvilinea s che si muove lungo le curve caratteristiche. Sto praticamente riscrivendo la soluzione in termini di un nuovo sistema di riferimento (s, t). Calcolo le curve caratteristiche soddisfacendo per ogni valore di t le seguenti equazioni.

Le prime due equazioni forniscono le equazioni geometriche delle curve caratteristiche, che sono equazioni che mi consentono di descrivere nel piano (x, y) la curva caratteristica. L’ultima equazione (du/ds = c) viene chiamata equazione di compatibilità e indica come varia l’incognita u lungo la curva caratteristica.

La trasformazione di coordinate da (x, y) a (s, t) ha avuto l’effetto di trasformare l’equazione differenziale del primo ordine alle derivate parziali in un sistema di equazioni differenziali ordinarie, in cui, da un lato calcolo le curve caratteristiche e dall’altro calcolo come varia l’incognita lungo la curva caratteristica. Se c fosse zero avrei un’equazione notevole, l’equazione di compatibilità mi direbbe che lungo la curva caratteristica la soluzione u non varia e rimane costante e pari al valore sulla curva iniziale.

Esempio: equazione del primo ordine a coefficienti costanti

Suppongo di conoscere la u ad y=0. Ho quindi già fornito la curva dei valori iniziali che sarà y=0. Se vogliamo vederla con la notazione usata nella teoria, t=0 sarà la nostra curva dei valori iniziali. Devo ora usare le formule scritte prima che nel nostro caso diventano:

Poiché du/ds è nullo, vuol dire che u(s) sarà costante e il suo valore dipendente dal valore sulla curva dei valori iniziali. Devo ora trovare l’espressione delle curve caratteristiche integrando le equazioni differenziali. Per valutare i valori delle costanti devo partire dai valori sulla curva dei valori iniziali dove s=0.

• Per s=0 sappiamo dalla traccia che y=0, quindi k₂=0
• Per s=0 x coincide proprio con la curva dei valori iniziali t

Dobbiamo valutare cosa sia t. Poiché c=cost le curve caratteristiche saranno un fascio di rette parallele di pendenza dy/dx=1/c. Abbiamo ricavato che il valore di u lungo le curve caratteristiche è costante; vediamo a quanto è uguale:

Cerchiamo di capire cosa rappresentano x e y in questo caso. Prendiamo un punto P(x, y) appartenente ad una curva caratteristica. In base a quello che ci dice l’equazione di compatibilità, la soluzione nel punto P(x, y) è proprio quella che abbiamo riportandoci sulla curva dei valori iniziali lungo la curva caratteristica che passa per questo punto. La curva caratteristica propaga l’informazione sulla curva dei valori iniziali.

Esempio 2: equazione quasi lineare

Conoscendo l’esempio precedente possiamo subito verificare che la famiglia delle curve caratteristiche rappresentano un fascio di rette che non sono parallele perché dipendono da u. du/ds = 0 quindi le curve caratteristiche sono delle rette; se fosse diverso da zero non sarebbero delle rette. Poiché la pendenza non è costante, ma dipende da u(x, y) allora le curve caratteristiche non sono parallele.

Possiamo ottenere la soluzione del nostro sistema se è possibile effettuare le trasformazioni scritte prima. Queste trasformazioni si hanno quando noi possiamo scrivere il seguente sistema:

Posso anche scrivere il forma matriciale:

Quello che sto facendo è praticamente un cambiamento di coordinate da (x, y) a (s, t) che è UNIVOCAMENTE DETERMINATO se per valori di s sufficientemente piccoli, ovvero in prossimità della curva dei valori iniziali, lo Jacobiano della trasformazione è diversa da zero. Se lo Jacobiano J della trasformazione è nullo accade che:

Sappiamo che x_s=a e che y_s=b, quindi:

Ci chiediamo se esista una curva lungo la quale lo Jacobiano della trasformazione sia nullo; Cerchiamo di valutare questa curva in maniera parametrica e il parametro di questa curva lo chiamo u_t. Utilizzo anche l’equazione differenziale. Divido membro a membro le due equazioni. Considero la condizione che abbiamo quando J=0.

Sostituiamo e otteniamo:

La trasformazione da (x, y) a (s, t) lungo la curva caratteristica non è univocamente determinata, quindi posso pensare di valutare l’equazione delle curve caratteristiche imponendo che lo Jacobiano della trasformazione sia nullo. Più in generale possiamo dire che le curve caratteristiche sono il luogo dei punti in cui le derivate di ordine massimo dell’equazione differenziale non sono univocamente determinate e quindi, a seconda dell’ordine dell’equazione definiremo le curve caratteristiche andando a verificare il luogo dei punti lungo cui la derivata di ordine massimo non è univocamente determinata.

Possiamo esprimere tutto questo anche in un altro modo. Supponiamo sempre che t sia il parametro che utilizziamo per valutare il luogo dei punti lungo cui la derivata di ordine massimo non è univocamente determinata. Dobbiamo capire se le derivate prime u_x(t) e u_y(t) possano essere univocamente determinate a fronte del cambio di riferimento. Scrivo l’equazione differenziale: u sarà la soluzione dell’equazione di partenza. Posso sempre scrivere la derivata di u rispetto a t calcolata con la chain rule (regola di derivazione a catena). Posso mettere a sistema le equazioni (4) e (5).

Posso valutare da questa trasforma le due derivate u_x e u_y, che in questo caso sono le derivate di ordine massimo, se riesco a invertire il sistema sopra. Come abbiamo detto prima, il luogo dei punti lungo cui le derivate di ordine massimo (in questo caso le prime) possono non essere univocamente determinate coincide con l’equazione della curva caratteristica. Questa è la maniera che useremo per definire le curve caratteristiche e anche le proprietà delle equazioni differenziali.

Classificazione delle EDDP del secondo ordine

Consideriamo un’equazione differenziale del secondo ordine nella variabile φ. Interessandoci principalmente i termini con la derivata di ordine massima, possiamo riscrivere l’equazione nella seguente forma:

Per semplicità chiamo pongo i termini a destra dell’equazione uguali a una funzione h. In maniera analoga a quanto visto per le EDDP del primo ordine posso costruire il seguente sistema:

Vediamo come si arriva a questo sistema. Ci interessa valutare una curva x(t) e y(t) lungo cui le derivate non sono definite in maniera univoca. Poiché la valutazione viene effettuata sulla derivata seconda, si ipotizza che le derivate prime siano valutabili e continue nel dominio di calcolo. Dobbiamo quindi valutare rispetto queste variabili la derivata in t.

Il sistema sarà quindi costituito dall’equazione differenziale di partenza (2.1) e dalle equazioni (2.2) e (2.3) scritte sopra. La trasformazione è possibile ed è unica soltanto se il determinate della matrice dei coefficienti è diverso da zero. Se, invece, il determinante è nullo otteniamo come risultato un’equazione che definisce, anche in questo caso, le famiglie di curve caratteristiche.

Calcoliamo il determinante:

dy/dx rappresenta la pendenza della curva, parametrizzata con t, nel piano x-y. Risolvendo l’equazione di secondo grado nell’incognita dy/dx otteniamo, come nel caso delle equazioni alle derivate parziali del primo ordine, la pendenza della curva caratteristica. Avremo due famiglie di curve caratteristiche che dipendono dal discriminante:

• Se il discriminante è maggiore di zero avremo due famiglie di curve caratteristiche reali. In questo caso l’equazione si dice iperbolica
• Se il discriminante è nullo avremo due radici coincidenti e quindi un’unica famiglia di curve caratteristiche nel piano fisico. L’equazione si dice parabolica
• Se il discriminante è minore di zero le radici saranno complesse e coniugate, non potremo tracciare curve caratteristiche nello spazio fisico ma in quello complesso. L’equazione si dice ellittica