FLUIDODINAMICA COMPUTAZIONALE - TERMOFLUIDODINAMICA COMPUTAZIONALE

Ingegneria Meccanica e Gestionale

· TERMOFLUIDODINAMICA COMPUTAZIONALE

· Valerio Marra

· INTRODUZIONE

· EQUAZIONE DI CONSERVAZIONE DELL’ENERGIA

· FLUSSI NON ISOTERMI

· APPROSSIMAZIONE DI BOUSSINESQ

· SCHEMA DI SOLUZIONE NUMERICA

· RISULTATI NUMERICI

INTRODUZIONE

Consideriamo un mezzo continuo costituito da un fluido Newtoniano (i.e. segue la legge di Newton) e un elemento infinitesimale di volume dV al suo interno. Il volume, per quanto piccolo, è comunque tale da avere una dimensione spaziale lineare molto più grande delle distanze intermolecolari e quindi il fluido al suo interno può essere considerato come continuo. Il campo di velocità in un sistema di riferimento cartesiano è rappresentato dal vettore u=ui+vj+wk. Oltre a tale campo nel fluido è definito un campo di pressione p e di temperatura T, il sistema è così completamente determinato dal punto di vista dinamico e termodinamico (le altre proprietà fisiche del fluido possono essere espresse come loro funzione). La conoscenza della distribuzione spazio-temporale dei campi u, p e T è resa possibile, relativamente a date condizioni iniziali e al contorno, dall’applicazione delle leggi fondamentali di conservazione, le quali consentono la definizione di un sistema composto da tre equazioni differenziali alle derivate parziali: una di tipo vettoriale relativa alla conservazione della quantità di moto e due di tipo scalare relative, rispettivamente, alla conservazione della massa e dell’energia.

Molti fluidi, in particolare l’aria e l’acqua, seguono in molte circostanze la legge di Newton (equazione costitutiva che caratterizza il fluido e ne descrive le reazioni ai carichi applicati, i.e. lega linearmente il tensore degli sforzi S al vettore velocità u) e la legge di Fourier (equazione costitutiva che lega linearmente il vettore flusso di calore q alla temperatura assoluta T). Le leggi citate sono descritte in Equazioni di Navier-Stokes, Ing. V. Marra (Ingegneria del Vento, Azioni del Vento sulle Strutture, Dinamica dei Fluidi).

La derivazione e la caratterizzazione teorica delle equazioni di conservazione della massa e della quantità di moto (per fluidi Newtoniani e incomprimibili) sono trattate nell’articolo Equazioni di Navier-Stokes, Ing. V. Marra, al quale si rimanda per maggiori dettagli. La soluzione numerica delle equazioni di Navier-Stokes, alla quale ci riferiremo nel seguito, è invece trattata in Metodi Numerici per le Equazioni di Navier-Stokes, Ing. V. Marra.

La derivazione, la caratterizzazione teorica e la soluzione numerica dell’equazione di conservazione dell’energia è invece oggetto del presente articolo.

EQUAZIONE DI CONSERVAZIONE DELL’ENERGIA

L’equazione di conservazione dell’energia descrive il bilancio energetico all’interno di un volume di riferimento ed è ricavata a partire dalla formulazione integrale, riferita all’unità di massa, del Primo Principio della Termodinamica:

(1)

dove t è il tempo, r la densità, u la velocità, E=e+1/2|u|² è l’energia specifica totale con e energia interna specifica, S è il tensore degli sforzi, q il vettore flusso di calore, a_e l’accelerazione dovuta all’azione di una qualsivoglia forza massica esterna (e.g. l’azione del campo gravitazionale è modellata come a_e=g), n il versore normale uscente dal contorno S del volume di integrazione V e D/Dt è la derivata convettiva o sostanziale (nelle formule analitiche le quantità vettoriali sono indicate come lettere soprassegnate da una freccetta) . In questa formulazione non sono presenti pozzi o sorgenti di energia (e massa). La formulazione lagrangiana della (1) porge:

(2)

Si possono ricavare varie formulazioni dell’equazione (1) in relazione alla grandezza fisica che si vuole assumere come incognita (e.g. la temperatura assoluta T, l’entalpia specifica h oppure l’energia interna specifica e). Scegliendo come incognita il campo scalare temperatura assoluta T e sostituendo nella (2) le rispettive espressioni della legge di Newton e di Fourier e la relazione fondamentale e=c_vT si ottiene la seguente forma euleriana dell’equazione di conservazione dell’energia:

(3)

dove k è il coefficiente di conducibilità termica, c_v la capacità termica a volume costante, p la pressione isotropa (i.e. in assenza di sforzi/deformazioni idrodinamiche) e F¢ è la funzione di dissipazione (esprimente il tasso al quale l’energia è dissipata irreversibilmente in ogni elemento di fluido a causa delle azioni viscose).

FLUSSI NON ISOTERMI

Per flussi non isotermi il set di equazioni di base dell’idrodinamica da risolvere è il seguente:

(4)

(5)

(6)

(7)

(8)

dove I è il tensore unitario, m la viscosità dinamica, D il tensore del tasso di incremento della deformazione nel fluido e l è il secondo coefficiente di viscosità (è l’analogo della prima costante di Lamé derivata dalle relazioni tensione-deformazione nell’ambito della meccanica dei solidi).

Dall’analisi del sistema di equazioni (4)-(8) si evince come si abbiano 9 incognite per 8 equazioni; per raggiungere la chiusura – e quindi la sua risolubilità teorica – il sistema deve essere completato per mezzo di un’equazione di stato. Per la gran parte dei fluidi di interesse ingegneristico possiamo scrivere:

(9)

dove T₀ è la temperatura alla quale r=r₀ e a il coefficiente di espansione volumetrica (o termica). Ricordiamo che:

(10)

(11)

APPROSSIMAZIONE DI BOUSSINESQ

Nella definizione del sistema di equazioni (4)-(11) non è stata fatta nessuna ipotesi riguardo la costanza o meno dei coefficienti m, l, a, k e c_v laddove sono stati introdotti, quindi questo sistema è di validità del tutto generale. Nella maggioranza dei casi di interesse ingegneristico, questo sistema è notevolmente semplificabile e ciò è dovuto in gran parte al fatto che le variabilità dei valori della densità e dei vari coefficienti è dovuta principalmente alle variazioni di temperatura, inoltre tale dipendenza è piuttosto debole. Analizziamo il comportamento della densità r: il coefficiente di espansione volumetrica a per gas e liquidi è dell’ordine di 10^-4 ¸ 10^-3, per variazioni della temperatura non superiori a 10 K la variazione della densità è al massimo dell’1%. La variazione dei coefficienti introdotti a seguito delle piccole variazioni di r è dello stesso ordine e quindi trascurabile.

Tuttavia c’è un’importante eccezione: le variabilità di r nel termine ra_e nell’equazione (6) non può essere trascurata. Questo poiché l’accelerazione risultante da:

(12)

può essere piuttosto grande: più grande, ad esempio, dell’accelerazione dovuta al termine inerziale (u×Ñ)u presente nell’equazione del moto (6). Di conseguenza, tratteremo la densità come constante in tutti i termini dell’equazione (6) eccetto che in quello relativo all’azione delle forze esterne: questa è la cosiddetta approssimazione di Boussinesq. Con DT=T-T₀ si è indicato il valore della variazione di temperatura che interessa il sistema studiato.

Equazioni della Termofluidodinamica nell’approssimazione di Boussinesq

Riscriviamo l’equazione di continuità (4) come:

(4’)

Sulla base delle precedenti osservazioni possiamo affermare che i termini del membro di sinistra sono dell’ordine di a se comparati con quelli del membro di destra. Per questa ragione sostituiamo l’equazione (4¢) con la seguente:

(4’’)

I termini in (5) che sono proporzionali a m e l definiscono gli sforzi dovuti alla viscosità, indicandoli con F possiamo definire il tensore degli sforzi viscosi:

(13)

che diventa in base alla condizione (4¢¢):

(13’)

dove, per le precedenti osservazioni, trattiamo m come una costante. L’espressione per F¢ porge come conseguenza della (4¢¢):

(7’)

In seguito a queste approssimazioni l’equazione del moto (6) diviene:

(6’)

dove n=m/r₀ è la viscosità cinematica. Considerando l’equazione di conduzione del calore (8), possiamo trattare c_v e k come costanti e portarli fuori dall’argomento degli operatori differenziali e ignorare il termine -pÑ×u presente nel membro di destra.

Anche il termine di dissipazione viscosa F¢ può essere ignorato. Come conseguenza dell’equazione (6¢) si ha che il campo di velocità dominante è dell’ordine di |aDTa_eL|^½, dove L è una misura della dimensione lineare del sistema in esame. Relativamente alla (8), il rapporto tra F¢ e il termine dovuto alla conduzione del calore porge:

(14)

dove tale rapporto è per la gran parte dei liquidi dell’ordine di 10^-4¸10^-3, per L»1 cm e |a_e|» g, con g pari all’accelerazione dovuta alla gravità terrestre. Effettuando lo stesso tipo di analisi sul termine che rappresenta la potenza generata dalle forze esterne si ha che il rapporto è dell’ordine di 10^-7¸10^-8. Come risultato di questa analisi e delle approssimazioni fin qui introdotte l’equazione (8) si riduce nella seguente:

(8’)

Riassumendo, il sistema di 5 equazioni per 5 incognite della termofluidodinamica nell’approssimazione di Boussinesq è il seguente:

(4’’)
(6’)
(8’)

I termini presenti in queste equazioni sono noti nella letteratura scientifica come:

·	forza volumetrica di galleggiamento

·	termine di conduzione di calore

·	termine di convezione del calore

Se la forza di galleggiamento è la sola causa del moto allora il termine (u×Ñ)T è detto di convezione libera, se invece essa è trascurabile allora è detto di convezione forzata.

SCHEMA DI SOLUZIONE NUMERICA

L’approssimazione di Bussinesq permette di disaccoppiare agevolmente le equazioni che costituiscono il sistema fondamentale di equazioni della termofluidodinamica (4)-(11). Dal punto di vista risolutivo il disaccoppiamento risulta nella soluzione non più in contemporanea ma in cascata delle equazioni in oggetto. In particolare, nell’ottica dell’approssimazione di Boussineq, lo schema numerico di soluzione – relativo alle equazioni (4¢¢), (6¢) e (8¢) – si esplica nei seguenti quattro passi:

1. inizializzazione dei campi u, p e T;

2. soluzione della (4¢¢) e della (6¢), considerando il campo di temperatura non come incognito ma come assegnato, secondo il metodo descritto in Metodi Numerici per le Equazioni di Navier-Stokes, Ing. V. Marra. Il campo vettoriale u è aggiornato;

3. soluzione della (8¢) secondo uno schema numerico, anch’esso alle differenze finite, nel quale il campo di velocità è considerato non come incognito ma come assegnato dal passo di soluzione precedente. Il campo scalare T è aggiornato;

4. i passi 2. e 3. sono ripetuti fino al soddisfacimento di un assegnato criterio di convergenza.

La variabile incognita temperatura è, in questa approssimazione, uno scalare passivo trasportato dal fluido in movimento (i.e. non esercita direttamente sulle equazioni di conservazione della massa e della quantità di moto alcuna influenza, la sua evoluzione dinamica è quindi completamente determinata da esse). Al fine di catturare questa peculiarità è indicato anche in questo caso l’uso di una tecnica upwind per il trattamento del termine di convezione libera. In questo schema di soluzione numerica, nel quale per la cattura della dinamica del campo di velocità è già utilizzato il metodo di Godunov, è sufficiente che la tecnica upwind utilizzata per la soluzione dell’equazione (8¢) sia un metodo upwind del primo ordine. Ciò consente una determinazione accurata del campo di temperatura evitando l’ulteriore complicazione dell’algoritmo di soluzione e l’aumentare dello sforzo computazionale richiesto.

RISULTATI NUMERICI

In Figura 1 sono mostrati i risultati relativi alla simulazione di una piastra riscaldata centralmente e sovrastata da un fluido caratterizzato da un campo di velocità nullo. In questa configurazione è riscontrabile il solo meccanismo di trasporto del calore per conduzione (diffusione).

Figura 1

In Figura 2 sono mostrati i risultati relativi alla simulazione di una piastra riscaldata centralmente e sovrastata da un fluido caratterizzato da un campo di velocità orizzontale diretto verso destra. In questa configurazione sono riscontrabili sia il meccanismo di trasporto del calore per conduzione (diffusione) che per convezione.

Figura 2

Nelle figure sono mostrate le isoaree del campo di temperatura.