Circolazione Naturale
Circolazione Naturale
Circolazione Naturale
Facoltà di Ingegneria
L’ obiettivo di questo esercizio è il calcolo della portata in massa di miscela liquido vapore
all’interno di una caldaia a circolazione naturale, ovvero il suo punto di funzionamento.
Figura 1
La circolazione dell’acqua fra i fasci vaporizzatori ed il corpo separatore può essere mantenuta per
convezione naturale solamente se esistono le condizioni idonee.
Indicando con q, che consideriamo costante, la potenza termica fornita all’unità di lunghezza del
fascio tubiero, abbiamo che la potenza termica assorbita alla generica ascissa l sarà q ⋅ l , mentre si
indicherà con Q = q ⋅ l la potenza termica totale assorbita nell’attraversamento di tutto il fascio di
altezza L .
L’insieme dei fasci vaporizzatori ed i tubi di discesa, costituiscono un sistema di vasi comunicanti;
in un ramo (fascio vaporizzatore) il fluido è di densità minore rispetto all’altro (tubi di discesa).
Ciò significa che la miscela acqua-vapore, più leggera, tenderà a salire, mentre l’acqua, nei tubi di
discesa, tenderà a scendere, determinando una portata di circolazione G.
Dall’ equazione d’equilibrio si può ricavare la forza risultante, e, di conseguenza, la differenza di
pressione che mantiene il fluido in movimento:
1 1
∆p m = g ⋅ ∫ − dl in cui vL è il volume specifico della fase liquida mentre v è il volume
L L
v v
specifico della miscela.
Il volume specifico della miscela è, ovviamente, dipendente dal titolo, e, quindi, aumenta
all’aumentare di l .
Indicando con vV il volume specifico del vapore si ha:
v = vV x + v L (1 − x ) = v L + (vV − v L )x
2
Il titolo aumenta man mano che si sale lungo il fascio tubero, in ragione del fatto che la potenza
termica assorbita determina l’evaporazione di parte della portata di circolazione. Possiamo scrivere:
.
ql
ql = G r x → x =
Gr
I dati di partenza che abbiamo utilizzato per il problema sono i seguenti.
Valori termodinamici:
Kg v
Gv = 50000 (portata di vapore in uscita dalla caldaia)
h
Pv = 50 bar (pressione di funzionamento)
. Kcal
Q IRR = 13828368 ≈ 16083 KW (potenza irradiata; si è considerato che all’interno di questa
h
parte di caldaia l’unica potenza in gioco sia quella irradiata)
3
Øi = 55 mm (diametro interno dei tubi)
Øe = 60 mm (diametro esterno dei tubi)
t
= 1 (rapporto tra passo e diametro dei tubi sulla superficie 1)
d
t
= 1,5 (rapporto tra passo e diametro dei tubi sulle superfici 2, 3 e 4)
d
SVOLGIMENTO:
Abbiamo considerato un range di valori plausibili per le velocità della miscela liquido vapore e del
liquido soltanto (dettati da considerazioni empiriche) ed un valore preciso di α (grado di vuoto) per
evitare il rischio di surriscaldamenti. I valori sono i seguenti:
m
Vplausibile = 1 ÷ 10 (velocità della miscela liquido-vapore)
s
m
Vplausibile = 1 ÷ 2,5 (velocità del liquido nel tratto di discesa)
s
α = 0,7
µH2O = 0,0001 Pa s
µVap = 0,000018 Pa s
L
1 1
∆pm = g ∫ ( − )dl Equazione 1
0
vl v
C’è una correlazione tra α ed x (titolo della miscela), ove il pedice v sta ad indicare vapore, il pedice
l sta ad indicare liquido ed il pedice f si riferisce alle condizioni finali della miscela:
x f ⋅ vv
α= Equazione 2
(1 − x f ) ⋅ vl + x f ⋅ v v
4
da cui si ottiene:
1
xf = = 0,070703 Equazione 3
0,3 vv
1+ ⋅
0,7 vl
L1 = B
L3 = 2a + b
Calcoliamo la lunghezza complessiva dei tubi sulle superfici 1 e 3 (le uniche che abbiamo
considerato nel nostro calcolo, semplificando così il lavoro):
− 1297,8
L = m = 6,69 m
194
Q = 16083 KW risulta:
Q W
q= _
= 2404048
m
L
q⋅l
x= Equazione 5
r ⋅G
q ⋅l ⋅ a
v = v L (1 + ) Equazione 6
G
Secondo quanto emerso a lezione, la formula per il calcolo della pressione motrice può essere scritta
come segue:
g⋅L G q ⋅ a ⋅ L
∆pm = 1− ln1 + Equazione 7
vl q ⋅ a ⋅ L G
5
Basandoci su questa abbiamo proceduto a riportare i valori delle pressioni motrici in funzione della
portata in massa (variabile da 1 a 500 kg/s); qui di seguito è possibile vedere le prime righe della
tabella da noi ottenuta.
G ∆p Motrice
(kg/s) (Pa)
1,00 50037,6
2,00 49283,9
3,00 48620,7
4,00 48016,0
5,00 47454,4
6,00 46926,9
7,00 46427,6
8,00 45952,4
9,00 45497,9
10,00 45061,8
11,00 44642,2
12,00 44237,4
13,00 43846,0
14,00 43467,1
15,00 43099,5
16,00 42742,6
17,00 42395,5
18,00 42057,6
19,00 41728,4
20,00 41407,3
Tabella 1
Procediamo ora con il calcolo delle ∆p resistenti, partendo dai valori di velocità nei tubi.
Per il calcolo della zona di discesa del liquido, abbiamo ipotizzato un numero di tubi N = 50 di
primo tentativo.
Nel tratto di salita, la velocità è:
G ⋅v d2
w= dove Ω1 tubo = π e Ω = 194 ⋅ Ω 1 tubo
Ω 4
Abbiamo quindi calcolato le velocità massime e minime in corrispondenza dei valori massimi e
minimi di titolo, rispettivamente allo sbocco e all’imbocco dei tubi di salita. Anche queste sono
state parametrizzate a seconda dei valori di portata ipotizzati nelle tabelle di Excel.
La velocità massima dipende dal volume specifico della miscela massimo, valore per il quale
abbiamo considerato il titolo finale:
1
xf = 0,070703 Equazione 3.
0,3 vv
1+ ⋅
0,7 vl =
m3
vd = vv − vl = 0,038154 in corrispondenza del quale si avrà:
Kg
G ⋅ v max
wmax =
Ω
6
Per quanto riguarda i valori minimi di velocità, si è considerato il volume specifico del liquido in
quelle condizioni di pressione e temperatura (si vedano tabelle).
Per i tubi di discesa, le velocità sono state calcolate sulla sezione totale di discesa, considerando il
volume specifico del liquido e sempre parametrizzate secondo i valori della portata in massa.
L’espressione delle perdite di pressione si può indicare come segue:
ove sono state trascurate le variazioni di pressione dovute all’accelerazione del fluido nel suo
cambiamento di fase (e quindi di volume specifico). Le perdite di pressione distribuite nel tratto in
salita sono dovute ad attriti interni al fluido e si possono calcolare con la seguente formula:
G2 f 1 q a L2
∆pdist = Z v L + Equazione 8
Ω 2 D 2 G 2
l
ε
ove f si ottiene dal diagramma di Moody e richiede la conoscenza dei valori di Re, di e Z,
D
fattore di Martinelli, che si ottiene in funzione del titolo finale e della pressione in caldaia, come si
può vedere dal diagramma seguente:
7
che ci ha portato ad un valore di scabrezza relativa di:
ε
= 0,0000008
D
quindi all’interno dei valori di tubo liscio. Per permettere poi al sistema di calcolare
automaticamente i valori di f in base alla portata G parametrizzata (i valori di portata influenzano la
velocità che influisce sul numero di Re il quale a sua volta modifica il valore di f) abbiamo
provveduto ad interpolare alcuni valori della curva di tubi lisci, come risulta dal
Grafico 1, e a calcolare la linea di tendenza mediante equazione esponenziale; basandoci su questa,
abbiamo calcolato la relazione tra Re ed f, ottenendo così un calcolo automatico da parte di Excel
dell’f che compete ad ogni valore di G.
Simulazione Curva di Moody
-0,2058
y = 0,2038x
2
0,045 R = 0,991
0,04
0,035
0,03
f 0,025
Curva di Moody
Potenza (Curva di Moody)
0,02
0,015
0,01
0,005
0
0 1000000 2000000 3000000 4000000 5000000 6000000 7000000 8000000
Re
Grafico 1
Per quanto riguarda Z, si veda il Diagramma 1, ove si può rilevare l’andamento di questo valore in
funzione del titolo e della pressione.
Le perdite di pressione distribuite per i tubi di discesa sono di più facile computazione e si possono
esprimere mediante l’equazione:
2
f TD 1 W
∆ptubo disc = ρL L L Equazione 9
DTD 2 2
8
tratto di salita e di discesa. I valori di ξ sono stati considerati uguali ad 1 e 0.5, rispettivamente per
lo sbocco e l’imbocco delle tubazioni. Le equazioni per le perdite concentrate, quindi, sono state
prese direttamente dalla meccanica dei fluidi e risultano:
Wmin
2
Wmax
2
Wmin
2
Wmin
2
∆pconc = ξ I ρ L + ξ S ρV + ξ I ρ L + ξS + υL Equazione 10
2 2 SALITA 2 2 DISCESA
Una volta calcolate le perdite concentrate per ogni accidentalità abbiamo sommato tutti i valori di
∆p concentrati, tutti quelli distribuiti e la loro somma; confrontando questa con il valore di ∆p
motrice calcolato precedentemente, ne abbiamo calcolato la differenza ed abbiamo ricavato il punto
di funzionamento (rimanendo nell’ordine di definizione di 1 kg/s). Nella Tabella 2 è possibile
vedere una parte della schermata del file Excel che ci ha portato alla soluzione del problema con i
dati di partenza e ad un valore di portata di funzionamento di 324 kg/s.
G ∆p Motrice ∆p Resistente
(kg/s) (Pa) (Pa) Differenza WSmin WSmax WD
1,00 50037,6 11,4 50026,2 0,00 0,01 0,01
2,00 49283,9 18,4 49265,5 0,01 0,02 0,02
3,00 48620,7 26,1 48594,6 0,01 0,03 0,03
0 4,00 48016,0 33,7 47982,3 0,01 0,03 0,04
0 5,00 47454,4 41,3 47413,1 0,01 0,04 0,05
0 6,00 46926,9 49,0 46877,9 0,02 0,05 0,06
0 7,00 46427,6 56,8 46370,9 0,02 0,06 0,08
0 8,00 45952,4 64,8 45887,6 0,02 0,07 0,09
0 9,00 45497,9 72,9 45425,0 0,03 0,08 0,10
0 10,00 45061,8 81,3 44980,5 0,03 0,09 0,11
0 320,00 14734,7 14374,4 360,3 0,89 2,77 3,46
0 321,00 14704,9 14457,3 247,6 0,90 2,77 3,48
0 322,00 14675,2 14540,4 134,8 0,90 2,78 3,49
0 323,00 14645,6 14623,7 21,9 0,90 2,79 3,50
324 324,00 14616,1 14707,2 -91,1 0,90 2,80 3,51
0 325,00 14586,8 14791,0 -204,2 0,91 2,81 3,52
0 326,00 14557,6 14875,0 -317,4 0,91 2,82 3,53
0 327,00 14528,4 14959,2 -430,8 0,91 2,83 3,54
0 328,00 14499,5 15043,7 -544,2 0,92 2,83 3,55
0 329,00 14470,6 15128,4 -657,8 0,92 2,84 3,56
0 330,00 14441,8 15213,3 -771,5 0,92 2,85 3,57
Tabella 2: Soluzione 1
Nella colonna all’estrema sinistra di Tabella 2 è possibile vedere un contatore che indica il valore di
portata in corrispondenza del quale si ha il cambiamento di segno della differenza tra la ∆p motrice
e quella resistente, ovvero il valore intero di portata immediatamente successivo a quello che
identifica il punto di funzionamento, questo valore ci è servito per il calcolo esatto del titolo
massimo raggiunto dalla miscela. Nelle colonne all’estrema destra di Tabella 2 è possibile vedere i
valori di velocità minima e massima nel tratto in salita, nonché la velocità del tratto di discesa. Il
risultato ottenuto è rappresentato nel
Grafico 2.
9
Curve delle Pressioni
20000
19000
18000
17000
16000
Pressioni (Pa)
15000
14000
13000
12000
11000
10000
250 260 270 280 290 300 310 320 330 340 350 360 370 380 390 400 410 420 430 440 450
Portate (kg/s)
Grafico 2: Soluzione 1
10
SOLUZIONE 2:
Da considerazioni sulle velocità raggiunte dal liquido nel tratto di discesa abbiamo pensato di
provare ad aumentare il diametro dei tubi in questa sezione per limitare questi valori. Utilizzando
quindi lo stesso foglio Excel (modificando soltanto i dati iniziali come visibile in
Tabella 3) siamo giunti alla soluzione visibile in Grafico 3.
Tabella 3: Soluzione 2
11
Curve delle Pressioni
17000
16000
15000
14000
Pressioni (Pa)
13000
12000
11000
10000
9000
8000
350 360 370 380 390 400 410 420 430 440 450 460 470 480 490 500
Portata (kg/s)
Grafico 3: Soluzione 2
12
SOLUZIONE 3:
Da ulteriori considerazioni sui valori di velocità abbiamo compreso che il singolo aumento di
diametro dei tubi di discesa non è sufficiente a ridurre questa entro valori accettabili; abbiamo
quindi provveduto ad aumentare il numero di tubi di discesa come visibile in
Tabella 4 e
Grafico 4.
U. di
Grandezza Valore misura
T= 264 °C
p= 50 bar
a= 0,000018 kg/J
r= 1640100 J/kg
vv = 0,03944 m^3/kg
vl = 0,001286 m^3/kg
g= 9,81 m/s^2
L= 6,69 m
q= 2404048 W/m
x = 0,0707 -
Ø interno tubo 0,055 m
Nsalita = 194 tubi
Sezione totsalita = 0,460911 m^2
Ndiscesa = 100 tubi
Sezione totdiscesa = 0,237583 m^2
vtit max = 0,003984 m^3/kg
ρl = 777,60 kg/m^3
ρtit max= 251,03 kg/m^3
ξ Imbocco = 0,5 -
ξ Sbocco = 1 -
µ acqua 0,0001 Pa *s
µ vapore 0,000018 Pa *s
µ miscela sbocco = 0,000094 Pa *s
ε acc comm 45*10^(-9) m
ε acc comm/D 0,000001 -
Z 4 -
Ø interno tubo discesa 0,055 m
x effettivo = 0,022439723
Tabella 4: Soluzione 3
13
Curve delle Pressioni
15000
14000
13000
Pressioni (Pa)
12000
11000
10000
9000
400,00 410,00 420,00 430,00 440,00 450,00 460,00 470,00 480,00 490,00 500,00
Portate (kg/s)
Grafico 4: Soluzione 3
Per verificare infine la bontà del modello da noi realizzato ed acquisire un po’ di dimestichezza con
le grandezze in gioco, abbiamo provato a calcolare i valori di portata all’interno di caldaie aventi
stessa geometria, ma pressioni funzionamento differenti. Si vedano a questo scopo i grafici e le
tabelle seguenti.
14
SOLUZIONE 4:
U. di
Grandezza Valore misura
T= 295,06 °C
p= 80 bar
a= 1,10971E-05 kg/J
r= 1441300 J/kg
vv = 0,02352 m^3/kg
vl = 0,001384 m^3/kg
g= 9,81 m/s^2
L= 6,69 m
q= 2404048 W/m
x = 0,12072575 -
Ø interno tubo salita 0,055 m
Nsalita = 194 tubi
Sezione totsalita = 0,460910523 m^2
Ndiscesa = 100 tubi
Sezione totdiscesa = 0,237582744 m^2
vtit max = 0,004056 m^3/kg
ρl = 722,54 kg/m^3
ρtit max= 246,52 kg/m^3
ξ Imbocco = 0,5 -
ξ Sbocco = 1 -
µ acqua 0,000088 Pa *s
µ vapore 0,000019 Pa *s
µ miscela sbocco = 0,000080 Pa *s
ε acc comm 45*10^(-9) m
ε acc comm/D 0,000001 -
Z 3 -
Ø interno tubo discesa 0,055 m
x effettivo = 0,029135072
Tabella 5: Soluzione 4
15
Curve delle Pressioni
14000,0
12000,0
10000,0
Pressioni (Pa)
8000,0
6000,0
4000,0
2000,0
0,0
300,00 310,00 320,00 330,00 340,00 350,00 360,00 370,00 380,00 390,00 400,00 410,00 420,00 430,00 440,00 450,00
Portate (kg/s)
Grafico 5: Soluzione 4
16
SOLUZIONE 5:
U. di
Grandezza Valore misura
T= 311,06 °C
p= 100 bar
a= 8,66646E-06 kg/J
r= 1317100 J/kg
vv = 0,018026 m^3/kg
vl = 0,001452 m^3/kg
g= 9,81 m/s^2
L= 6,69 m
q= 2404048 W/m
x = 0,158214252 -
Ø interno tubo salita 0,055 m
Nsalita = 194 tubi
Sezione totsalita = 0,460910523 m^2
Ndiscesa = 100 tubi
Sezione totdiscesa = 0,237582744 m^2
vtit max = 0,004074 m^3/kg
ρl = 688,71 kg/m^3
ρtit max= 245,44 kg/m^3
ξ Imbocco = 0,5 -
ξ Sbocco = 1 -
µ acqua 0,000081 Pa *s
µ vapore 0,00002 Pa *s
µ miscela sbocco = 0,000071 Pa *s
ε acc comm 45*10^(-9) m
ε acc comm/D 0,000001 -
Z 3 -
Ø interno tubo discesa 0,055 m
x effettivo = 0,035190143
Portata calcolata
= 347 kg/s
Tabella 6: Soluzione 5
17
Curve delle Pressioni
12000
11000
10000
9000
Pressioni (Pa)
8000
7000
6000
5000
4000
300 310 320 330 340 350 360 370 380 390 400
Portate (kg/s)
Grafico 6: Soluzione 5
18
CONCLUSIONI:
Il file Excel da noi ideato si è rivelato utile e sufficientemente versatile per il lavoro che ci eravamo
prefissi e per le successive modifiche. I valori di titolo che abbiamo calcolato, in base al dato di
input q, ottenuto da un esercizio precedente svolto in classe, si sono sempre dimostrati al di sotto
dei valori imposti dal grado di vuoto α desiderato. Dalle simulazioni ottenute nella soluzione 1,
abbiamo verificato un valore di velocità di discesa del liquido troppo elevato, che non siamo riusciti
a contenere neppure con un aumento del diametro dei tubi di discesa (soluzione 2) ed al quale
abbiamo dovuto ovviare con un aumento del numero di tubi di discesa (da 50 a 100, soluzione 3).
Come era prevedibile, abbiamo verificato che sia un aumento del diametro dei tubi, sia un aumento
del loro numero, portano il punto di funzionamento a valori di portata più elevati, proprio a causa
della diminuzione di velocità nel tratto di discesa, quindi delle perdite da calcolarsi nel suddetto
ramo. Un aumento di pressione in caldaia (come verificato nelle soluzioni 4 e 5), come era anche in
questo caso prevedibile, ci ha portato a valori di titolo più elevati (dovuti alla diminuzione del
calore latente di vaporizzazione con l’aumentare della pressione e temperatura in caldaia) e a valori
di portata più ridotti (in questo caso, la diminuzione di ∆pm è dovuta al calo del valore di volume
specifico differenziale, con l’aumentare di p e T, congiuntamente alla diminuzione di r, calore
latente di vaporizzazione).
Possibili ulteriori sviluppi di questo lavoro potrebbero comprendere un calcolo dei valori
termodinamici della miscela liquido vapore considerando non il titolo dedotto dal valore di α
accettabile in caldaia, ma calcolato nel punto più alto della stessa mediante l’equazione 5, sviluppo
che porterebbe a valori termodinamici della miscela più attendibili. Inoltre, risultati più precisi di
portata nel punto di funzionamento si otterrebbero se non si trascurassero la variazione di pressione
dovuta all’accelerazione del fluido lungo i tubi di salita per variazione di volume specifico e le
perdite concentrate situate in corrispondenza delle curve nelle tubazioni.
Non dimentichiamo, infine, che la nostra semplificazione non ha tenuto conto dei tubi disposti sulle
superfici 2 e 4 di Figura 2, e che una loro computazione avrebbe portato a risultati più attendibili.
Nonostante questo, i valori numerici da noi ottenuti appaiono più che plausibili, e le semplificazioni
adottate più che ragionevoli nel campo di studio nel quale ci siamo trovati ad operare.
19