Simulazione ed analisi dei fenomeni di trasporto

Simulazione ed analisi dei
fenomeni di trasporto per
moto intorno ad un disco
Autori:
Ubaldo Fierro 0622200082
Fabio Fulchini 0622200110
Mirko Granese 0622200084
Alessandro Visconti 0622200085
Indice
Indice delle figure ......................................................................................................... I
Indice tabelle ............................................................................................................... II
Introduzione ................................................................................................................. 1
Teoria del moto intorno ad oggetti ............................................................................ 2
1.1.
Generalità sul moto dei fluidi ................................................................................................ 2
1.2.
Moto intorno ad oggetti ......................................................................................................... 4
1.3.
Case study: “moto intorno ad un disco” ................................................................................ 6
Uso di FlexPDE per la risoluzione del moto attorno ad un disco ......................... 10
2.1.
Introduzione......................................................................................................................... 10
2.2.
FlexPDE .............................................................................................................................. 11
2.3.
Costruzione dominio del problema ..................................................................................... 13
2.3.1.
Il sistema di coordinate ................................................................................................ 13
2.3.2.
Il dominio ..................................................................................................................... 13
2.4.
Individuazione variabili del problema ................................................................................. 15
2.5.
Definizione principali parametri del problema ................................................................... 16
2.6.
Introduzione equazioni che descrivono i fenomeni di trasporto ......................................... 20
2.7.
Scelta delle condizioni al contorno ..................................................................................... 21
2.8.
Visualizzazione dei risultati ottenuti ................................................................................... 23
Analisi fluidodinamica .............................................................................................. 24
3.1.
Introduzione......................................................................................................................... 24
3.2.1.
Analisi del profilo di velocità: caso standard .................................................................. 24
3.2.2.
Analisi del profilo di velocità: componente vz ............................................................. 25
3.2.3.
Analisi del profilo di vz su diverse sezioni................................................................... 27
3.2.4.
Analisi relativa a vr ...................................................................................................... 30
3.3.
Analisi della pressione......................................................................................................... 31
3.4.
Introduzione allo studio delle tensioni e del fattore di attrito.............................................. 32
3.5.
Caratterizzazione del coefficiente di attrito per variazione di Re ....................................... 33
3.6.
Caratterizzazione del coefficiente di attrito per variazione del raggio e velocità ............... 35
3.7. Caratterizzazione del coefficiente di attrito per variazione del solo raggio a parità di
velocità v0 ...................................................................................................................................... 37
Trasporto energia ...................................................................................................... 39
4.1.
Introduzione......................................................................................................................... 39
4.2.
Analisi del profilo di temperatura ....................................................................................... 40
4.3.
Analisi del flusso di calore: qz ............................................................................................. 43
4.4.
Analisi del flusso di calore: qr ............................................................................................. 44
4.5.
Analisi del flusso di calore .................................................................................................. 45
4.6.
Effetti della variazione del numero di Reynolds ................................................................. 45
4.7.
Effetti della variazione del raggio del disco ........................................................................ 48
Conclusioni ................................................................................................................. 50
Bibliografia................................................................................................................. 51
Indice delle figure
Figura 1 Diagramma Cd-Re. ................................................................................................................ 5
Figura 2 Correlazione JH-Re per una lastra investita. .......................................................................... 6
Figura 3 Moto intorno ad un disco ....................................................................................................... 6
Figura 4 Tipi di mesh: strutturata-strutturata a blocchi- non strutturata. ........................................... 10
Figura 5 Sezioni FlexPDE.................................................................................................................. 12
Figura 6 Disco come sodido di rotazione........................................................................................... 13
Figura 7 Sezione COORDINATES. .................................................................................................. 13
Figura 8 Dominio del problema. ........................................................................................................ 14
Figura 9 Mesh disco investito. ........................................................................................................... 15
Figura 10 Sezione VARIABLES. ...................................................................................................... 15
Figura 11 Sezione DEFINITIONS. .................................................................................................... 16
Figura 12 Feature. .............................................................................................................................. 18
Figura 13 Sezione EQUATIONS. ...................................................................................................... 20
Figura 14 Sezione BOUNDARIES. ................................................................................................... 22
Figura 15 Grafico tipo VECTOR. ...................................................................................................... 25
Figura 16 Grafico tipo CONTOUR di vz. .......................................................................................... 26
Figura 17 Effetto geometria sulla scia a valle del disco. ................................................................... 27
Figura 18 ELEVATION a monte del disco. ...................................................................................... 28
Figura 19 ELEVATION in corrispondenza del disco. ....................................................................... 28
Figura 20 ELEVATION a valle del disco. ......................................................................................... 29
Figura 21 ELEVATION lungo la coordinata assiale. ........................................................................ 30
Figura 22 CONTOUR vr. .................................................................................................................. 31
Figura 23 CONTOUR pressione. ....................................................................................................... 32
Figura 24 Coefficiente d’attrito-Re. ................................................................................................... 34
Figura 25 Coefficiente d’attrito-Re: Analisi di letteratura comparata ai risultati ottenuti con
FlexPDE. ............................................................................................................................................ 35
Figura 26 Coefficiente d’attrito-Re per raggi differenti..................................................................... 37
Figura 27 h e Nu nella sezione DEFINITIONS ................................................................................. 39
Figura 28 Flussi termici nella sezione DEFINITIONS. ..................................................................... 39
Figura 29 Equazione dell’energia nella sezione EQUATIONS......................................................... 40
Figura 30 Grafico CONTOUR temperatura....................................................................................... 40
Figura 31 ELEVATION temperatura sezione ingresso dominio. ...................................................... 41
Figura 32 ELEVATION T lungo la coordinata assiale...................................................................... 41
Figura 33 ELEVATION T lungo l’asse z. ......................................................................................... 42
Figura 34 CONTOUR qz. ................................................................................................................... 43
Figura 35 CONTOUR qr. ................................................................................................................... 44
Figura 36 Grafico VECTOR del flusso di calore. .............................................................................. 45
Figura 37 Q-Re................................................................................................................................... 46
Figura 38 Effetto Re sul profilo di T.................................................................................................. 46
Figura 39 Correlazione disco investito. ............................................................................................. 47
Figura 40 h-raggio .............................................................................................................................. 49
I
Indice tabelle
Tabella 1 Moto intorno ad un disco: caso standard. .......................................................................... 24
Tabella 2 Misure del coefficiente d’attrito e forza al variare di Re. .................................................. 34
Tabella 3 Misure con raggio = 1.5 volte raggio standard. ................................................................. 36
Tabella 4 Misure con raggio = 2 volte raggio standard. .................................................................... 36
Tabella 5 Coefficiente d’attrito al variare del raggio a parità di velocità v0. ..................................... 37
Tabella 6 Studi al variare del numero di Reynolds. ........................................................................... 45
Tabella 7 Effetto raggio del disco ...................................................................................................... 48
II
Introduzione
Numerosi fenomeni fisici, chimici e biologici, sono caratterizzati dal trasporto di specifiche
grandezze. La disciplina dei fenomeni di trasporto ha identificato i meccanismi di trasporto di
quantità fisiche quali:



Quantità di moto
Calore
Materia
cercando di unificare l’approccio fisico e matematico dei tre trasporti.
Le leggi di bilancio per queste grandezze conducono ad equazioni differenziali alle derivate parziali
(PDE), che risultano spesso non risolvibili analiticamente.
Con lo sviluppo tecnologico degli ultimi decenni, sono nate nuove discipline che sfruttano software
per la risoluzione delle PDE, vedasi in ambito del trasporto di quantità di moto la fluidodinamica
computazionale.
Uno dei programmi utilizzati a questi scopi è FlexPDE, sviluppato e venduto dalla PDE solution
Inc., che consente la risoluzione di equazioni differenziali alle derivate parziali di ordine superiore
al primo, sfruttando il metodo numerico alle differenze finite.
Il moto attorno ad oggetti sommersi è un processo di studio ingegneristico, che data la natura delle
equazioni che lo descrivono può essere ottimizzato con l’uso di software. Di particolare interesse è
la valutazione di:


Coefficiente d’attrito (f)
Numero di Nusselt (Nu)
Lo studio di questi due parametri, al variare del numero di Reynolds e della geometria del corpo,
consente di caratterizzare il sistema in termini fluidodinamici e termici. Infatti entrambi i parametri
sono legati alla valutazione del profilo di velocità e temperatura, nonché ai flussi di quantità di moto
ed energia.
Lo scopo di questo progetto è lo studio di un disco investito da un fluido, a partire da una analisi
teorica del sistema in esame, caratterizzando il sistema nell’ambiente di calcolo del programma
FlexPDE.
1
Capitolo 1
Teoria del moto intorno ad oggetti
1.1.
Generalità sul moto dei fluidi
Lo studio del moto dei fluidi è uno dei problemi più frequenti e rilevanti che si possono incontrare
nell’analisi dei fenomeni di trasporto. Possibili esempi sono: il moto di un fluido attorno a tubi nei
quali scorre un altro fluido per ottenere un certo scambio termico; il flusso di una miscela reagente
attraverso un letto fisso o fluidizzato di particelle quando si ha a che fare con reazioni chimiche
eterogenee; il moto di un fluido in torri di assorbimento riempite con opportuno materiale di
riempimento, come gli anelli o i corpi a sella.
In tutte queste situazioni vi sono importanti quantità da valutare per via sperimentale, sia in sede di
progetto che di conduzione dell’unità di interesse. In particolare, in sede di progetto è fondamentale
avere dati sperimentali in quanto una progettazione fatta puramente su una base teorica non potrà
mai realizzare un’unità operativa funzionante e funzionale.
Tra le quantità di maggior interesse nei casi di moto intorno ad oggetti sommersi vi sono la forza di
trascinamento, il calore e la materia scambiata. Queste quantità possono essere valutate solo tramite
la conoscenza dei coefficienti di scambio rispettivi, cioè il coefficiente d’attrito (f), il coefficiente di
scambio termico (h) e quello di scambio di materia (Kc). Un importante annotazione riguarda il
fatto che f e h non sono coefficienti di scambio in senso stretto perché lo sono più propriamente le
quantità
e h/
. Infatti, moltiplicando questi coefficienti (che racchiudono in loro tutta
l’ignoranza sul meccanismo di trasporto) per le rispettive forza motrici si ottengono le quantità
macroscopiche desiderate.
I coefficienti di scambio possono essere valutati con correlazioni empiriche o semi-empiriche in cui
compaiono dei numeri adimensionali, che contengono alcune grandezze caratteristiche del sistema,
e proprietà fisiche, come la densità, la viscosità, il calore specifico, la conducibilità termica e la
diffusività di materia (per condizioni di laminarità si possono anche usare equazioni derivate per via
teorica). Per il coefficiente di attrito è necessario conoscere il numero di Reynolds (Re); il
coefficiente di scambio termico viene valutato attraverso il numero di Nusselt (Nu), per il quale
bisogna conoscere sia il numero di Reynolds che quello di Prandtl (Pr), oltre alla conducibilità
termica e una lunghezza caratteristica; il coefficiente di scambio di materia si trova tramite il
numero di Sherwood (NuAB) e per questo necessita del numero di Reynolds e di quello di Schmidt
(Sc), oltre alla diffusività di materia e una lunghezza caratteristica. Le correlazioni suddette valgono
sia per moto intorno ad oggetti che per quello intubato, con l’unica differenza, però, che in
quest’ultimo, in regime laminare, si deve tener conto nella correlazione anche della geometria
(rapporto diametro/lunghezza). Gran parte delle correlazioni viene presentata in forma grafica, dal
momento che una sola correlazione difficilmente riesce a descrivere la relazione coefficiente di
scambio-numeri adimensionali per tutto il campo di interesse.
I numeri adimensionali sono definiti come segue:
2
Dove μ,
k, DAB, Cp, , vc, L, sono, rispettivamente, la viscosità, la densità, la conducibilità
termica, la diffusività di materia, il calore specifico del fluido, lo sforzo alla parete, la velocità e la
lunghezza caratteristica dello specifico sistema.
Tutti questi numeri adimensionali hanno un preciso significato fisico di cui, per brevità, non
discorreremo.
In letteratura vi sono molte correlazioni per diversi numeri di Reynolds, cioè per diversi regimi di
moto.
I fenomeni di trasporto della quantità di moto, energia e materia sono fondati su meccanismi fisici
molto analoghi, soprattutto gli ultimi due e questo comporta che le correlazioni per ciascun
fenomeno sono molto simili. Inoltre sono state trovate delle analogie tra i numeri adimensionali
caratteristici dei tre fenomeni che mettono in relazione i coefficienti di scambio tra loro.
A volte si può aver la necessità di dover conoscere anche i profili delle variabili principali, come la
velocità, lo sforzo, la pressione, la temperatura; addirittura, per scopi di ottimizzazione si può anche
aver bisogno di conoscere i coefficienti di scambio puntualmente, e non solo quelli medi.
L’unico modo per avere un controllo completo di tutte le variabili e per calcolare le quantità di
interesse, sia globalmente che localmente, è di usare un calcolatore impostando in modo corretto le
equazioni di bilancio, che sono equazioni differenziali alle derivate parziali non lineari di grado
superiore al primo. Purtroppo, non sempre il calcolatore riesce a trovare una soluzione delle
equazioni, anche se il problema è ben posto. Questa eventualità ha una probabilità di accadimento
che aumenta all’aumentare del numero di Re, cioè all’aumentare della turbolenza. Infatti il regime
turbolento completamente sviluppato non può essere spiegato in modo esaustivo, nemmeno per via
semi-empirica. Di seguito riportiamo le suddette equazioni nella loro forma generale (è stata omessa
l’equazione per la materia):
Equazione di continuità (scalare)
Equazione del moto (vettoriale)
Equazione dell’energia (scalare)
3
dove a primo membro di tutte le equazioni è presente l’operatore di derivata sostanziale, il quale
racchiude il termine di accumulo e quelli di trasporto convettivo della variabile di interesse mentre a
secondo membro ci sono i termini di trasporto molecolare e le forze motrici.
Purtroppo, anche nel caso in cui esista una soluzione esatta, gli output di un calcolatore che risolve
tali equazioni complesse hanno comunque un certo margine di errore.
1.2.
Moto intorno ad oggetti
Tra i vari sistemi fluidodinamici esistenti, una particolare importanza assumono quelli in cui si ha
moto intorno ad oggetti. Questo tipo di moto è molto difficile da studiare a causa della forma non
rettilinea delle linee di flusso e, quindi, non si riesce a descrivere analiticamente gli andamenti delle
velocità e della pressioni. Nonostante queste difficoltà, per valori del numero di Reynolds minori di
0.1 si è riusciti a trovare espressioni analitiche per le velocità, la pressione e gli sforzi per il moto
intorno ad una sfera: il problema di Stokes. A parte questo esempio, comunque, per studiare i profili
delle variabili fluidodinamiche principali attorno ad un oggetto investito bisogna ricorrere ad
esperimenti.
Il corpo che viene investito subisce una forza, detta di trascinamento o anche resistenza del mezzo,
che è dovuta sia all’impatto del fluido contro di esso, sia a causa di un gradiente di pressione tra
monte e valle dello stesso. La resistenza del mezzo interviene in molte applicazioni dell’ingegneria
chimica; vi sono strumenti di regolazione che possono rendere massima o minima la resistenza del
mezzo, secondo le necessità; altri esempi di resistenza del mezzo sono la caduta di una particella
solida in un fluido e le sollecitazioni impresse alle strutture per l’azione del vento o dell’acqua
corrente.
Parametro importante per valutare la forza di trascinamento è il coefficiente di attrito, la cui
espressione è:
dove A è la superficie d’impatto, è la densità del fluido,
è la velocità di libero flusso, detta
anche di avvicinamento, ed F è la forza di trascinamento. Il numero di Reynolds, invece, ha la
seguente espressione
dove D è la lunghezza caratteristica dell’oggetto, data da volume/superficie, e è la viscosità del
fluido.
Come già detto f è funzione di Reynolds. La sua correlazione empirica viene presentata
graficamente con un diagramma logaritmico mostrato in figura 1 (specifico per sfera, disco e
cilindro) [1]. Come si può ben vedere dalla figura, per valori di Re minori di circa 0.1 gli sforzi
viscosi predominano e si hanno andamenti lineari (in un diagramma logaritmico); in questo
intervallo il coefficiente di attrito ha espressione:
4
Dove C è una costante diversa a seconda dell’oggetto.
Per valori di Reynolds tra 0.1 e 105 si ha una zona di transizione in cui l’andamento di f non è
lineare.
Per valori di Reynolds superiori a circa 105 si ha la prevalenza degli sforzi inerziali e una zona dove
f è più o meno costante; ad un certo punto però, si ha un brusco abbassamento della curva e di
nuovo un andamento all’incirca costante. Quest’ultimo fenomeno è spiegato dal fatto che ad un
certo punto lo strato limite, causato dall’attrito del fluido sull’oggetto, diminuisce bruscamente a
causa del distacco di una sua parte, facendo così diminuire l’ingombro dovuto all’oggetto. Il
distacco dello strato limite avviene ad un determinato Reynolds per una determinata geometria; per
il disco il Reynolds critico è all’incirca 200000.
Figura 1 Diagramma Cd-Re.
Per quanto riguarda il trasporto del calore parametro importante, come già detto, è il coefficiente di
scambio termico h, valutabile tramite il numero di Nusselt; la lunghezza caratteristica da assumere
per valutare questo numero è la stessa usata per il coefficiente di attrito. Anche per il coefficiente di
scambio h vi sono in letteratura molte correlazioni, soprattutto per il moto intorno a sfere e cilindri,
ciascuna valida in un range di valori di Reynolds definito. In genere queste correlazioni vengono
presentate in forma grafica, tramite la funzionalità tra il numero di Nusselt e quelli di Reynolds e di
Prandtl o tra il fattore di Colburn JH e il numero di Reynolds. La definizione di JH è la seguente:
5
La figura 2 [2] mostra una correlazione per il caso della lastra piana investita. In questo diagramma,
effettivamente, viene rappresentato l’andamento del fattore di Colburn contro Reynolds.
Figura 2 Correlazione JH-Re per una lastra investita.
1.3.
Case study: “moto intorno ad un disco”
Nel nostro studio, cioè il “moto attorno ad un disco” (fig. 3), abbiamo studiato gli effetti relativi ai
fenomeni di trasporto della quantità di moto e dell’energia.
Più dettagliatamente il sistema studiato è costituito da un disco che viene impattato assialmente da
una corrente d’aria. Abbiamo assunto che il nostro fluido abbia proprietà costanti ( , μ e ) e che il
sistema sia in condizioni di stazionarietà.
Figura 3 Moto intorno ad un disco
Le equazioni utilizzate sono state l’equazione di continuità, del moto e dell’energia in forma
differenziale in coordinate cilindriche:
6
Equazione di continuità
L’equazione di continuità è l’espressione matematica del bilancio microscopico di materia; è
costituita dal termine di accumulo
, e dai termini che esprimono la velocità di ingresso e di uscita
della massa. Il termine di generazione è nullo, in accordo alla legge di Lavoisier.
Equazione del moto, componente r
Equazione del moto, componente z
L’equazione del moto si origina dal bilancio microscopico della quantità di moto ed ha natura
vettoriale. La componente
dell’equazione non è stata riportata perché trascurabile nel nostro
sistema. In queste equazioni abbiamo al primo membro l’accumulo e i termini legati al trasporto
convettivo di quantità di moto; a secondo membro si ha la generazione della quantità di moto per
gradiente di pressione e per gravità e il suo trasporto molecolare (termine moltiplicante la μ).
Equazione dell’ energia
L’equazione dell’energia ha la stessa struttura dell’equazione del moto; al primo membro accumulo
e trasporto convettivo, mentre, al secondo, generazione per dissipazione viscosa (termine vicino a
μ) e trasporto molecolare (termine vicino a k).
Questa struttura simile delle equazioni testimonia l’analogia dei fenomeni, la quale, comunque,
viene a mancare un po’ a causa dei diversi tipi di generazione.
Queste equazioni, tuttavia, sono state opportunamente semplificate per il caso di studio in
questione.
La prima semplificazione riguarda la stazionarietà del processo che ha permesso di eliminare i
termini di accumulo:
,
,
e
. Inoltre, considerando la direzione assiale del moto
7
possiamo trascurare la componente
della velocità e anche la dipendenza di tutte le altre variabili
da questa coordinata; i termini eliminati sono:
dall’equazione del moto lungo r;
,
z;
e
dall’equazione di continuità;
e
e
dall’equazione del moto lungo
dall’equazione dell’energia. In quest’ultima è stata anche trascurata la generazione
per dissipazione viscosa.
Equazione di continuità semplificata
Equazione del moto lungo z semplificata (assiale)
Equazione del moto lungo r semplificata (radiale)
Equazione dell’energia semplificata
Come detto la risoluzione al calcolatore di queste complesse equazioni può fornire i profili della
pressione (P), delle componenti della velocità lungo r e z (vr e vz) e della T; a partire da questi è poi
possibile ottenere anche i profili degli sforzi e dei flussi di energia.
In coordinate cilindriche, con le nostre condizioni, le espressioni degli sforzi viscosi e dei flussi
termici (leggi generalizzate di Newton della viscosità e di Fourier della conducibilità) sono:
Tramite l’integrazione sulla superficie del disco dei profili degli sforzi e dei flussi di energia è
possibile arrivare alla forza di trascinamento e alla portata termica. In particolare la forza di
trascinamento totale è data dalla somma di 3 termini che sono: integrale dello sforzo
sulle due
superfici circolari, l’integrale dello sforzo
sulla superficie laterale e l’integrale della pressione
sulle due superfici circolari. La portata termica totale, invece, è data dalla somma dei seguenti 2
termini: integrale di sulle due superfici circolari e l’integrale di sulla superficie laterale.
8
Infine, una volta ottenuta la forza di trascinamento e la portata di calore complessive, possiamo
calcolare sia il fattore di attrito f che il coefficiente di scambio termico h con le equazioni:
dove
è la differenza tra la temperatura del disco e quella ambiente e
disco. Trovato h è immediato il calcolo del numero di Nusselt.
9
è l’area totale del
Capitolo 2
Uso di FlexPDE per la risoluzione del moto attorno ad un disco
2.1.
Introduzione
I fenomeni di trasporto sono descritti da equazioni differenziali alle derivare parziali. La funzione è
cosi rappresentata indirettamente attraverso una relazione fra se stessa e le sue derivate parziali.
Non esistono metodi generali analitici per risolvere tali equazioni, tuttavia, sono adottate delle
semplificazioni per specifici casi (campo velocità irrotazionale, strato limite, creeping flow,
problema del Graetz, etc.).
L’evoluzione informatica ha permesso lo sviluppo di software per lo studio dei problemi dei
fenomeni di trasporto, potendo cosi andare oltre al solo studio di problemi in flusso laminare e di
semplici geometrie.
L’approccio tipico di tali software è di discretizzare il sistema da analizzare in celle elementari, per
ottenere una griglia di calcolo, detta mesh. In questo modo il dominio iniziale è suddiviso in
sottodomini più piccoli e di geometria più semplice (esaedri e tetraedri in 3D, quadrilateri e
triangoli in 2D ecc.), così da risolvere le equazioni all’interno di ciascun sottodominio e di
raccordare le soluzioni ai bordi in maniera tale da rispettare le condizioni iniziali fissate sul bordo
del dominio originale.
La mesh (fig. 4) può presentarsi in forma:



Strutturata: costituito da elementi finiti di tipo uguale, lo svantaggio principale sta nella
scarsa adattabilità della griglia a geometrie complesse ma tale mesh consente la possibilità
di sviluppare algoritmi più efficienti
Strutturata a blocchi: costituito da elementi finiti misti (uguali e differenti), consente di
superare i limiti geometrici della mesh strutturata. Gli svantaggi sono legati ad una difficile
gestione dei sottodomini e nelle loro interconnessioni
Non strutturata: costituito da elementi finiti differenti, presenta un elevata flessibilità
geometrica
Figura 4 Tipi di mesh: strutturata-strutturata a blocchi- non strutturata.
10
Esistono 3 tecniche di discretizzazione:



Metodo ai volumi finiti: le equazioni vengono risolte in un volume di controllo discreto e
sono integrate definendo le condizioni al contorno sui confini del dominio
Metodo agli elementi finiti: Si tratta di una tecnica numerica che ha lo scopo di cercare
soluzioni approssimate attraverso la risoluzione di equazioni differenziali alle derivate
parziali le quali vengono ridotte ad un sistema di semplici equazioni algebriche
Metodo alle differenze finite: tale metodo si basa sull’approssimazione diretta delle
equazioni differenziali alle derivate parziali ottenute sostituendo alle derivate parziali delle
differenze definite sul dominio del problema
Per ogni singolo metodo, la procedura di analisi risulta simile:
1. Viene definita la geometria del problema da analizzare
2. Il sistema in esame viene discretizzato
3. Viene definito il modello fisico-matematico (equazioni del moto, equazione dell'energia)
e quindi quello numerico (metodo di discretizzazione delle equazioni, algoritmi per la
risoluzione delle equazioni)
4. Vengono definite le condizioni al contorno e quelle iniziali
5. Vengono risolte le equazioni in maniera iterativa. Il calcolo viene interrotto una volta
che sia stato raggiunto il grado di accuratezza desiderato.
2.2.
FlexPDE
Per lo sviluppo di questo progetto, è stato utilizzato il software FlexPDE 6.20, sviluppato e venduto
dalla PDE solution Inc. Il programma esegue le operazioni necessarie alla risoluzione delle PDE
con il metodo agli elementi finiti, presentando risultati grafici e numerici [3].
Uno script di descrizione del problema è un file di testo leggibile, nel quale l’utilizzatore specifica
l’equazione differenziale alle derivate parziali da risolvere, la geometria del problema, le condizioni
al contorno e la visualizzazione degli output. Il programma legge tale file, costruisce la mesh e
risolve l’equazione, presentando gli output che sono stati specificati dall’utente. FlexPDE agisce
come una “scatola nera” (black box), infatti l’utente non deve fare i conti con i dettagli degli
algoritmi del programma.
I contenuti del file sono costituiti da un numero di sezioni, ciascuna identificata da un nome. Le
sezioni fondamentali sono (fig.5):
 TITLE: nome file
 COORDINATES: sistema di coordinate adottato
 SELECT: selezionare accuratezza nei calcoli
 VARIABLES: variabili del problema
 DEFINITIONS: parametri del problema (densità, viscosità, diffusività, etc.)
 EQUATIONS: equazioni che descrivono il sistema in esame
 BOUNDARIES: dominio problema
 MONITORS e PLOTS: risultati grafici
 END
11
Figura 5 Sezioni FlexPDE.
Lo scopo del progetto è stato di analizzare il moto attorno ad un disco di piccolo spessore, in
termini fluidodinamici e termici.
Le operazioni compiute sono state cronologicamente:
1.
2.
3.
4.
5.
Costruzione dominio del problema
Individuazione variabili del problema
Definizione principali parametri del problema
Introduzione equazioni che descrivono i fenomeni di trasporto
Scelta delle condizioni al contorno
6. Visualizzazione dei risultati ottenuti
12
2.3.
Costruzione dominio del problema
2.3.1. Il sistema di coordinate
La prima operazione compiuta è stata la scelta del sistema di coordinate. Poiché l’oggetto investito
da un fluido è un disco, è stato utilizzato un sistema di coordinate cilindriche. Il disco è stato
ottenuto dalla rotazione di un rettangolo attorno ad un asse, detto asse di rotazione (l’asse z, ad
esempio, rappresenta l’asse di rotazione della regione piana tratteggiata mostrata in figura 6).
Figura 6 Disco come sodido di rotazione.
Tale funzione è stata impostata nel programma attraverso il comando XCylinder ('z','r') nella
sezione COORDINATES (fig. 7), che consente di definire l’asse delle ascisse (individuato dal
programma con X) come asse di rotazione ed è stato scelto di attribuire a tale asse la lettera z e
all’asse ortogonale la lettera r.
Figura 7 Sezione COORDINATES.
2.3.2. Il dominio
Il dominio del problema è stato definito nella sezione BOUNDARIES.
Il progetto ha previsto la descrizione di un disco investito da un fluido in un ambiente aperto; però
nessun software è in grado di simulare un tale sistema (moto attorno ad oggetti).
Data l’impossibilità di descrivere un dominio aperto con il programma, è stata definita una regione
esterna, come se il disco fosse contenuto in un condotto esterno di dimensioni tali da simulare un
sistema aperto. Tale dominio esterno è stato caratterizzato dai parametri introdotti nella sezione
DEFINITIONS:


raggio_condotto: rappresenta il raggio dell’ipotetico dominio esterno di forma cilindrica; è
stato fissato con una dimensione di 0.5 m, pari a 10 volte il raggio del disco così da
simulare un sistema aperto.
L: rappresenta la lunghezza del dominio esterno; è stata fissata a 0.75 m tale che il flusso di
velocità potesse completamente svilupparsi dopo l’impatto con il disco.
13
Il disco è stato caratterizzato dai parametri geometrici:


raggio: individua il raggio del disco, è stato fissato a 0.05 m
Lz: il suo doppio rappresenta lo spessore del disco; consente di individuare la posizione
rispetto all’origine degli assi pari a (0.20-Lz) m
Il disco è stato tracciato nel dominio del problema, seguendo il seguente percorso (come si vede in
figura):
(0,0) → (z-Lz, 0) → (z-Lz, raggio) → (z+Lz, raggio) → (z+Lz, 0) → (L, 0)
Figura 8 Dominio del problema.
Come è possibile notare FlexPDE lavora con una mesh di tipo non strutturato. Il programma
focalizza la mesh nelle regioni dove sono presenti elevati gradienti, consentendo così un uso
efficace della griglia.
14
Figura 9 Mesh disco investito.
2.4.
Individuazione variabili del problema
Le variabili del sistema studiato sono il vettore velocità, costituito dalla componente radiale vr e da
quella parallela al flusso vz, la pressione p e la temperatura T. Lo studio del processo è stato
condotto allo stato stazionario; in questo modo, il tempo non è stato inserito nelle variabili e
soprattutto non è stato necessario inserire condizioni iniziali per la risoluzione delle equazioni.
Figura 10 Sezione VARIABLES.
15
2.5.
Definizione principali parametri del problema
I principali parametri del sistema sono stati introdotti nella sezione DEFINITIONS come mostrato
dalla figura seguente.
Figura 11 Sezione DEFINITIONS.
Nel sistema come fluido-test è stata impiegata aria; nella sezione DEFINITIONS ne sono state
introdotte le principali proprietà fisiche (densità, viscosità, conducibilità termica, calore specifico).
Tali parametri sono stati valutati alla temperatura di film del sistema, ed è stato fatto uso di unità di
misura del SI (sistema internazionale).
16
I parametri fisici dell’aria sono stati definiti con la seguente nomenclatura:
→ densità
→ viscosità
→ calore specifico
→ conducibilità
In particolare il prodotto tra la densità (dens) ed il calore specifico (Cp) è stato individuato col
parametro rCp, per così facilitare la scrittura dell’equazione dell’energia.
Le grandezze fluidodinamiche introdotte in questa sezione sono state:




La velocità di libero flusso (v0)
Gli sforzi viscosi ( )
La forza di drag (Fz)
Il coefficiente d’attrito (f)
La velocità di libero flusso o velocità di avvicinamento (free stream velocity, v0) è stata impostata a
partire dalla definizione del numero di Reynolds. È stato necessario lavorare con un Re molto basso
e pari a 0.0001. A questo valore il sistema è rimasto in regime laminare, così da consentire al
programma una corretta risoluzione delle equazioni.
Con Re che rappresenta il numero di Reynolds ed il doppio del raggio che individua la lunghezza
caratteristica del caso di studio (il diametro del disco).
Sulla variabile velocità (v) è stata definita la divergenza (div_v), così da poter essere utilizzata
nell’equazione del moto nella sezione EQUATIONS:
Per poter ottenere output grafici specifici, ad esempio per visualizzare il campo vettoriale della
velocità), è stato impiegato il comando vector. Con questo comando, la variabile v stata identificata
come un vettore costituito da componenti vz e vr.; anche il modulo del vettore v è stato introdotto
(speed):
Per la valutazione della forza agente sul disco in direzione parallela alla corrente di fluido (Fz, forza
di drag), è stato necessario definire gli sforzi viscosi impiegando l’equazione di Newton della
viscosità, opportunamente semplificata per il caso in questione:
17
;
;
;
Con il comando surf_integral (si veda la figura 12) è stato possibile calcolare l’integrale degli sforzi
viscosi e della pressione sulle superfici su cui agiscono. Tale comando consente di valutare
l’integrale superficiale su di una superficie specificata e descritta nella sezione BOUNDARIES con
il comando Feature (figura sottostante).
Figura 12 Feature.
Si riporta la nomenclatura impiegata per la valutazione di Fz:
con:
e:
, forza dovuta allo sforzo
agente sulla superficie S1
, forza dovuta allo sforzo
agente sulla superficie S2
18
, forza dovuta alla pressione agente sulla superficie S1
, forza dovuta alla pressione agente sulla superficie S2
, forza dovuta allo sforzo
agente sulla superficie S3
La forza Fz è stata utilizzata nella definizione del coefficiente d’attrito f:
con:
, superficie ottenuta proiettando il disco su un piano perpendicolare alla
velocità di avvicinamento del fluido.
Le grandezze legate allo scambio termico introdotte nello script sono state:




Temperatura disco e ambiente
Portate termiche
Coefficiente di scambio termico (h)
Numero di Nusselt (Nu)
La temperatura del disco (Tdisco) e dell’ambiente circostante (T0) sono state fissate rispettivamente a
400
e 300 .
Utilizzando il comando surf_integral sono state valutate le portate termiche attraverso le superfici
S1, S2, S3. Tali grandezze sono state impiegate per il calcolo del coefficiente di scambio h e del
numero di Nusselt (Nu):
dove:
;
19
, portata termica attraverso S1 dovuta a
, portata termica attraverso S2 dovuta a
, portata termica attraverso S3 dovuta a
2.6.
Introduzione equazioni che descrivono i fenomeni di trasporto
Il software consente di inserire le equazioni che descrivono i fenomeni di trasporto, digitando ogni
termine come in un file di testo come si nota dalla figura sottostante.
Figura 13 Sezione EQUATIONS.
Il problema è descritto da 4 equazioni in 4 incognite (vr, vz,
, T):
→
→
→
I termini di trasporto convettivo nelle equazioni (
e
) rendono l’equazione di NavierStokes e dell’energia equazioni non lineari. Dato che le incognite eguagliano le equazioni il
problema è matematicamente ben posto, ma FlexPDE non è in grado di risolvere le equazioni
differenziali del primo ordine come è l’equazione di continuità. Tale limite è superato effettuando la
divergenza dell’equazione di Navier-Stokes; combinando questa con l’equazione di continuità si
arriva all’equazione di Poisson della pressione [4]:
A questa equazione è sottratto il termine
per la stabilità numerica. Tale termine aggiunge
una pressione addizionale tale da ridurre la divergenza della quantità di moto e così soddisfare la
legge di conservazione della materia. Otteniamo così:
Con c che deve essere un numero sufficientemente grande tale da far rispettare l’equazione di
conservazione della materia ed L che rappresenta una grandezza caratteristica del sistema.
20
Sperimentalmente è stato deciso di operare con un c pari a 10^-6, tale valore infatti ha consentito di
ottenere risultati attendibili.
2.7.
Scelta delle condizioni al contorno
Nella sezione BOUNDARIES (fig. 14), la definizione del dominio del problema è stata
accompagnata dalla scelta delle condizioni al contorno per le quattro variabili del sistema (vz, vr, p,
T).
FlexPde consente di utilizzare due comandi:


Value: definisce il valore che la variabile assume sul limite del domino
Natural: definisce il flusso sul limite del dominio
I limiti su cui sono state fissate le condizioni al contorno sono stati:





Ingresso dominio: sono stati fissati il valore di velocità di libero flusso (v0, con vz unica
velocità non nulla) e la temperatura dell’ambiente (T0)
Uscita del dominio: per ottenere un flusso completamente sviluppato si è adoperato il
comando Natural fissato a zero per tutte le variabili
Parete superiore dominio: la simulazione di un dominio aperto ha richiesto di fissare la
velocità vz al valore di libero flusso
Asse di simmetria: uso del comando Natural fissato a zero per la simmetria del sistema
Bordi disco: per soddisfare la condizione di no slip alle pareti del disco la velocità è stata
fissato a zero.
21
Figura 14 Sezione BOUNDARIES.
22
2.8.
Visualizzazione dei risultati ottenuti
FlexPDE consente di visualizzare i grafici relativi alla risoluzione delle equazioni
differenziali implementate sia durante i calcoli che a fine processo (una volta
raggiunto il grado di accuratezza desiderato) [3]:
 Monitors: soluzione in tempo reale (output intermedi)
 Plots: soluzione finale (output finale)
I grafici utilizzati nel progetto sono stati:
 Contour: contour variabile/parametro indicato
 Vector: campo vettoriale variabile specificata, che deve essere specificata
come vettore col comando vector nelle DEFINITIONS
 Elevation: andamento parametro/variabile indicato lungo una linea specificata
Particolarmente utile è stata la sezione aggiuntiva SUMMARY; ha consentito di
estrapolare degli output numerici mediante il comando report, così da poter inserire
tali valori in un foglio di calcolo e valutare gli andamenti delle variabili del sistema al
variare di specifici parametri.
23
Capitolo 3
Analisi fluidodinamica
3.1.
Introduzione
In questa sezione è stata analizzata la fluidodinamica del sistema, e sono stati discussi i risultati
ottenuti facendo rifermento a casi di studio presenti in letteratura.
In particolare sono stati valutati:



Il profilo di velocità: è stata effettuata l’analisi delle componenti di velocità
e . Gli
studi di letteratura indicano la necessità di lavorare a basso numero di Reynolds, che è stato
imposto <<1 (“moto di scorrimento”) così da evitare tutte le problematiche relative
all’instabilità risolutiva apportata da sistemi di tipo turbolento.
Forze agenti sul disco: a partire dall’integrazione degli sforzi sulle superfici di azione è stato
possibile calcolare il valore della forza di trascinamento sul disco
Il coefficiente d’attrito
Tali valutazioni sono state condotte oltre che per il caso standard (caratteristiche riassunte in
tabella), anche al variare del numero di Reynolds e del raggio del disco.
Re
Raggio disco, m
Spessore disco,m
v0 ,m/s
0.0001
0.05
0.01
2.06
Tabella 1 Moto intorno ad un disco: caso standard.
3.2.1. Analisi del profilo di velocità: caso standard
Il modo più intuitivo di visualizzare il profilo di velocità con FlexPDE è quello di
sfruttare il grafico vector, come mostrato in figura:
24
Figura 15 Grafico tipo VECTOR.
Questo tipo di rappresentazione consente di visualizzare in una maniera molto chiara la direzione
del fluido, ed il suo comportamento in prossimità del disco. Il vettore velocità è, infatti,
rappresentato come una freccia la cui punta ne specifica il verso.
Come è possibile notare dalla legenda associata alla figura, la velocità è compresa in uno specifico
range dipendente quasi esclusivamente dalla velocità lungo l’asse del disco, essendo vz e vr di ordini
di grandezza ben differenti:
0
< |v| < 2.15, scala
m/s
3.2.2. Analisi del profilo di velocità: componente vz
La componente rappresenta il contributo della velocità nella direzione parallela al movimento del
fluido.
Tale componente è quella più importante per lo studio del moto attorno al disco, essendo di due
ordini di grandezza superiore rispetto alla componente radiale (come sarà illustrato nei prossimi
paragrafi).
Il profilo di velocità è stato studiato utilizzando come grafici:
25


Contour
Elevation
Dall’immagine tipo “contour” sottostante è possibile constatare che il fluido in corrispondenza
dell’ostacolo rallenta, subendo un cambiamento di velocità progressivo (sulla destra dell’immagine
sono presenti i risultati numerici della simulazione e questi sono associati alla gradazione di colori),
ed essa diventa nulla sulle superfici del disco per la condizione di non scivolamento.
Superato l’oggetto, il fluido ritorna allo stesso valore della velocità di avvicinamento ma più
gradualmente essendo il fluido reale. Infatti un fluido a viscosità nulla avrebbe aggirato l’ostacolo
restituendo esattamente lo stesso profilo anche dopo il superamento dell’oggetto. La condizione di
fluido reale è anche la causa del riscontro di velocità negative a valle del disco; tali valori mostrano,
dunque, un ricircolo del fluido.
Figura 16 Grafico tipo CONTOUR di vz.
Dal grafico è possibile notare un aumento di velocità soprattutto nella zona superiore al disco. Tale
zona è individuata dal programma con una “x”, ed è dovuta al sommarsi del contributo di velocità
imposto alla parete superiore del dominio esterno tracciato, ed al profilo di velocità caratteristico
per la presenza del disco. È possibile costatare che lo sviluppo del profilo di velocità è dipendente
dalle caratteristiche geometriche del disco, come mostrato dal confronto dei tre grafici seguenti a
raggio del disco differente (rispettivamente 0.05, 0.075 e 0.1 m).
26
Figura 17 Effetto geometria sulla scia a valle del disco.
3.2.3. Analisi del profilo di vz su diverse sezioni
Il profilo delle vz del fluido può essere analizzato in diverse sezioni del dominio, sia trasversali al
suo movimento, che parallele. Sono state così monitorate le zone che possono destare maggiore
interesse.
Si è voluto evidenziare nei grafici seguenti, come i valori di v z tendono a raggiungere in maniera
diversa il valore v0 imposto alla parete superiore del dominio esterno, a seconda della posizione
lungo sull’asse z. In particolare sono stati tracciati dei grafici di tipo elevation che mostrano
l’andamento di vz al variare della coordinata radiale.
Sono stati monitorati lungo z:



Elevation pre-disco
Elevation con z in prossimità del disco
Elevation dopo il disco
Prima del disco è possibile notare che in prossimità dell’asse la velocità è vicina ai 3*10^-10 m/s, e
quindi raggiunge un massimo intorno a R=0,2m per poi portarsi al valore di velocita .
27
Figura 18 ELEVATION a monte del disco.
In corrispondenza della coordinata assiale dove è presente il disco, per effetto della restrizione della
sezione di passaggio vi è un forte aumento della velocità. Il valore è 0 alla parete e raggiunge
“velocemente” il valore massimo in corrispondenza di un terzo dell’altezza di passaggio per poi
ritornare a quello di avvicinamento.
Figura 19 ELEVATION in corrispondenza del disco.
28
Superato il disco, ad una certa distanza dall’ostacolo il fluido risente ancora della sua presenza. Il
valore di massimo è più marcato rispetto alla sezione iniziale ed è ad un R di circa 0.35m.
Figura 20 ELEVATION a valle del disco.
Infine, in direzione parallela al fluido può essere osservato un comportamento della componente z
della velocità “ad ali di gabbiano”. Il valore di plateau è praticamente lo stesso. A contatto con
l’oggetto raggiunge il valore nullo che è stato imposto.
29
Figura 21 ELEVATION lungo la coordinata assiale.
3.2.4. Analisi relativa a vr
La componente radiale di velocità presenta due ordini di grandezza inferiori rispetto alla
componente assiale:
-3.14 < vr < 3.77, scala
m/s
Tale componente è massima in prossimità delle estremità del disco, come mostrato in figura 22;
come ci si aspetta a sinistra è positiva (verso l’alto) e a destra è negativa (verso il basso). I valori
minori di zero sono giustificati dalla scelta del sistema di riferimento.
Superato l’ostacolo il valore ritorna gradualmente nullo.
30
Figura 22 CONTOUR vr.
3.3.
Analisi della pressione
La pressione, attenendosi alla legge di conservazione dell’energia meccanica, subisce un
innalzamento in corrispondenza dei punti in cui la velocità diminuisce. Questo accade in prossimità
del disco, dove è massima; poi, man mano che ci si allontana da questo si porta fino a valore nullo.
La pressione risulta essere quasi simmetrica rispetto al disco (in valore assoluto) e non precisamente
simmetrica, essendo il fluido reale (fig. 23).
31
Figura 23 CONTOUR pressione.
3.4.
Introduzione allo studio delle tensioni e del fattore di attrito
Il quadro completo degli sforzi presenti all’interno di un fluido è definito dal tensore degli sforzi.
Nel caso di un sistema in coordinate cilindriche le componenti diagonali
,
rappresentano
gli sforzi normali, nel senso che ciascuna di esse fornisce la componente ortogonale della forza che
agisce sulla superficie piana elementare normale ad ogni asse coordinato. Le sei componenti non
diagonali
,
,
sono invece gli sforzi di taglio tangenti alle superfici elementari, e
si manifestano quando strati paralleli adiacenti di materia scorrono uno rispetto all'altro. Tenendo
poi in considerazione che gli elementi relativi agli sforzi tangenziali al di fuori della diagonale sono
uguali per simmetria, per completare la perfetta descrizione del tensore sono necessari 6
componenti complessive. Le componenti τ vengono sfruttate per il calcolo delle componenti della
forza ricevuta dal disco per effetto del suo movimento. Le componenti della forza si ottengono
dall’integrazione degli sforzi sulla superficie di azione.
Il tensore degli sforzi in questo caso ha la seguente struttura:
T=
32
Ma diverse componenti sono nulle, per cui non hanno effetto sulla forza di trascinamento.
e
sono entrambe nulle per simmetria. Anche le componenti
e
,
e
danno effetto nullo
sul disco. Il tensore finale, quindi risulta:
T=
La forza lungo la direzione del moto agente sul disco è data dalla somma delle forze di pressione,
più quelle dovute alle τ.
Per il moto intorno ad oggetti sommersi, tale forza è proporzionale ad A, che è l’area che deriva
dalla proiezione del solido in un piano perpendicolare alla velocità del fluido in avvicinamento. Per
un disco è l’area del cerchio A=
.
Nel nostro caso sono noti (calcolata integrando gli sforzi e le pressioni sulle superfici su cui
agiscono) e
(parametro fissato), quindi è possibile ricavare il fattore di attrito . Il problema è
stato oggetto di studio già in alcuni testi [5].
Lo studio della forza di trascinamento agente sul disco e del coefficiente d’attrito è stato effettuato
variando i seguenti parametri:


Numero di Reynolds
Raggio disco
Attraverso tali variazioni è stato potuto generalizzare il caso di studio.
3.5.
Caratterizzazione del coefficiente di attrito per variazione di Re
Un primo studio del coefficiente d’attrito ha riguardato la sua variazione con Re.
Fissati i parametri geometrici del sistema (caso standard), è stato fatto variare il numero di
Reynolds, che per come è stato definito lo script ha implicato anche una variazione della velocità di
libero flusso v0.
I dati ottenuti dalla simulazione di f (tabella 2), presi nella schermata finale del programma nella
SUMMARY, sono stati diagrammati in un grafico log-log contro il numero di Reynolds, si può
notare un andamento praticamente lineare (fig. 24). I dati sono sovrapponibili al diagramma
estrapolato dal Perry’s Chemical Engineers Handbook (fig. 25). Di buono c’è comunque la
possibilità di esplorare il comportamento nel moto di Stokes con una precisione elevatissima.
Sebbene questo tipo di simulazioni siano teoricamente attendibili soltanto nel regime di Stokes,
limitato a Re<0.1, è stato possibile estendere i risultati ottenuti anche a parte della regione
intermedia. Il coefficiente di attrito fino a Re≈10 mostra ancora proporzionalità inversa; invece, per
valori superiori a tale soglia FlexPDE non riesce a portare a convergenza i calcoli essendo quel
regime di flusso non simulabile.
33
Re
Forza [N]
Superficie di Impatto
[m2]
Coefficiente
di Attrito f
0.00001
4.83765*10^-14
7.853982*10^-3
2874481
0.0001
4.8371*10^-13
7.853982*10^-3
287415.2
0.001
4.83704*10^-12
7.853982*10^-3
28741.2
0.01
4.83718*10^-11
7.853982*10^-3
2874.2
0.1
4.83853*10^-10
7.853982*10^-3
287.5
1
4.93109*10^-9
7.853982*10^-3
29.3
10
7.40505*10^-8
7.853982*10^-3
4.4
Tabella 2 Misure del coefficiente d’attrito e forza al variare di Re.
f vs Re
10000000
2874481
1000000
287415.2
100000
28741.15
10000
2874.171
f
1000
287.5077
100
29.32876
10
1
0.000001
0.1
4.419023
0.00001
0.0001
0.001
0.01
Re
Figura 24 Coefficiente d’attrito-Re.
34
0.1
1
10
Figura 25 Coefficiente d’attrito-Re: Analisi di letteratura comparata ai risultati ottenuti con
FlexPDE.
3.6.
Caratterizzazione del coefficiente di attrito per variazione del raggio e
velocità
La versatilità del codice può essere mostrata andando a modificare anche gli altri parametri che
caratterizzano il sistema in esame. Un diverso tipo di indagine può essere effettuato fissando la
velocità di avvicinamento del fluido e facendo variare, invece, il raggio del disco. In tale modo il
numero di Reynolds risulta univocamente determinato.
Sono state effettuate simulazioni con il raggio del disco incrementato rispetto al caso standard di:


Raggio disco = 1.5 volte raggio standard (tabella 3)
Raggio disco = 2 volte raggio standard (tabella 4)
È possibile notare che anche in questo caso il coefficiente d’attrito che si ottiene cade sulla retta. I
diversi risultati sono ancora correlabili al grafico dei coefficienti di attrito e sono evidenziati con
colori diversi (fig. 26).
35
v0 [m/s]
Re
Forza [N]
Superficie
di Coefficiente
2
Impatto [m ]
di Attrito f
2.06*10^-8
0.00015
6.944196*10^-13
0.017671
195020.8
2.06*10^-7
0.001499
7.377487*10^-12
0.017671
19500.3
2.06*10^-6
0.014989
7.378004*10^-11
0.017671
1950.2
2.06*10^-5
0.149895
7.384005*10^-10
0.017671
195.2
2.06*10^-4
1.498947
7.550484*10^-9
0.017671
19.9
2.06*10^-3
14.98947
1.152603*10^-7
0.017671
3.0
Tabella 3 Misure con raggio = 1.5 volte raggio standard.
v0 [m/s]
Re
Forza [N]
Superficie di Coefficiente
Impatto [m2]
di Attrito f
1*10^-8
0.0000970192 5.075280*10^-13 0.031416
320220.4
1*10^-7
0.000970192
5.075299*10^-12 0.031416
32022.2
1*10^-6
0.009701923
5.075166*10^-11 0.031416
3202.1
1*10^-5
0.097019231
5.074174*10^-10 0.031416
320.2
1*10^-4
0.970192308
5.095524*10^-9
0.031416
32.1
1*10^-3
9.701923077
6.551359*10^-8
0.031416
4.1
Tabella 4 Misure con raggio = 2 volte raggio standard.
36
Figura 26 Coefficiente d’attrito-Re per raggi differenti.
3.7.
Caratterizzazione del coefficiente di attrito per variazione del solo raggio a
parità di velocità v0
Un possibile confronto può essere fatto utilizzando lo stesso valore di velocità (quella standard),
facendo variare poi soltanto il valore del raggio del disco.
Raggio
[m]
Disco v0 [m/s]
Re
Forza [N]
S impatto [m2] f
0.05
1*10^-6
0.004851
2.47006*10^-11
0.007853982
6233.8
0.0625
1*10^-6
0.006064
2.8488*10^-11
0.012272
4675.9
0.075
1*10^-6
0.007276
3.58141*10^-11
0.017671
4017.2
0.0875
1*10^-6
0.008489
4.63382*10^-11
0.024053
3818.7
0.1
1*10^-6
0.009702
5.07517*10^-11
0.031416
3202.1
Tabella 5 Coefficiente d’attrito al variare del raggio a parità di velocità v0.
Dai risultati ottenuti è possibile notare che la forza di trascinamento totale aumenta all’incrementare
del raggio, ma la superficie di impatto cresce più velocemente di quanto aumenti la forza. Ne
consegue che il denominatore nell’espressione di f prevale e quindi esso si riduce progressivamente.
37
Possono essere effettuate delle modifiche allo script, così da consentire di incrementare i rapporti
dimensionali del disco rispetto all’ambiente. Tuttavia lo script non ha validità universale; infatti
ampliando le dimensioni del raggio eccessivamente il flusso non riesce a svilupparsi nel dominio
descritto (per esempio quando il raggio è 3 volte quello originale) e quindi la velocità di
avvicinamento non viene più raggiunta.
38
Capitolo 4
Trasporto energia
4.1.
Introduzione
L’analisi relativa al trasporto di energia termica ha permesso di studiare la capacità del sistema in
esame di scambiare calore; a questo proposito si è valutato il coefficiente di scambio termico h
(rappresentante l’attitudine di un conduttore termico ad essere percorso da corrente termica) ed il
numero di Nusselt Nu ( il quale esprime il rapporto tra il flusso di calore scambiato per convezione
e il flusso di calore scambiato per conduzione). Ovviamente si tratta di valori relativi al fluido, nel
nostro caso aria.
Rispettivamente h e Nu sono definiti come:
in FlexPDE:
Figura 27 h e Nu nella sezione DEFINITIONS
Come già visto nel paragrafo 2.5, FlexPDE ha permesso di definire il flusso di calore q come un
vettore di componenti qz e qr, restituendo in output i relativi diagrammi attraverso i comandi
contour e vector.
Figura 28 Flussi termici nella sezione DEFINITIONS.
Uno degli scopi di questo progetto è risalire al profilo di temperatura del fluido, la cui temperatura
di bulk è diversa (in particolare inferiore di 100 °C) da quella del disco. La temperatura del disco è
un dato del problema ed è fissata a 400 °C. Per come è definito il problema, unicamente nel fluido
(caso ideale) si instaurerà un profilo di temperatura legato all’equazione dell’energia (la cui
variabile incognita è appunto la temperatura). L’aria subirà infatti un rilevante riscaldamento nei
pressi del disco.
39
Figura 29 Equazione dell’energia nella sezione EQUATIONS.
Nei prossimi paragrafi sono riportati i risultati ottenuti.
4.2.
Analisi del profilo di temperatura
Figura 30 Grafico CONTOUR temperatura.
Il profilo di temperatura è influenzato dall’andamento della velocità del flusso d’aria, ciò è evidente
dalla scia a valle del disco; infatti, ci sono porzioni di spazio più grandi dove la temperatura è
omogenea. È da specificare che nel caso di fluido in quiete il profilo di temperatura risulterebbe
perfettamente simmetrico con il raggio. Dal grafico successivo (fig. 31) si può notare come la
temperatura del fluido nella sezione d’ingresso al sistema sia costante e valga 300°C (condizione al
contorno); questi tipi di grafici sono stati ottenuti tramite il comando elevation.
40
Figura 31 ELEVATION temperatura sezione ingresso dominio.
Figura 32 ELEVATION T lungo la coordinata assiale.
In prossimità del disco la temperatura ha un picco: essa sale molto rapidamente da 300°C fino a
circa 400°C per poi scendere in maniera più graduale una volta che il flusso ha oltrepassato il disco
(fig. 32).
41
Il grafico seguente (fig. 33) mostra che l’effetto riscaldante del disco è rilevante anche a distanza
raggio_condotto/2 dal centro del disco.
Figura 33 ELEVATION T lungo l’asse z.
42
4.3.
Analisi del flusso di calore: qz
Figura 34 CONTOUR qz.
Tramite il comando contour è stato possibile conoscere il valore della componente parallela all’asse
z del flusso di calore.
Come si può evincere dal grafico la superficie laterale del disco, ovviamente, non da contributi; il
fatto che il flusso a monte del disco assuma valori minori di zero è giustificato dalla scelta del
sistema di riferimento.
43
4.4.
Analisi del flusso di calore: qr
Figura 35 CONTOUR qr.
Con lo stesso criterio si è giunti alla visualizzazione dei valori che assume la componente radiale
del flusso di calore; ovviamente le superfici frontale e posteriore non danno alcun contributo.
44
4.5.
Analisi del flusso di calore
Figura 36 Grafico VECTOR del flusso di calore.
Tramite il comando vector si è giunti al grafico che mostra l’ andamento del flusso di calore q in
qualità di vettore di componenti qz e qr. È intuitivo affermare che il flusso è uscente dal disco e
massimo alla parete.
4.6.
Effetti della variazione del numero di Reynolds
Per generalizzare i risultati ottenuti, si è proceduto facendo variare il numero di Reynolds (tabella 6)
dal quale è stata ricavata la velocità di bulk del fluido.
Re
Nu
Nu_ corr
h [=] W/m2
h_corr
Nu_corr/Nu
Pr
Re^0.5*Pr^(1/3)
1
2.250061
2.534245
0.668268
0.752671
1.126301
0.705939
0.890408
0.1
2.202884
2.168943
0.654256
0.644176
0.984592
0.705939
0.281572
0.01
2.202572
2.053424
0.654164
0.609867
0.932285
0.705939
0.089041
0.001
2.202595
2.016894
0.654171
0.599018
0.91569
0.705939
0.028157
0.0001
2.202598
2.005342
0.654172
0.595587
0.910444
0.705939
0.008904
0.00001
2.202598
2.001689
0.654172
0.594502
0.908785
0.705939
0.002816
Tabella 6 Studi al variare del numero di Reynolds.
45
Come si può evincere dalla tabella, al variare del numero di Reynolds il numero di Nu aumenta
leggermente ma resta sempre dell’ ordine di 2.2; allo stesso tempo aumenta il coefficiente di
scambio h, ma in maniera non molto rilevante. Ciò è in accordo con la fisica del sistema.
Dal grafico successivo si può notare che la portata di calore dipende linearmente da Re.
Q= 0.0269Re+ 1.2326
R² = 0.9916
1.23311422
Q [=] W
1.23311422
1.25968518
1.23309914
1.233112335
1.23327256
Q
Lineare (Q)
1
0.000001 0.00001
0.0001
0.001
0.01
0.1
1
Re
Figura 37 Q-Re.
All’aumentare del numero di Reynolds il profilo di temperatura è sempre più influenzato dalla
fluidodinamica del sistema. Infatti, si sviluppano scie che impediscono al fluido di riacquisire la
temperatura di bulk.
Sono di seguito riportati i profili di temperatura rispettivamente a Re uguale a 0.00001, 0.01 e 1.
Re
Figura 38 Effetto Re sul profilo di T.
46
È stata sfruttata inoltre l’analogia fluidodinamica che esiste tra disco e sfera, in regime laminare,
per collegare il coefficiente di scambio termico del disco con quello della sfera, la cui correlazione
presente in letteratura [2] è:
Per un fluido in quiete (Re uguale a 0) il numero di Nusselt è uguale a 2. In questo caso Re è molto
basso ed è presente il coefficiente correttivo 0.6 che fa in modo che Nu converga a 2. Valutando
quindi il Nu della correlazione per la sfera (Nu_ corr), sfruttando il Re e il Pr del sistema in esame,
si è visto che i valori discostavano di poco da quelli ottenuti da FlexPDE (vedi 6a colonna della
tabella).
Lasciando la dipendenza di Nusselt da
, l’equazione è stata riarrangiata sostituendo i
fattori correttivi 2 e 0.6 con altri coefficienti, ottenendo una correlazione per disco investito valida
per 0.00001<Re<1.
Successivamente è riportato un grafico logaritmico che mostra l’andamento di Nu contro
:
10
Nu= 0.0534Re1/2 Pr1/3 + 2.199
R² = 0.9124
correlazione disco investito
Nu
Lineare (correlazione disco
investito)
1
0.001
0.01
0.1
1
Re1/2 Pr1/3
Figura 39 Correlazione disco investito.
Il coefficiente di determinazione R² è una misura della bontà dell'adattamento (in inglese fitting)
della regressione lineare stimata ai dati osservati.
Per il caso studiato la linea di tendenza presenta un R2=0.9124, il che vuol dire che al massimo la
correlazione può dare valori che discostano da quelli reali di circa il 9%.
47
4.7.
Effetti della variazione del raggio del disco
Cambiamenti del raggio del disco permettono di comprendere al meglio il comportamento del
sistema quando se ne varia la geometria. In particolare sono stati studiati i casi in cui il raggio del
disco è 1.5 e 2 volte quello di riferimento, tenendo costante la velocità del fluido. Ovviamente il
sistema risente del cambiamento di geometria sia dal punto di vista fluidodinamico che energetico.
Raggio
[=]m
Re
v0 [=]
m/s
Nu
h[=] W/m2K
Superficie
[=] m2
Q[=] W
q[=] W/m2
r=0.05
0.0001
2.06E-08
2.202598
0.654172
0.01885
1.23308
65.41716
r*1.5
0.000149895
2.06E-08
2.44495
0.4841
0.040055
1.93907
48.41001
r*2
0.00019986
2.06E-08
2.44875
0.363639
0.069115
2.51329
36.36394
Tabella 7 Effetto raggio del disco
Per esempio Reynolds è direttamente proporzionale al raggio, ma ciò non va ad influenzare il
coefficiente di scambio in quanto la velocità resta costante.
Andando ad analizzare i risultati ottenuti da FlexPDE si è visto che l’ aumento della superficie
disponibile allo scambio è preponderante rispetto l’ aumento della quantità scambiata o portata di
calore Q. . Il fattore che tiene conto della proporzionalità tra la portata di calore e la superficie di
scambio è proprio il coefficiente di scambio, definito appunto come
con
costante; esso quindi diminuisce al crescere del raggio, da un punto di vista matematico il
denominatore predomina sul numeratore, da un punto di vista fisico invece possiamo dire che dischi
con raggi maggiori scambiano di più ma peggio (Q aumenta, h diminuisce).
Ovviamente il flusso di calore q decresce con lo stesso andamento di h essendo i due uguali a meno
di una costante ΔT.
Il tutto è riassunto nel grafico successivo (fig. 40).
48
h(raggio)
0.7
0.65
0.6
h
0.55
0.5
0.45
0.4
0.35
0.3
0.05
0.06
0.07
0.08
raggio
Figura 40 h-raggio
49
0.09
0.1
Conclusioni
Il presente lavoro ha mirato ad analizzare gli effetti fluidodinamici e termici che si hanno quando un
disco rigido viene investito da una corrente laminare di aria. Le proprietà dell’aria sono state prese
alla temperatura di film (media aritmetica tra quella di bulk e quella del disco).
Questo studio è stato portato avanti grazie alla simulazione di tale sistema con un software
specializzato nella risoluzione di equazioni differenziali alle derivate parziali. Per effettuare la
simulazione con questo programma è stato necessario costruire uno script strutturato in sezioni ben
precise.
Il punto di partenza per questa analisi è rappresentato dalla corretta scrittura delle equazioni di
bilancio e delle condizioni al contorno associate; questo sia nell’elaborazione del modello
matematico del sistema, sia nell’implementazione nel software. La risoluzione di dette equazioni ha
restituito i profili e gli andamenti delle principali variabili fluidodinamiche e termiche: P, v e T. Da
queste è stato possibile risalire alle quantità macroscopiche salienti, cioè la forza di trascinamento
Fz e la portata termica Q; da queste ultime, infine, sono stati ricavati il coefficiente di attrito e di
scambio termico. Soprattutto per i coefficienti di scambio sono stati ricercati andamenti
caratteristici al variare del numero di Reynolds e del raggio del disco, cioè delle condizioni di flusso
e della geometria. I risultati, confrontati con quelli di letteratura, hanno mostrato che sia il
coefficiente d’attrito che quello di scambio termico hanno andamenti correlabili, rispettivamente,
con il numero di Reynolds e con il numero di Reynolds e quello di Prandtl. I valori di f trovati sono
sovrapponibili a quelli trovati per via sperimentale; quelli di h invece, si discostano di poco da
quelli ottenuti con la correlazione di letteratura. L’andamento col raggio ha dato risultati spiegabili
sia matematicamente che fisicamente.
50
Bibliografia
[1] J.H. Perry-Chemical Engineering’s Handbook
[2] Fenomeni di trasporto, Bird R.B., Steward W.E., Lightfoot E.N. (1970)
[3] Manuale FlexPDE
[4] Simple fields of physics by finite element analysis, B.Gunnar Backstrom
[5] Fluid Mechanics and Heat Transfer, Knudsen e Katz (1958)]
51

Download Report