83
Modelação de caldeira de condensação João Pedro Vasconcelos Barros Dissertação para obtenção do Grau de Mestre em Engenharia Mecânica Júri Presidente: Prof. Luís Rego da Cunha de Eça Orientador: Prof. João Luís Toste de Azevedo Co-Orientador: Eng. Luís Monteiro Vogal: Prof. Edgar Caetano Fernandes Outubro de 2012

Modelação de caldeira de condensação - Autenticação · potência perdida foram de 0,16 kW a 2,86 kW para a caldeira de condensação. Já para as caldeiras sem condensação

  • Upload
    lamnhan

  • View
    229

  • Download
    0

Embed Size (px)

Citation preview

Modelação de caldeira de condensação

João Pedro Vasconcelos Barros

Dissertação para obtenção do Grau de Mestre em

Engenharia Mecânica

Júri

Presidente: Prof. Luís Rego da Cunha de Eça

Orientador: Prof. João Luís Toste de Azevedo

Co-Orientador: Eng. Luís Monteiro

Vogal: Prof. Edgar Caetano Fernandes

Outubro de 2012

ii

“Scientists dream about doing great things. Engineers do them. “

James A. Michener, in Space, 1982

iii

Resumo

O presente trabalho visou a criação de um modelo matemático simplificado de uma caldeira de

condensação para análise da influência da humidade relativa com base em uma caldeira de

condensação da Bosch e posterior validação através da análise de ensaios experimentais já

realizados fornecidos, bem como generalização com recurso a outros modelos de caldeiras sem

condensação. Este modelo foi criado em VBATM .

O modelo criado assenta em quatro zonas: queimador, câmara de combustão, permutador de

calor primário e zona de condensação. O modelo considera combustão completa instantânea com as

perdas por radiação na chama. Os produtos de combustão são considerados cinzentos. O

permutador de calor primário foi modelado utilizando correlações existentes na literatura para alhetas

onduladas tipo herringbone. A zona de condensação foi considerada utilizando dois modelos, um que

considera condensação pura e outro que tem em conta a existência de gases não condensáveis.

O modelo matemático foi resolvido com recurso ao método Broyden, melhorado com recurso à

fórmula de Sherman-Morrison. Os resultados, obtidos com base nas diferentes caldeiras, mostram

que o modelo prevê o rendimento, para situações e humidade relativa nula, com um desvio padrão de

3,63%, enquanto em situações de humidade relativa de 100% o desvio-padrão é de 4,46%.

O efeito da humidade relativa varia consoante a potência da caldeira. As diferenças para a

potência perdida foram de 0,16 kW a 2,86 kW para a caldeira de condensação. Já para as caldeiras

sem condensação as diferenças para a potência perdida situam-se entre 0,72 kW a 0,99 kW.

Palavras-chave: caldeira de condensação, gases não condensáveis, humidade relativa, eficiência

iv

Abstract

In the present work a simplified mathematical model for a condensing boiler, was built using the

VBATM programming language. The model was configured to a Bosch condensing boiler and validated

through the analysis of experimental tests already carried out. In particular the influence of the air

humidity in the boiler efficiency is addressed. The model was also applied to other boilers without

condensation.

The model was created based on four zones: the burner, the combustion chamber, the primary

heat exchanger and the condensation zone. The energy balance in the combustion zone considers

complete conversion with heat losses from the flame. Combustion products are considered a gray gas

mixture. The primary heat exchanger was modeled using correlations for the convection coefficient for

herringbone corrugated fins. For the condensation zone two models were built. The first model

considers pure condensation while the second took into account the existence of non-condensable

gases.

The mathematical model was solved using the Broyden method improved with Sherman-Morrison.

The results obtained from the different boilers, show that the model predicts the boiler efficiency, in

null relative humidity, with a standard deviation of 3.63% and with 100% of relative humidity, the

standard deviation was of 4.46%.

The effect of the relative humidity varies depending on input power. The difference on power loss

was of 0.16kW to 2.86 kW for the condensing boiler. For the non-condensing boilers, a difference of

0.72 kW to 0.99 kW was obtained.

Keywords: condensing boiler, non-condensable gases, relative humidity, efficiency

v

Agradecimentos

Ao Professor Toste Azevedo pela orientação, paciência e ensinamentos que revelaram-se indispensáveis.

Ao Engenheiro Luís Monteiro pelos dados disponibilizados, sem os quais não seria possível elaborar a dissertação e pela disponibilidade em tirar qualquer dúvida.

Aos meus pais e à minha irmã, por tudo. Sem o apoio incondicional nunca teria sido possível chegar onde estou hoje.

Ao resto da minha família por todos os esforços realizados, preocupação e expectativa demonstrados nos maus e bons momentos.

Aos meus amigos madeirenses que sempre estiveram presentes nos bons e maus momentos.

Aos amigos daqui que sempre mostraram disponibilidade para ajudar a partilhar fardos e alegrias.

Aos meus colegas por todo o apoio prestado ao longo destes anos.

vi

Índice

1 Introdução ................................................................................................................................ 1

1.1 Objetivo .................................................................................................................................... 1

1.2 Trabalhos Anteriores ............................................................................................................... 1

1.3 Enquadramento ....................................................................................................................... 2

1.4 Estrutura da dissertação .......................................................................................................... 3

2 Revisão Bibliográfica de Correlações para Convecção .......................................................... 4

2.1 Interior dos tubos ..................................................................................................................... 4

2.2 Exterior dos tubos .................................................................................................................... 6

2.3 Convecção no condensador .................................................................................................. 11

2.3.1 Modelo 1 – Condensação de vapor ...................................................................................... 11

2.3.2 Modelo 2 - Condensação de vapor em mistura .................................................................... 13

3 Modelação de uma caldeira de condensação ....................................................................... 17

3.1 Modelo ................................................................................................................................... 17

3.2 Análise da Combustão .......................................................................................................... 18

3.3 Balanço de energia à câmara de combustão ........................................................................ 20

3.4 Balanço ao permutador de calor primário ............................................................................. 23

3.4.1 Método ε-NTU-θ ..................................................................................................................... 24

3.5 Balanço ao condensador ....................................................................................................... 24

3.6 Perdas na envolvente da caldeira ......................................................................................... 26

3.7 Propriedades dos fluídos ....................................................................................................... 27

3.7.1 Propriedades dos gases ........................................................................................................ 27

3.7.2 Propriedades da água ........................................................................................................... 28

3.8 Metodologia ........................................................................................................................... 30

3.8.1 Esquemática de resolução do modelo .................................................................................. 30

3.8.2 Método de resolução de sistemas de equações não lineares .............................................. 31

4 Resultados ............................................................................................................................. 36

4.1 Análise do método numérico ................................................................................................. 37

4.2 Análise da Combustão .......................................................................................................... 38

4.3 Análise das Perdas ................................................................................................................ 41

4.4 Análise do lado dos gases de combustão ............................................................................. 42

4.5 Análise do lado da água ........................................................................................................ 48

4.6 Outras caldeiras .................................................................................................................... 48

5 Conclusões ............................................................................................................................ 51

vii

6 Referências Bibliográficas ..................................................................................................... 53

Anexo 1. Outras correlações ................................................................................................................ 57

Anexo 2. Propriedades ......................................................................................................................... 59

Anexo 3. Instruções .............................................................................................................................. 65

Anexo 4. Ensaios Caldeiras .................................................................................................................. 67

Listafiguras

Figura 2.1 - Exemplo da influência dos tubos elípticos a diferentes ângulos de ataque para uma razão

de eixos de 0,5 (esq.) e 0,25 (dir.) - Ibrahim et al. (2009) ....................................................................... 7

Figura 2.2 – Tubos alhetados com alhetas onduladas do tipo herringbone com indicação de alguns

parâmetros geométricos – Kim et al. (2008) ........................................................................................... 8

Figura 2.3 - Evolução do número de Nusselt em função de Reynolds para Kim et al. (2008) (azul),

com parâmetros geométricos correspondentes ao permutador analisado e Xie (2009) (vermelho) .... 11

Figura 2.4 - Exemplo de excesso de condensado criado em um tubo alhetado – Webb et al. (1983) 12

Figura 2.5 - Exemplo de uma camada limite de difusão na presença de gases não condensáveis (NC)

– Shi et al. (2011) .................................................................................................................................. 14

Figura 2.6 - Ilustração da variação do coeficiente de convecção em função do caudal e fracção

mássica de ar a diferentes temperaturas do gás – Collier et al. (1994) ............................................... 16

Figura 3.1 - Modelo completo com ilustração esquemática das trocas de calor e dos fluxos de entalpia

............................................................................................................................................................... 17

Figura 3.2 - Exemplo de caldeira de condensação - Cortesia Vulcano ................................................ 18

Figura 3.3 - Zona convectiva - permutador de calor ............................................................................. 21

Figura 3.4 – Exemplo do balanço às paredes laterais + tubos ............................................................. 22

Figura 3.5 - Ilustração da zona de condensação .................................................................................. 25

Figura 3.7 - Representação das entradas e saídas das folhas de cálculo representadas a negrito .... 31

Figura 3.8 - Representação das entradas e saídas das folhas de cálculo representadas a negrito .... 31

Figura 3.9 - Ilustração do método de Broyden ...................................................................................... 32

Figura 4.1 - Número de iterações utilizando várias combinações do método de Broyden disponíveis 38

Figura 4.2 - Comparação do número de iterações necessário para tolerância de 108 ........................ 38

Figura 4.3 - Gráfico exemplificativo da evolução da temperatura adiabática com o excesso de ar para

HR=0% .................................................................................................................................................. 39

Figura 4.4 - Gráfico da evolução das frações molares com o excesso de ar para HR=0% ................. 39

viii

Figura 4.5 - Temperatura adiabática de chama vs. Humidade relativa ................................................ 40

Figura 4.6 - Comparação da influência dos combustíveis nos mecanismos de transmissão de calor. 41

Figura 4.7 - Perdas pela envolvente da zona de combustão vs. Humidade relativa ............................ 42

Figura 4.8 - Variação da temperatura na zona de combustão vs. Emissividade da parede em

situações extremas de humidade relativa ............................................................................................. 43

Figura 4.9 - Variação de temperatura no permutador de calor primário vs. Correlação usada ............ 43

Figura 4.10 - Gráfico da temperatura de saída do permutador primário e de saída dos gases de

combustão para os ensaios #3 (linhas vermelhas) e #4 (linhas amarelas) e a tracejado os valores

medidos ................................................................................................................................................. 46

Figura 4.11 -Gráfico com o rendimento da caldeira de condensação vs. HR para os ensaios #3 (linha

vermelha) e #4 (linha amarela) e a tracejado os rendimentos medidos ............................................... 47

Figura 4.12 - Comparação de rendimentos experimentais com rendimentos calculados para +-10%

para HR=0% .......................................................................................................................................... 49

Figura 4.13 - Comparação de rendimentos experimentais com rendimentos calculados para +-10%

para HR=100% ...................................................................................................................................... 50

Figura A.1 - Exemplo da variação da massa específica com temperatura e humidade relativa -

Tsilingiris (2007) .................................................................................................................................... 64

Figura A.2 - Exemplo da folha "Solve" em configuração 0 ................................................................... 66

Listatabelas

Tabela 2.1 - Limite de validade das correlações para escoamentos internos ........................................ 5

Tabela 2.2 - Correlações com o uso de turbulators ................................................................................ 5

Tabela 2.3 - Limite de validade da correlação de Kim et al. (2008) ........................................................ 9

Tabela 2.4 - Limites de validade da correlação de Xie et al. (2009) ..................................................... 11

Tabela 3.1 - Limite de validade das expressões para a água............................................................... 30

Tabela 4.1- Várias configurações na análise do gás ............................................................................ 36

Tabela 4.2 - Resultados dos ensaios #1, #2 e #5 ................................................................................. 44

Tabela 4.3 - Temperatura dos gases de combustão em diferentes posições para o ensaio #1

(H.R.=50%) ............................................................................................................................................ 45

Tabela 4.4 - Comparação calor total trocado e calor latente para ensaio #3 ....................................... 47

Tabela 4.5 – Entalpias de saída dos gases de combustão e diferença obtida em unidades de potência

para diferentes ensaios ......................................................................................................................... 48

ix

Tabela 4.6 - Temperatura de saída da água da caldeira de condensação para diferentes ensaios em

situações extremas de humidade .......................................................................................................... 48

Tabela 4.7 - Diferença obtida para os gases de combustão na saída para os diferentes ensaios ...... 50

Tabela A.1 - Limite de validade das correlações para escoamento exterior ........................................ 57

Tabela A.2 - Tabela com dados da correlação de Leckner .................................................................. 58

Tabela A.3 - Coeficientes de �� ............................................................................................................ 59

Tabela A.4 - Coeficientes de � ........................................................................................................... 59

Tabela A.5 - Coeficientes de �� ............................................................................................................ 59

Tabela A.6 - Coeficientes de �� ............................................................................................................ 60

Tabela A.7 - Coeficientes da equação de Shomate .............................................................................. 60

Tabela A.8 - Coeficientes da equação de Shomate .............................................................................. 60

Tabela A.9 - Coeficientes dos polinómios dos gases de combustão para uma dada gama de

temperaturas ......................................................................................................................................... 61

Tabela A.10 - Coeficientes dos polinómios dos combustíveis para uma dada gama de temperaturas 61

Tabela A.11 - Correlações usadas no ensaio #1 .................................................................................. 67

Tabela A.12 - Correlações alteradas para o ensaio #3 ......................................................................... 67

Tabela A.13 - Correlações alteradas para o ensaio #4 ......................................................................... 68

Tabela A.14 - Correlações alteradas para o ensaio #5 ......................................................................... 68

x

Nomenclatura

� – Maior diâmetro de uma elipse / 2 (m)

� – Área exposta dos tubos (sem alhetas)

(m2)

� – Área total exterior do tubo (m2)

Af – Área da alheta (m2)

� – Área do interior do tubo (m2)

��� – Área da região conectiva do modelo 2

(m2)

��� – Área da envolvente da zona combustão

(m2)

���� – Área passagem pelo permutador de

calor (m2)

������� – Área frontal do permutador de calor

(m2)

���� – Área da região da raiz do tubo do

modelo 2 (m2)

�� – Área da superfície radiativa i (m2)

� – Matriz ligada ao método de Broyden

� – Menor diâmetro de uma elipse /2 (m)

� – Razão do passo da mola (!/#���) �%���� – Concentração molar dos gases

(kmol/m3)

�� – Razão das capacidades caloríficas

c – Parâmetro do modelo 2 (ad.)

&�,� – Calor específico do ar seco (kJ/kgK)

&�,� – Calor específico da água (kJ/kgK)

&�( – Calor específico do ar húmido (kJ/kgK)

&�,� – Calor específico do vapor de água

(kJ/KgK)

&�,%���� – Calor específico dos gases (kJ/kgK)

D – Coeficiente de difusividade térmica (m2/s)

)* – Desvio-padrão

#� – Diâmetro do colar do tubo, #� = #,� + 2/�

(m)

#��� – Diâmetro da espiral (m)

#,� – Diâmetro hidráulico do tubo (m)

#� – Diâmetro exterior do tubo (m)

#0 – Diâmetro exterior da alheta circular (m)

#� – Diâmetro do interior do tubo (m)

#�� – Diâmetro do interior da alheta anular (m)

#1 – Direcção de Broyden

2� – Poder emissivo da superfície i (W)

3��� – Comprimento da região convectiva ao

longo da alheta (modelo 2) (m)

3� – Comprimento da região da alheta ao

longo da alheta (modelo 2) (m)

3���� – Comprimento da região da raiz ao

longo da alheta (modelo 2) (m)

3� – Razão do passo da asa em delta (*/4) 35 – Razão da largura da asa em delta (6/4) 7 – Vector das funções do método de Broyden

8�9 – Factor de forma entre as superfícies i e j

8� – Espaçamento das alhetas (m)

8� – Passo das alhetas (m)

: – Fator de atrito (ad.)

;� – Irradiação da superfície i (W)

g – Aceleração da gravidade =9,81A !B⁄ )

xi

D , �9 – Parâmetros ad. relativos a

propriedades

– Humidade relativa (%)

ℎ�� – Entalpia do ar húmido (kJ/kg ar seco)

ℎF – Coeficiente de convecção do excesso de

condensado (W/m2K)

ℎ� – Entalpia de condensação (relativa a vapor

de água a G��� = 25℃) (kJ/kg)

ℎ��(F – Coeficiente de convecção da câmara

de combustão (W/m2K)

ℎ��� – Coeficiente de convecção da região

convectiva (modelo 2) (W/m2K)

ℎ���J – Coeficiente de convecção por

condensação (modelo 2) (W/m2K)

ℎKL – Incremento usado pelo método das

diferenças finitas

ℎ� – Coeficiente de convecção no exterior

(W/m2K)

ℎ� – Coeficiente de convecção das alhetas

(W/m2K)

ℎ�� – Coeficiente de convecção da zona

inundada (W/m2K)

ℎ��( – Coeficiente de convecção do filme de

líquido (W/m2K)

ℎ�% – Entalpia de vaporização (kJ/Kg)

ℎ%� – Entalpia dos gases de combustão no

ponto i (W/m2K)

ℎ, - Coeficiente de convecção do tubo

horizontal (W/m2K)

ℎ� – Coeficiente de convecção no interior dos

tubos (W/m2K)

ℎ�% – Coeficiente de convecção dos gases não

condensáveis e vapor de água (W/m2K)

ℎ��� – Entalpia de referência (kJ/kg)

ℎ���� – Coeficiente de convecção da região da

raiz do tubo (modelo 2) (W/m2K)

ℎM�� – Coeficiente de convecção da zona não

inundada (W/m2K)

ℎ� – Entalpia do fluido gasoso (kJ/kg)

ℎ5 – Entalpia da água líquida (ℎ5 = ℎ�) (kJ/kg)

ℎN – Coeficiente de convecção médio exterior

do condensador (W/m2K)

ℎN��� – Coeficiente de convecção médio

exterior da câmara combustão (W/m2K)

ℎN� =G) – Entalpia de formação à temperatura T

e pressão p (KJ/Kmol)

ℎN�0=G) – Entalpia de formação à temperatura T

e pressão atmosférica (KJ/Kmol)

O� – Número de Jakob (ad.)

O� – Radiosidade da superfície i (W)

� – Factor de Colburn P� = ,�QR Pr

UV = WM

X�YZ[ V⁄ \ �� - Condutibilidade térmica do ar seco

(W/mK)

�% - Condutibilidade térmica do gás (W/mK)

�� - Condutibilidade térmica da água (W/mK)

�( - Condutibilidade térmica do ar húmido

(W/mK)

�� - Condutibilidade térmica do vapor de água

(W/mK)

�� – Condutibilidade térmica da alheta (W/mK)

�� – Condutibilidade térmica da superfície i

(W/mK)

xii

�� – Condutibilidade térmica do tubo (W/mK)

�N- Condutibilidade térmica normalizada (ad.)

]0, 0, ]_, _ – Funções modificadas

(hiperbólicas) de Bessel

� – Comprimento da parede lateral (m)

�� – Comprimento da alheta (altura) (m)

�� , ��9 – Parâmetros ad. relativos às

propriedades

�( – Espessura média da camada de gás na

câmara de combustão (m)

���� – Meia distância entre dois tubos da zona

radiativa (m)

�� – Comprimento do tubo (m)

�3 – Número de Lewis (ad.)

` – Massa molar (kg/kmol)

A – Parâmetro da alheta seca (m-1)

Aa �� – Caudal mássico de ar de entrada (kg/s)

Aa � – Caudal mássico de condensado (kg/s)

Aa �M – Caudal mássico de gás

natural/combustível (kg/s)

Aa %���� – Caudal mássico dos gases de

combustão (kg/s)

Aa %����,� – Caudal mássico dos gases de

combustão após condensador i (kg/s)

Aa %����,�� – Caudal dos gases sem conteúdo de

vapor de água (kg/s)

Aa ����� – Caudal mássico de vapor de água à

temperatura do gás (kg/s)

Aa 5 – Caudal mássico de água nos tubos

(kg/s)

bGc – Número de unidades de transferência

(ad.)

bd – Números de tubos em cada fila

(transversalmente) (ad.)

be - Número de filas de tubos (ad.)

b� – Número de alhetas para o comprimento

analisado, b� = f/8� (ad.)

bg – Número de Nusselt médio hbg = ,eD i (L –

dimensão característica) (ad.)

j – Número de passagens no permutador (ad.)

j� – Número de moles da espécie i (Kmol)

j� – Número de moles de vapor de água

(Kmol)

* – Passo da fita a cada 180º (m)

*�^ – Poder calorifico in ferior (MJ/kg comb.)

*�k – Poder calorifico inferior (MJ/kg comb.)

* – Razão do passo (1 + !/#�) (ad.)

*3 – Número de Péclet

*l – Número de Prandtl

*K – Altura da forma da alheta excluindo /� (m)

*d – Passo transversal dos tubos (m)

*e – Passo longitudinal dos tubos (m)

*e∗ – Passo longitudinal equivalente dos tubos

(m)

*( – Perímetro molhado (m)

n – Pressão atmosférica (MPa)

n���@d – Pressão de saturação à temperatura

T (MPa)

p� – Taxa de transferência de calor

condensador n (W)

xiii

pqq – Fluxo de transferência de calor (W)

p���,����� – Taxa de transferência da

envolvente à câmara combustão (W)

0 – Constante dos gases ( 0 = 8,1344O/Auv])

� – Número de Rayleigh (ad.)

3Jw – Número de Reynolds baseado no

diâmetro #� (ad.)

3Jx – Número de Reynolds baseado no

diâmetro #, (ad.)

3Jy – Número de Reynolds baseado no

diâmetro #� (ad.)

3e – Número de Reynolds baseado no

comprimento da placa (ad.)

l_ – Raio interior do tubo (m)

lB – Raio exterior do tubo (m)

lB� – Raio exterior do tubo corrigido (m)

l�,�z – Raio do círculo equivalente (m)

! – Passo (m)

k& – Número de Schmidt hk& = {|Ki

!1 – Vector de Δf do método de Broyden

k�J – Comprimento do condensado convexo

(m)

k& – Número de Schmidt (ad.)

kℎ – Número de Sherwood (ad.)

k~ – Número de Stanton (ad.)

GN� – Temperatura média da interface (K)

GN��� – Temperatura média da envolvente (K)

GN� – Temperatura média do ponto i (K)

GN�� – Temperatura média da superfície i (K)

GN5 – Temperatura média da água (K)

G�J – Temperatura adiabática chama (K)

G���J – Temperatura do condensado (K)

G%� – Temperatura dos gases no ponto i (K)

G��� – Temperatura de saturação (K)

G50� – Temperatura da água à saída do

condensador i (se i=0, não entra) (K)

G5� – Temperatura água à saída de uma das 3

zonas (K)

G� – Temperatura ambiente (K)

c – Coeficiente de convecção global (W/m2K)

6 – Largura da fita (m)

�� – Fracção molar da substância i (ad.)

�� – Comprimento projectado da alheta

(*J/~��) (m)

�d – Parâmetro geométrico do permutador de

calor (m)

f� – Vetor aproximação da solução do sistema

de equações do método de Broyden

f�� – Comprimento da zona de combustão (m)

f����� – Comprimento da zona de convecção

(m)

f����J – Comprimento da zona de

condensação (m)

� – Razão da dobra da fita (*/6) �� – Fração mássica da espécie (ad.)

�1 – Vector de Δ7 do método de Broyden

��� – Largura da zona de combustão (m)

������ – Largura da zona de convecção (m)

�����J – Largura da zona de convecção (m)

xiv

��� – Altura da zona de combustão (m)

������ – Altura da zona de convecção (m)

�����J – Altura da zona de condensação (m)

Símbolos gregos

�� – Absorvidade da superfície i (ad.)

� – Ângulo entre a vertical e o início de

excesso de condensado (°)

� – Razão de convergência (ad.)

�– Rendimento (ad.)

��– Rendimento global das alhetas (ad.)

ΔGeW – Diferença de logarítmica de

temperatura (K)

/� – Espessura da alheta (m)

���� – Emissividade da envolvente (ad.)

�� – Emissividade da superfície i (ad.)

� – Eficiência do permutador de calor (ad.)

� – Tensão superficial (mN/m)

Υ – Parâmetro da alheta molhada (ad.)

Λ – Temperatura adimensional

� – Difusividade térmica (m2/s)

�( – Difusividade térmica ar húmido (m2/s)

λ – Coeficiente de excesso de ar (ad.)

� – Viscosidade dinâmica normalizada (ad.)

�� – Viscosidade dinâmica do ar seco (kg/sm)

�� – Viscosidade dinâmica da água (kg/sm)

�( – Viscosidade dinâmica do ar húmido

(kg/sm)

�� - Viscosidade dinâmica do vapor de água

(kg/sm)

� – Viscosidade cinemática (m2/s)

Φ – Ângulo da ondulação das alhetas (°)

� – Propriedade calorífica especial (ad.)

�� – Massa específica água (kg/m3)

�( – Massa específica ar húmido (kg/m3)

�� – Massa específica vapor de água (Kg/m3)

�� – Reflectividade da superfície i (ad.)

σ – Constante de Stefan-Boltzmann

(5,67040×10-8Js-1m-2K-4) ¢£ – Razão entre a área frontal ao permutador

com a área mínima de escoamento livre

(ambas perpendiculares ao escoamento) (ad.)

¤� – Transmissividade da superfície i (ad.)

� – Coeficiente do método ε-NTU-θ (ad.)

Θ – Factor de sucção (ad.)

¦ – Relação com temperatura adimensional

(1 − ¤) ¨ – Parâmetro das alhetas (ad.)

© – Parâmetro da correlação de Gnielinski

(ad.)

1

1 Introdução

1.1 Objetivo

Este trabalho tem como objetivo a modelação de uma caldeira de condensação em MS ExcelTM

com recurso à linguagem de programação VBATM (Visual Basic for Applications) de modo a avaliar o

desempenho da instalação de uma caldeira de condensação, quando existe a variação de um ou

mais parâmetros da caldeira, nomeadamente o comportamento da caldeira de condensação para

diferentes valores de humidade relativa do ar ambiente e o impacto que esta variação gera no

rendimento e funcionamento da caldeira. Para tal é necessário proceder à caracterização do sistema

envolvido e recorrer a algumas simplificações e/ou premissas de modo a tornar o modelo matemático

o mais simples possível.

Deste modo a modelação da caldeira de condensação irá basear-se nas trocas de calor por

zonas como a zona da câmara de combustão (queimador incluído) e permutadores de calor (sem e

com condensação). Com base nestes balanços e usando o método melhorado de Broyden para a

resolução de sistemas de equações não lineares será possível obter as temperaturas ao longo da

caldeira de condensação. A validação destes resultados pode ser feita com recurso aos dados

experimentais já adquiridos pela Bosch.

O trabalho também pretende generalizar o programa criado de maneira a comportar-se de

maneira aceitável para vários modelos e configurações de caldeiras de condensação ou não, de

modo a dar utilidade à mesma no desenvolvimento ou aperfeiçoamento de outros modelos. Esta

generalização irá implicar a escolha de variadas correlações bem como a adoção de uma geometria

standard e simplificada que pode ser influenciada por um fator corretivo sempre que necessário uma

vez que as geometrias das caldeiras de condensação, nomeadamente nos permutadores podem

variar significativamente. Disto isto o programa será o mais abrangente possível apesar de existirem

situações em que modificações mais profundas no programa possam ocorrer.

Em suma, é esperado que o modelo criado recrie o melhor possível, as condições reais de

funcionamento de uma caldeira de condensação.

1.2 Trabalhos Anteriores

O desenvolvimento de modelos de balanços de energia com transferência de calor é efetuado em

geral em paralelo com o desenvolvimento de equipamentos e não são muito divulgados. Makaire et

al. (2010) refere um modelo para caldeiras de condensação previamente proposto por Hanby (2007)

que utiliza um modelo geral de caldeira com modelo de tubo em espiral para a zona de condensação.

Makaire et al. (2010) propõe então um modelo baseado no modelo de Hanby (2007) mas com

uma diferença. Na zona de condensação o tubo em espiral é dividido em 5 partes e será utilizado o

modelo de Morisot (2000) em cada uma. Foram utilizadas duas caldeiras das quais as suas

propriedades, nomeadamente coeficientes globais de convecção dos permutadores e caudais foram

determinados a priori, pelo que a incerteza máxima obtida em relação ao rendimento calculado foi de

2

2,5% ao variar a temperatura da água de entrada, potência e excesso de ar. Anthony, S. (1986)

também criou um modelo para caldeiras de condensação que necessita de dados experimentais

sobre os coeficientes globais de convecção mas no entanto propõe igualmente correlações baseadas

no fator de Colburn para os permutadores de calor.

Miranda et al. (2003) apresentam um modelo simples para caldeiras sem condensação, incluindo

um modelo para a combustão com escoamento tampão que lhes permitiu calcular a altura da zona de

combustão e ainda a produção de NOx e de CO. Para a transferência de calor consideraram a zona

da combustão com equações de radiosidade e no permutador ambos os mecanismos de convecção e

radiação. A contribuição da convecção nesse permutador é preponderante e aumenta com o caudal

de água a aquecer.

Já ao nível da análise do efeito da quantidade de vapor de água presente no ar à entrada Kuck

(1995) conclui que, à pressão atmosférica, o ponto de orvalho não excede os 60°C pelo que o uso de

ar pré-aquecido com um aumento de vapor de água à entrada aumenta de maneira significativa a

eficiência da caldeira.

1.3 Enquadramento

O elevado rendimento da caldeira de condensação (> 90%) e a subida alarmante dos preços do

gás fazem esperar que estas se tornem cada vez mais comuns nas habitações residenciais uma vez

que, os governos estão a pressionar o consumidor para uma redução de consumos energéticos

através de incentivos fiscais para a troca de equipamentos de “baixo” rendimento por outros de

elevado rendimento. Até é possível encontrar países onde é obrigatório possuir uma caldeira de

condensação ao invés de um esquentador normal segundo o relatório The Warren Report (2007). Os

incentivos fazem inferir que uma caldeira de condensação é mais dispendiosa que uma caldeira

normal (i.e. atmosférica).

A certificação das caldeiras requer a verificação em ensaios do seu elevado rendimento sendo

este valor determinante para poderem ser alvo de benefícios na sua aquisição. As normas para

condução de ensaios (10 CFR Part 430, Point 6.2.2.2 – Test Procedure Rule – U.S. D.O.E.) apesar

de extensas, não preveem qualquer controle da humidade relativa aquando da realização dos

ensaios experimentais, nem métodos para corrigir o efeito deste parâmetro para condições de

referência. Por outro lado a temperatura de entrada da água e do ar são fixas pelo que não é

necessário introduzir correções ao valor aferido para a eficiência da caldeira de condensação. A

norma também prevê que se reporte a percentagem de dióxido de carbono em base seca e não

húmida. Pretende-se então que o modelo a desenvolver neste trabalho indique o efeito da humidade

relativa no rendimento. O rendimento tradicional de uma caldeira de condensação pode ultrapassar

os 100%, uma vez que na Europa é utilizado como referência o poder calorífico inferior. Na legislação

Americana é utilizado o poder calorífico superior para o cálculo do rendimento dos equipamentos pelo

que este não pode ultrapassar os 100%. As caldeiras de condensação apresentam um elevado

rendimento (> 90%) mas, requerem mais manutenção que outros equipamentos equivalentes de

aquecimento de águas sanitárias (AQS). Estes elevados rendimentos são atingidos a potências

3

reduzidas a moderadas, para uma dada caldeira, uma vez que a estas potências existe condensação

pois a temperatura dos gases de combustão é inferior ao ponto de orvalho. Deste modo as caldeiras

de condensação são mais eficientes em regimes de carga parcial do que em regimes de carga total.

As caldeiras de condensação têm sido usadas em conjunto para AQS e para aquecimento central

dado o elevado rendimento das mesmas. Desta forma, é de notar que dada a complexidade da

maioria das caldeiras, se torna necessário recorrer a funcionários especializados durante a instalação

de modo não prejudicar o rendimento efetivo com uma instalação deficiente (Inman, 2005).

Estes equipamentos também têm outras particularidades que incrementam o preço,

nomeadamente a baixa temperatura de saída dos gases de combustão da caldeira que faz necessitar

de ventilação mecânica e pelo facto de o condensado que saí por uma tubagem ser ácido, o que faz

com que o(s) condensador(es) tenham que ter materiais resistentes à corrosão (e.g. alumínio).

Também o simples facto de existir uma tubagem para escoamento do condensado faz aumentar o

custo de aquisição de uma caldeira de condensação embora de forma insignificante.

1.4 Estrutura da dissertação

Esta introdução pretende explicar de maneira sucinta o conteúdo e âmbito da dissertação pelo

que está dividia em 5 capítulos incluindo a introdução.

O segundo capítulo é uma revisão bibliográfica de correlações para convecção para o permutador

de calor primário e zona de condensação. O terceiro capítulo é a apresentação do modelo tal como

foi idealizado e inclui a caracterização do método numérico usado e caracterização das propriedades

dos fluídos.

O quarto capítulo trata os resultados incluindo a discussão de resultados. O quinto capítulo

apresenta as conclusões e recomendações para futuros trabalhos.

4

2 Revisão Bibliográfica de Correlações para Convecção

A revisão bibliográfica apresentada tem como objetivo identificar correlações de transmissão de

calor específicas para as configurações observadas na caldeira de condensação. Deste modo na

secção 2.1 apresentam-se correlações e informações para o escoamento no interior de tubos onde é

prática usual incluir elementos para o aumento da transferência de calor. A secção 2.2 apresenta

correlações para o escoamento no exterior de tubos alhetados. Na secção 2.3 apresentam-se as

correlações para a transferência de calor em permutador com condensação de vapor de água dos

produtos de combustão.

2.1 Interior dos tubos

O coeficiente de convecção para o interior dos tubos pode ser calculado com recurso a uma

correlação existente na literatura que se aproxime do desejado. Com o objetivo de aumentar o

coeficiente de transmissão de calor por convecção no interior dos tubos usam-se superfícies que

promovam turbulência designadas por turbulators que apresentam várias configurações de modo que

é necessário escolher correlações da literatura apropriadas. Estes turbulators podem criar ou

aumentar a turbulência do escoamento e em alguns casos este efeito pode ser previsto com recurso

a correlações de escoamento turbulento, sem ser necessário correlações complexas. Também é

preciso ter em conta que o tubo pode não ser circular, existindo correlações específicas para estes

casos como discutido em Khan et al. (2004). Para o escoamento no interior de tubos existem as

correlações propostas por Sieder e Tate (1936) e Shah et al. (1987) para regime laminar que

satisfazem e preveem relativamente bem o teste de tubos elípticos usando o diâmetro hidráulico em

regime laminar e de transição (Khan et al., 2004). Para a zona de transição a correlação de Gnielinski

(1976) para tubos circulares e em regime turbulento estima o número de Nusselt por excesso. Meng

et al. (2005) verificam que a correlação de Gnielinski (1976) satisfaz o cálculo do número de Nusselt

com o escoamento de transição e turbulento em tubos elípticos com eixos alternados com um erro de

6% em relação aos valores experimentais e que a correlação de Sieder e Tate (1936) também o

satisfaz com um erro de 10%. Desta maneira, para os tubos elípticos existe a possibilidade de usar

correlações de tubos circulares desde que usando o diâmetro hidráulico garantindo assim erros

dentro do admissível segundo Khan et al. (2004) e Meng et al. (2005). É apresentada, de seguida, a

correlação de Gnielinski (1995) para o número de Nusselt médio para escoamento de transição e

turbulento. (Huber et al., 2010)

bgd = =: 8⁄ ) 3Jy*l1 + 12,7ª: 8⁄ = *l=B «⁄ ) − 1) ¬1 + P#���\

B «⁄ ­ (1)

: = ®1,8 − log_0® 3Jy± − 1,5±²B (2)

Para a zona de transição Gnielinski (1995) introduziu uma interpolação entre os números de

Nusselt calculados para os regimes laminar e turbulento para os números de Reynolds limite

indicados no número de Nusselt.

bg( = =1 − ©)bge,B«00 + ©bgd,_0000 (3)

5

© = 3Jy − 230010³ − 2300 , 0 ≤ © ≤ 1 (4)

bge,B«00 = µ49,371 + ®bge,_,B«00 − 0,7±« + ®bge,B,B«00±«¶_ «⁄ (5)

bge,_,B«00 = 1,615 P2300*l #���\_ «⁄

(6)

bge,B,B«00 = P 21 + 22*l\ P2300*l

#���\_ B⁄

(7)

Tabela 2.1 - Limite de validade das correlações para escoamentos internos

Correlação Limites de validade

Gnielinski (1995) 2300 ≤ Re ≤106, 0,1 ≤ Pr ≤1000

As seguintes correlações são já para escoamento turbulento na presença de turbulators, o que

faz naturalmente aumentar o coeficiente de convecção no interior dos tubos, pois o diâmetro

hidráulico diminui. A transição para o regime turbulento também ocorre para um número de Reynolds

superior ao de tubos circulares, como sugerido por Incropera e De Witt (2006).

Para a modelação selecionaram-se correlações de configurações de turbulators semelhantes à

do modelo em causa mas também para outras configurações diferentes, de modo a permitir uma

posterior generalização caso seja necessário.

Estas correlações são também válidas para o escoamento interior na zona de condensação e a

tabela seguinte indica os limites de validade das correlações para os erros que os autores referem.

Quando a aplicação é efetuada fora da gama referida é gerado um aviso não sendo conhecido a

margem de erro introduzida pelas correlações. As correlações apresentadas de seguida dizem

respeito às geometrias e referências apresentadas na tabela 2.2.

Tabela 2.2 - Correlações com o uso de turbulators

Designação, Referência e Correlação Ilustração da geometria

Delta-Wing

Eiamsa-ard e Promvonge (2011)

bg = 0,112 30,·«_*l0,³3�²0,B¸¹350,Bº·

Caso de escoamento na direção das asas

(Esquerda para direita)

Limites de validade

4000 ≤ 3 ≤ 20000, 0,75 ≤ 3� ≤ 1,25,

0,5 ≤ 35 ≤ 0,83,*l = 0,7

6

Twisted tape and wire coil

Promvonge (2008)

bg = 4,47 30,¹ *l0,³ � ²0,«ºB�²0,«º Limites de validade

3000 ≤ 3 ≤ 18000, d/D½0,1, 4 ≤ � ≤ 8, 4 ≤ � ≤ 8, �/) ½32,*l = 0,7

V-Nozzle

Promvonge e Eiamsa-ard (2007)

bg = 0,37 30,¸«¹*l0,³* ²0,BB Limites de validade 8000 ≤ 3 ≤ 18000, 2 ≤ * ≤ 7, 4 ≤ � ≤ 8, �/) ½ 32,*l = 0,7

Wire Coil

Eiamsa-ard, Kongkaitpaiboon e

Promvonge (2011)

bg = 0,0545 30,·¸·� 0,0B Limites de validade 4500 ≤ 3 ≤ 20000, d/D½0,1, 4 ≤ � ≤ 8

2.2 Exterior dos tubos

O coeficiente de convecção para o exterior dos tubos requer que se usem correlações que

tentem aproximar ao máximo a configuração corrente na indústria, de maneira a obter o menor erro

possível uma vez que cada fabricante e/ou modelo tem o seu conjunto especifico de alhetas

desenhadas de tal modo que, recorrendo unicamente à literatura é apenas possível apresentar uma

aproximação dos resultados que se devem obter experimentalmente para o coeficiente de convecção.

Deste modo, os desenhos da indústria atualmente querem maximizar a transferências de calor

por isso é possível ter tubos alhetados com alhetas não regulares bem como tubos não regulares (por

exemplo tubos elípticos). Deste modo a escolha mais comum para obter uma correlação aproximada

da realidade é pensando nas alhetas como alhetas onduladas do tipo herringbone (wavy fins). Esta

aproximação pode no entanto ser conservadora para alguns modelos. Segundo Goldstein e Sparrow

(1976) o uso de alhetas onduladas do tipo herringbone (não planas) aumentou em 45% a

7

transferência de calor local, comparado com alhetas planas, para um valor de Reynolds baseado em

dc de 1000.

Quanto ao facto dos tubos poderem ser elípticos é necessário ter em conta que segundo Brauer

(1964) o uso de um tubo elíptico permite diminuir a perda de pressão à passagem do escoamento

perpendicular aos tubos, bem como um aumento da transferência de calor que só se torna

significativa quando a razão de eixos da elipse é maior que 0,30, algo que foi verificado por

Goldschmidt et al. (2002).

Deste modo para filas de tubos em linha, Nishiyama et al. (1988) verificou que quanto mais juntos

estiverem os grupos de tubos em linha e quanto menor o angulo de ataque menos arrasto irão

provocar melhorando a transferência de calor e tornando o permutador mais compacto. Então,

segundo Khan et al. (2004) ter-se-ão aumentos na transmissão de calor de 40% do lado exterior

quando comparados com tubos circulares alinhados. Note-se que no entanto este aumento foi

realizado para regime laminar ou de transição.

Em regime turbulento para tubos desalinhados e com recurso à figura 2.1 de Ibrahim et al. (2009)

tem-se que o ângulo de ataque melhora a transmissão de calor em que o máximo está nos 90° para

diferentes razões de eixos. No entanto isto cria uma queda de pressão intensa o que leva a que piore

a performance térmica sendo esta performance maior a 0°. A performance térmica foi definida como

sendo a razão entre o número de Stanton e o coeficiente de atrito normalizados a um caudal mássico

fixo hh St

St0i h f

f0i¾ i. Logo existe uma melhoria na transmissão de calor para ângulos a rondar os 90° para

tubos desalinhados mas no entanto a perda de pressão torna-se intensa.

Figura 2.1 - Exemplo da influência dos tubos elípticos a diferentes ângulos de ataque para uma razão

de eixos de 0,5 (esq.) e 0,25 (dir.) - Ibrahim et al. (2009)

8

Tao et al. (2007) estudou filas de tubos elípticos alhetados desalinhados com alhetas onduladas

do tipo herringbone com ângulos de ataque de valor zero e verificou um aumento de 30% na

transmissão de calor e um aumento de 10% no fator de atrito quando comparados com tubos

circulares. Este valor de aumento da transmissão de calor está em linha com os aumentos verificados

por Khan et al. (2004) para filas de tubos elípticos alinhados. Lin et al. (2011) confirma que tubos

elípticos comportam-se melhor ao atingirem coeficientes de convecção maiores.

Apresentam-se uma correlação da literatura para conjuntos de tubos desfasados com alhetas

onduladas do tipo herringbone. Existem diversos trabalhos para esta configuração geométrica

conforme revisão apresentada por Kim et al. (2008). Adicionalmente Youn e Kim (2007) aumentam a

gama de parâmetros geométricos.

Figura 2.2 – Tubos alhetados com alhetas onduladas do tipo herringbone com indicação de alguns

parâmetros geométricos – Kim et al. (2008)

Desta maneira a correlação para o coeficiente de convecção exterior deve ter em conta a figura

2.2. A correlação é para um número de uma a quatro filas de tubos desfasados, sem apresentar

correções para a variação do número de filas, que apresenta um efeito importante para baixos

números de Reynolds segundo Youn e Kim (2007). Na tabela 2.3 é apresentado os limites de

validade dos parâmetros da correlação usada. Outras duas correlações estão disponíveis no anexo 2

a título de apoio para validação de resultados. No anexo 2 são também apresentados os valores

escolhidos para a caldeira em questão, demonstrando assim que “aproximações” foram feitas.

Φ

Φ

9

�W�« = ¿1,04 + 0,0067bn�l� 3Kw > 1000, b = 1,21,69 − 0,254bn�l� 3Kw ≤ 1000, b = 1,2 (8)

�« = 0,170 3Jw²0,««Á P*d*e\0,_Á« P8�#�\²_,º« P��*K\²B,·¹ P*K8� \²_,º³ (9)

Tabela 2.3 - Limite de validade da correlação de Kim et al. (2008)

3Jw *d *e⁄ * 8�/#� ��/*K *K/8� NL Φ (°) 500-6000 1,15-1,33 0,13-0,44 3,01-4,82 0,30-1,03 1-4 11,9-19 * - Indicativo de que esta gama não é cumprida para o modelo em causa.

De acordo com o modelo da Bosch os tubos na zona radiativa não são circulares mas sim

elípticos. É então necessário aplicar uma de duas possíveis correções. Usar o diâmetro hidráulico

e/ou então introduzir um fator corretivo na correlação. Optou-se por utilizar apenas o diâmetro

hidráulico. Deste modo as seguintes equações indicam o diâmetro hidráulico do tubo, (exterior) bem

como o fluxo mássico e o número de Reynolds em função do diâmetro do colar aqui também

indicado.

#,� = 4 Ä��ÄÅ2(�B + �B) − (� − �)B2

(10)

#� = #,� + 2/� (11)

3Jw = Aa %����#��%����¢£ ������� (12)

Para o cálculo de ¢£, i.e. a razão entre a área de passagem e a área frontal do permutador de

calor, foi considerada a expressão apresentada de seguida em que a mesma terá em conta a

geometria padrão de uma fila de tubos alhetados ondulados do tipo herringbone. Desta maneira a

espessura da alheta pode variar em relação ao valor de espessura da alheta real, uma vez que é

preciso acomodar o facto de o(s) modelo(s) a testar poderem não ter alhetas onduladas do tipo

herringbone mas semelhantes.

¢£ = (f�)����� − #�f�����bd − b� P� − bd#��� \ P ��&u!Φ\ /�(f�)�����

(13)

De forma a calcular o rendimento global é necessário conhecer o rendimento da alheta bem como

as suas áreas. Deste modo, o rendimento da alheta é dado então por Schmidt (1949):

�� = tanh (Al¨)Al¨ (14)

Onde A é calculado a partir da seguinte expressão:

A = Ê 2ℎ���/� (15)

10

Considerando alheta hexagonal (tubos desalinhados, ou NT>1) é possível saber o rendimento da

alheta ondulada do tipo herringbone introduzindo uma dimensão característica para o comprimento

da alheta em relação ao raio do tubo na forma:

¨ = hl�,�zl − 1i µ1 + 0,35 ln hl�,�zl i¶ (16)

Nesta expressão rf,eq é o raio do círculo equivalente em torno de cada tubo que tem a mesma

área que a alheta completa. Esta dimensão para tubos desalinhados é calculada a partir de:

l�,�zl = 1,27 *dl ˪(*d 2⁄ )B + *eB*d − 0,3Ì0,¹ (17)

Enquanto para tubos alinhados ou para uma única fila de tubos é dada por:

l�,�zl = 1,28 *dl ˪(*d 2⁄ )B + *eB*d − 0,2Ì0,¹ (18)

Deste modo é possível obter o rendimento global que pode ser calculado da seguinte forma.

Como a geometria das alhetas é complexa e normalmente é dada a área total do permutador de calor

primário e o número de alhetas no permutador criou-se PL* que, não é mais que o valor que PL teria,

i.e. PL equivalente, se a área total fosse a indicada. Note que para as correlações utilizou-se à

mesmaPL.

� = 2be Í*e∗� − Ä#�B4 bdÎ (19)

� = b� � + � = b� � + ®Ä#�fbdbe − bebd/�b�#�± (20)

�� = 1 − � � (1 − ��) (21)

O coeficiente global de convecção que pode ser descrito pela equação seguinte para tubos

alhetados.

c = P 1 �ℎ� + ln (#�/#�)2Ä���� + 1�� �ℎ�\²_ (22)

Os coeficientes de convecção para o interior e para o exterior podem ser calculados recorrendo

às correlações propostas na secção 2 tal como o rendimento global pelo fato do permutador ser

constituído por tubos alhetados.

É possível ainda verificar as vantagens do uso de alhetas onduladas ao comparar a correlação

proposta por Kim et al. (2008) com Xie et al. (2009). Apresenta-se então a correlação de Xie et al.

(2009) para alheta planas em que a tabela 2.4 contém os limites de validade para os vários

parâmetros.

bg = 1,565 3Jw0,«³_³ Pbe *�#�\²0,_¸¹ P*d*e \0,0¹¹º

11

Tabela 2.4 - Limites de validade da correlação de Xie et al. (2009)

3Jw #� *� *d *e NL∗ 1000-6000 16-20 mm 2-4 mm 38-46 mm 32-36 2-6

* - Indicativo de que esta gama não é limitada pelo autor

A figura 2.3 revela então o aumento significativo do número de Nusselt para as alhetas onduladas

se for considerado um número de Prandtl constante de 0,7 e duas filas de tubos (N=2). A diferença é

mais dramática ao considerar-se tubos alinhados, como será o caso, uma vez que Kim e Kim (2005)

concluíram que existe uma diminuição de 10% relativamente ao número de Nusselt obtido para tubos

desalinhados.

Figura 2.3 - Evolução do número de Nusselt em função de Reynolds para Kim et al. (2008) (azul),

com parâmetros geométricos correspondentes ao permutador analisado e Xie (2009) (vermelho)

2.3 Convecção no condensador

As correlações a serem usadas no condensador devem levar em conta a influência da

condensação do vapor de água que se encontra misturado com outros gases de produtos de

combustão não condensáveis. São apresentados então dois modelos, um de condensação de fluido

puro e outro para a condensação de um gás misturado com outros gases não condensáveis, que é a

situação do vapor de água nos produtos de combustão. Apesar de o modelo de condensação de

fluido puro não considerar gases condensáveis, o mesmo é apresentado de modo a avaliar a

discrepância entre os modelos e se é aceitável usá-lo no modelo.

2.3.1 Modelo 1 – Condensação de vapor

O primeiro modelo considerado baseia-se na teoria de formação de filme de condensado na

presença de vapor. As correlações a utilizar dependem da configuração geométrica do permutador e

0

10

20

30

40

50

0 1000 2000 3000 4000 5000 6000

Nu

Re

12

do escoamento do vapor. Considera-se o caso do escoamento em torno de um tubo com alhetas

anulares de perfil retangular. Neste caso o condensado acumula-se entre as alhetas, cobrindo o tubo

(e alheta) ao longo de um ângulo β ilustrado na figura 2.4.

Figura 2.4 - Exemplo de excesso de condensado criado em um tubo alhetado – Webb et al. (1983)

A figura 2.4 representa o condensado a acumular-se ao longo do tubo alhetado e deste

condensado é possível tirar a primeira propriedade da correlação, o ângulo de excesso de

condensado desenvolvido por Webb et al. (1983).

� = cos²_ P1 − 4�#����8�\ (23)

Na equação acima não é considerado o efeito da velocidade do vapor. O efeito da velocidade do

vapor para tubos alhetados no coeficiente de convecção é limitado de acordo com Gogonin e

Dorokhov (1981) mas Bella et al. (1993) concluiu que o efeito começa a ser significativo para um

Reynolds (Redo) de vapor com caudal de vapor superior a 10. Também verificou-se que para

velocidades de vapor superiores a 30 m/s o aumento da transferência de calor aumentava em 50%

em relação ao vapor estagnado.

Ainda segundo Fitzgerald et al. (2012) a velocidade de vapor tem uma relação com β como já foi

dito mas que, no entanto, não é semelhante para todos os fluidos mas que no caso da água e do

etilenoglicol existe uma coerência quando ao sentido da evolução uma vez que, quando a velocidade

de aproximação aumenta β diminui.

De acordo com o modelo de Webb et al. (1985) tem-se que o coeficiente de convecção exterior

distingue a parte inundada da parte em que está presente o filme de condensado. Isto deve-se ao

fato de o coeficiente da zona inundada corresponder a apenas 2% do total do coeficiente médio

calculado. Isto leva a que o efeito da zona inundada para efeitos do filme de condensado seja

desprezado.

De acordo com Webb et al. (1983) obtém-se a expressão de hf. A expressão para hh mantém-se

inalterada e consiste no coeficiente de convecção para um tubo horizontal presente em Incropera e

De Witt (2006). Note que Ts do condensador será igual à temperatura água à saída do primeiro

condensador, Tw01.

13

ℎ� = 0,943 ¬ ����«ℎ�%q��=G��� − GN�$­_³ ¬ 8�=#� − #�$B Í 18� + 1/�έ

_³ (24)

ℎ, = 0,729 ¬���=�� − ��$��«ℎ�%q��=G��� − GN�$#� ­_³ (25)

Desta forma, e para o caso de não existir condensação, a seguinte correlação de Briggs e Young

(1963) satisfaz os requisitos mas, no entanto, é necessário saber o diâmetro hidráulico do

condensador que para os autores da correlação considera-se como o diâmetro exterior do tubo

alhetado e σA.

� = 0,134 3Kx²0,«_Á Í8���Î0,B Í8�/�Î

0,__ (26)

3K, = Aa %����#��%����¢£ ������� (27)

¢£ = =�#0$����J − #������J − /�=#0 − #�$b�=�#0$����J (28)

2.3.2 Modelo 2 - Condensação de vapor em mistura

O segundo modelo já prevê o facto de existir gases não condensáveis como ilustra a figura 2.5. O

modelo usado para prever esta existência é baseado em Shi et. al. (2011), Herranz et al. (2000) e

Webb et al. (1983,1985). O modelo utiliza mecanismos de difusão para caracterizar a camada de

gases e a abordagem de Nusselt para a película de condensado que serão separados por zona não

inundada e zona inundada. Devido à existência de gases não condensáveis existe então uma

resistência adicional correspondente aos gases não condensáveis pelo que se introduz uma

temperatura G� que corresponde à temperatura da interface entre o filme de condensado e a camada

dos gases não condensáveis. Também é ilustrado na figura 2.5 o fato de considerar-se dois fluxos de

calor, um latente (i.e. condensação) e outro sensível (i.e. convecção). Como nota, repare-se que a

fração mássica de vapor de água junto à interface será sempre menor que fora das camadas em

discussão devido ao fato de, ao condensar, a sua fração diminui junto ao filme de condensado.

14

Figura 2.5 - Exemplo de uma camada limite de difusão na presença de gases não condensáveis (NC)

– Shi et al. (2011)

A taxa de transferência de calor por parte dos gases para a camada de filme é quantificada para

a área M�� que é a área pertencente à zona não inundada na forma:

p = M��ℎ�%®GN% − GN�± (29)

Em que a área não inundada é dada pela expressão:

M�� = P1 − �Ä\ � (30)

Deste modo, e de acordo com Shi et al. (2011) e Herranz et. al. (2000), o coeficiente de

convecção dos gases é dado por uma componente de convecção (calor sensível) e outra de

condensação (calor latente).

ℎ�% = ℎ����Θ���� + ℎ���JΘ���J = Íbg�%#� ÎΘ���� + Pkℎ)�#� \ Í �ℎ�%qGN% − GN�ÎΘ���J (31)

Em que é ℎ�%q a entalpia de vaporização modificada dada por:

ℎ�%q = ℎ�% + &n�®G�,� − 25℃± − &n�®G�,� − 25℃± (32)

Onde, G�,� é a temperatura do vapor de água à entrada e G�,� é a temperatura do condensado à

saída. Já � é a driving force do mecanismo de difusão e é dada por:

� = ��,F − ��,�1 − ��,� (33)

O número de Sherwood este pode ser retirado usando as mesmas equações para o número de

Nusselt desde que o número de Prandtl seja substituído pelo número de Schmidt segundo a analogia

de transferência de calor-transferência de massa. Note-se que o número de Reynolds é retirado

através da equação 27.

15

bg = 0,134 3Kx0,¸º_ Í8���Î0,B Í8�/�Î

0,__ *l_/« (34)

kℎ = 0,134 3Kx0,¸º_ Í8���Î0,B Í8�/�Î

0,__ k&_/« (35)

Juntamente, Θ é considerado o efeito de sucção que é causado pela diminuição da camada limite

do filme de condensador, devido à absorção de vapor de água pela camada que resulta em

gradientes crescentes. Desta forma o fator de sucção para a condensação Θcond é dado por Bird et.

al. (1960) e para a convecção Θconv por Herranz (1996). Note-se que Xv,avg e Xnc,avg são médias

logarítmicas da concentração do vapor de água e dos gases não condensáveis.

Θ���J = lnЮ1 − ��,F± ®1 − ��,�±¾ Ñ®��,� − ��,F± ®1 − ��,�±¾ (36)

Θ���� = &1 − 3²� (37)

& = Í�%���� %����&�,%����)�% ÎÍ���,� − ���,F���,��% Î (38)

Utilizando ainda a figura 2.5 é possível retirar a transferência de calor no filme de condensado:

p = M��ℎ��(=GN� − GN�$ (39)

Em que o coeficiente de convecção do filme de condensado é dado por:

ℎ��( = ℎ, � � + ��ℎ� � � (40)

Para o qual ℎ, e ℎ� apresentam a diferença de temperatura da interface com a temperatura da

superfície (ao invés da saturação com a superfície). O rendimento da alheta pode ser retirado

utilizando as expressões presentes para o cálculo do coeficiente global do permutador primário.

É necessário então encontrar a temperatura da interface TIÒ . Mais uma vez, com recurso ao

método de newton simplificado, apresentado por Kelly (2003), é possível retirar TIÒ a partir do balanço

de energia à camada limite total. Repare que ambos os coeficientes de convecção dependem de TIÒ .

ℎ��(=GN� − GN�$ = ℎ�%®GN% − GN�± (41)

Com esta formulação é possível determinar o coeficiente de convecção médio, quando na

presença de gases não condensáveis durante a condensação, para a parte não inundada. É de

salientar que para o modelo 1 apenas existe ℎ��(.

ℎM�� = Í 1ℎ��( + 1ℎ�%β_

(42)

Quer para o modelo 1 quer para o modelo 2, o coeficiente global de convecção pode ser expresso, separadamente, pela zona não inundada e a zona inundada.

c = = c$M�� + = c$�� = � ÓP1 − �Ä\cM�� + P�Ä\c��Ô (43)

16

Em que cada uma das zonas é a soma de resistências térmicas em série em que estão por metro quadrado de área não inundada e inundada, respetivamente.

cM�� = Õ 1ℎM�� + P � �\ 1ℎ� + � Öln P#�#��\2Ä���� + ln P#��#� \2Ä����×ز_

(44)

c�� = Õ 1ℎF + P � �\ 1ℎ� + � Öln P#�#��\2Ä���� + ln P#��#� \2Ä���� ×ز_

(45)

Note-se que quer para o modelo 1 quer para o modelo 2, o termo em que ℎF aparece é ignorado pelas razões expostas no modelo 1.

A figura 2.7 mostra que para convecção forçada (quanto maior o caudal melhor) e pequenas

concentrações mássicas de ar é possível minimizar este efeito. Deste modo, a figura também mostra

que grandes concentrações de ar fazem diminuir significativamente o calor transferido para a parede.

Figura 2.6 - Ilustração da variação do coeficiente de convecção em função do caudal e fracção

mássica de ar a diferentes temperaturas do gás – Collier et al. (1994)

17

3 Modelação de uma caldeira de condensação

Este capítulo apresenta a formulação do modelo que foi utilizado para descrever o funcionamento

da caldeira com condensação do vapor de água dos produtos de combustão. O modelo é constituído

por diversos módulos correspondentes a zonas da caldeira. A secção 3.1 apresenta uma descrição

geral do modelo implementado sendo detalhado enquanto a secção 3.2 apresenta uma análise à

zona de combustão. Já a secção 3.3 descreve os balanços de energia à câmara de combustão

enquanto as secções 3.4 e 3.5 contêm os balanços ao permutador primário e ao(s) condensador(es).

A secção 3.6 apresenta o balanço das perdas. Por fim a secção 3.7 apresenta correlações para

caracterizar as propriedades dos fluidos que compõem o modelo e a secção 3.8 indica a metodologia

de resolução do modelo matemático com recurso ao método melhorado de Broyden.

3.1 Modelo

Figura 3.1 - Modelo completo com ilustração esquemática das trocas de calor e dos fluxos de

entalpia

Para a modelação da caldeira de condensação, esta foi conceptualmente dividida em quatro

zonas para o cálculo dos balanços de energia como ilustrado na figura 3.1 que também indica os

balanços idealizados para o modelo utilizando uma caldeira semelhante à ilustrada na figura 3.2 que,

indica a localização de cada um dos componentes bem como dos dispositivos auxiliares de ventilação

e sensores. Estas áreas são a zona da câmara de combustão (radiação), a zona dos tubos alhetados

(convecção) e a zona de condensação. Assim sendo os pontos que se seguem dizem respeito a cada

uma destas zonas separadamente mas, que garantem que as fronteiras possam interagir.

18

Figura 3.2 - Exemplo de caldeira de condensação - Cortesia Vulcano

A caldeira de condensação foi então modelada da seguinte maneira: a zona da câmara de

combustão é considerada com trocas de calor por radiação entre as superfícies que a delimitam e

que são a zona da chama na base, as paredes laterais que transferem calor para a água e a

superfície superior que corresponde à saída dos gases de combustão e entrada no permutador. A

base e topo são considerados como corpos negros e a emissividade das paredes são determinadas a

partir de resultados experimentais do calor transferido. Esta emissividade é também importante pois

influencia a temperatura das paredes e assim as perdas de calor para o exterior. Estas perdas são

calculadas com base em correlações para convecção natural no exterior da zona de combustão e nas

trocas por radiação com a superfície do invólucro da caldeira.

A zona dos tubos alhetados é considerada como tubos alhetados com alhetas onduladas do tipo

herringbone pelo que adaptações terão que ser implementadas para obter valores equivalentes do

coeficiente global de convecção uma vez que as alhetas não terão uma geometria simples. Já a zona

de condensação é modelada de acordo com uma geometria simples, ignorando geometrias

complicadas nesta zona. Também não é considerada transmissão de calor por radiação nesta zona

uma vez que é praticamente desprezável.

3.2 Análise da Combustão

Para a combustão no queimador é necessário saber a composição do combustível gasoso, o

coeficiente de excesso de ar e a humidade absoluta. Sabendo a equação química é possível então

calcular a composição química dos produtos de combustão em termos de massa ou moles sem

considerar dissociação. A equação química considerada é a seguinte:

®jÙÚÛ�³ + jÙUÚU�BB + jÙUÚÛ�B³ + jÙUÚÜ�B¸ + jÙÚVÝÚ�«Þ + jÙVÚß�«º + jÙ[àÚU[�_0B_+ jLÙÝU�ÞB + jLWUbB±�M�� + ájÝU��=ÞB + 3,76bB + j�BÞ$→ ®ájÝU��j� + jÚUݱBÞ + jÙÝU�ÞB + jÙÝ�Þ + jÝUÞB + jWUbB

(46)

19

A equação estequiométrica permite definir o número de moles de oxigénio para a combustão

completa do combustível definido acima na forma:

jÝUã� = 3jÙÚÛ + 3jÙUÚU + 4jÙUÚÛ + 5jÙUÚÜ + 2,5jÙÚVÝÚ + 7jÙVÚß + 15,5jÙ[àÚU[ − jLÙÝU (47)

Na equação acima indica-se a presença de CO e O2 que ocorre na realidade mas, como o

coeficiente de excesso de ar usado evita a formação de CO, considera-se apenas combustão

completa e jÙÝ = 0. Deste modo o sistema de equações a resolver é dado por:

äåæåç�: jÙÚÛ + 2jÙUÚU + 2jÙUÚÛ + 2jÙUÚÜ + jÙÚVÝÚ + 3jÙVÚß + 10jÙ[àÚU[ + jÙÝU = jÙÝU + jÙÝ: 4jÙÚÛ + 2jÙUÚU + 4jÙUÚÛ + 6jÙUÚÜ + 4jÙÚVÝÚ + 8jÙVÚß + 21jÙ[àÚU[ = 2jÚUÝÞ: jÙÚVÝÚ + 2jLÙÝU + 2ájÝUã�ÞB = jÚUÝ + 2jÙÝU + jÙÝ + 2jÝUb: jLWU + 3.76ájÝUã� = jWU

(48)

O sistema apesar de extenso é facilmente resolvido ao saber-se a composição do combustível. A

humidade, i.e. vapor de água, é adicionada aos reagentes e produtos na mesma quantidade

sabendo-se a humidade relativa e a temperatura de entrada do ar.

j� = n���@dn − n���@d (4,76jÝU��$ (49)

Consequentemente, desta forma, é possível saber a quantidade de vapor de água à entrada dos

condensadores. Esta situação revela que a temperatura adiabática de chama irá naturalmente baixar

com o aumento da humidade dado que o vapor de água não interage com as restantes espécies

limitando-se a absorver parte do calor sensível fornecido pela reação.

Em suma, o modelo está preparado para lidar com vários combustíveis assim como coeficientes

de excesso de ar maiores ou menores que 1. O excesso de ar pode ser calculado a partir da fração

molar de dióxido de carbono obtida nos ensaios através de uma análise seca (logo o excesso de ar é

da análise seca). A partir da equação 50, em base seca (i.e. excluindo o vapor de água), pode-se

obter a equação para o coeficiente de excesso de ar quando xCO=0, uma vez que o excesso de ar

será sempre maior que 1 e desprezam-se concentrações de monóxido de carbono.

á = ®jÙÝU±�� − fÙÝU®jÙÝU±��4,76fÙÝU®jÝU±�� + 14,76 (50)

Para determinar PCI, que é o simétrico da entalpia da reação/combustão com reagentes e

produtos à temperatura de referência (25°C) temos a seguinte equação:

*�^ = −∆Ò����çã�®G���± = ℎN���%®G���± − ℎN���J®G���± = í j9W9,���J ℎN�,9� ®G���± − í j�î

�,���% ℎN�,�� ®G���± (51)

Em que a entalpia absoluta é dada pela entalpia de formação nas condições de referência. Também

é possível calcular o poder calorifico superior (PCS) da seguinte forma:

*�k = *�^ + ÍjÚUÝj�M��Îℎ�%@dyï (52)

Para determinar a temperatura adiabática utiliza-se uma função com o método de Newton-

Raphson, apresentado por Kelly (2003), que com a composição dos produtos, permite obter a

20

temperatura correspondente a determinado valor da entalpia conhecendo a composição. A entalpia

dos produtos é igual à dos reagentes:

ð��J(G�J$ = X��%(GX��%$ (53)

∆Ò����çã�®G���± + í j9ñ &�,9dòódôõö

W9,���J (G$#G = í j�î

�,���% ñ &�,�(G$#Gdòôdôõö (54)

A temperatura pode, assim, ser facilmente calculada com base nos coeficientes do balanço de

massa. Note-se que a temperatura adiabática de chama depende da temperatura dos reagentes e

também da razão de equivalência ∅=1/λ, o que mostra que G�J pode ser controlada pelo excesso de

ar sendo que é máxima para uma mistura estequiométrica (na ausência de dissociação) segundo

Coelho e Costa (2007). A quantidade de vapor de água é importante pois faz aumentar a proporção

de gases inertes em relação ao oxigénio e combustível, diminuindo a temperatura adiabática para um

valor fixo de excesso de ar. Para manter as temperaturas adiabáticas o coeficiente de excesso de ar

tem de ser diminuído correspondendo a um aumento do caudal de combustível de modo a manter o

caudal total constante.

A temperatura adiabática de combustão (considerando eventualmente a formação de CO) é uma

temperatura que não se observa na caldeira, no entanto é bastante útil para formular os balanços de

energia à superfície do queimador que troca calor com as outras superfícies da câmara de combustão

como se analisa de seguida.

3.3 Balanço de energia à câmara de combustão

Para o balanço à câmara de combustão foi necessário desenvolver um sistema de equações para

a radiosidade das superfícies de modo a determinar os fluxos radiativos existentes na câmara de

combustão. Nos balanços à câmara de combustão considerada a transferência de calor por radiação

e por convecção.

As trocas por radiação entre as superfícies são calculadas considerando a participação dos gases

como gases cinzentos com a emissividade calculada pela correlação de Leckner apresentada e.g. por

Modest (2003). A correlação permite um bom resultado com erros de 5% para vapor de água e 10%

para o dióxido de carbono para temperaturas do gás superiores as 400K quando comparada com a

emissividade total calculada por integração da emissividade espectral. Esta correlação é apresentada

no anexo 1.

As trocas por convecção podem por sua vez ser calculadas utilizando correlações para placa

plana presentes em Incropera e DeWitt (2006) para regime laminar, turbulento ou misto. Para estes

balanços são assumidas várias hipóteses que são apresentadas de seguida. Com base nas figuras

3.1 e 3.3, na zona da câmara de combustão consideram-se três superfícies. A base da câmara que

corresponde ao queimador é considerado uma superfície negra e as paredes laterais são

consideradas como superfícies cinzentas. Quanto à superfície fictícia 2 (topo da câmara de

combustão – permutador de calor), como a radiação que entra no permutador pode sofrer múltiplas

reflexões nas alhetas, considera-se como aproximação que se comporta como uma superfície negra

21

(Azevedo, 2000). Note que 0 e 1 na figura 3.3 refere-se à irradiação incidente do queimador e

paredes laterais. Considerando então a superfície superior da câmara de combustão como corpo

negro as equações de radiosidade para esta zona são:

ø O0 = σTù0³J_ = ε_σTù_³ + �_¤%(F__J_ + F_BJB + F_0J0$ + �_εûσT%«³JB = σTùB³ (55)

Note que os fatores de forma podem ser retirados usando a relação de reciprocidade, a regra do

somatório e a expressão para o cálculo de duas superfícies de iguais dimensões a uma determinada

distância.

A partir das radiosidades calculadas pode-se definir a taxa de transferência de calor para cada

uma das superfícies como:

üzM��(�J�� = �0(;0 − O0$ (56)ü����J�� = �_(;_ − O_$ + �_ℎ��(F(G%« − GN�_$ (57)ü�M����í���B = �B(;B − OB$ (58)

Note que para a equação 57 o termo da troca de calor por convecção tem que entrar. Então, as

expressões das irradiações são dadas por:

;0 = ¤%(80_O_ + 80BOB$ + �%σT%«³ (59);_ = ¤%(8_0O0 + 8__O_ + 8_BOB$ + �%σT%«³ (60)

;B = ¤%(8B_O_ + 8B0O0$ + �%σT%«³ (61)

Figura 3.3 - Zona convectiva - permutador de calor

A conclusão de que a superfície 2 (i.e. superfície fictícia) se comporta como um corpo negro à

temperatura da superfície 4 (i.e. superfície do permutador primário) pode ser analisada a partir da

equação 62 da radiosidade desta última superfície com F44=1, desprezando F

42. A temperatura da

superfície fictícia 2 pode ser assim determinada como sendo aproximadamente a média das

temperaturas de entrada e saída do permutador de calor, Ts2=Tw1+Tw2

2.

22

O³ = �³1 − �³8³³ 2�³ + �³8³B1 − �³8³³ ;��� (62)

Balanço de energia ao queimador

Como já foi referido na secção 3.2, na análise da combustão define-se a temperatura adiabática

de chama que não é atingida devido às perdas por radiação para a câmara de combustão.

Designando a temperatura adiabática de chama por Tg1 e a temperatura dos gases após a perda de

calor por Tg2 e considerando ainda que este é o valor da temperatura da superfície correspondente

ao queimador para a câmara de combustão Ts0, obtém-se:

Aa %����®ℎ%_ − ℎ%B± = �0(O0 − ;0$ (63)

Balanço aos gases de combustão na câmara de combustão

Como o modelo criado assenta na premissa de que os gases na câmara de combustão são

cinzentos que absorvem e emitem energia pode-se efetuar um balanço a estes, para definir a

temperatura de saída desta zona Tg4. A temperatura para o cálculo das propriedades dos gases de

combustão na câmara é dada pela média geométrica da temperatura dos gases na câmara de

combustão Tg3=ªTg2Tg4 e, é de salientar que por se considerar os gases como cinzentos obtém-se a

seguinte relação para a emissividade em que εg=αg=1-τg. Desta forma o balanço exprime que a

variação de energia nos gases é a diferença entre a energia emitida e absorvida. Adicionalmente, a

troca de calor por convecção pode ser facilmente adicionada.

Aa %����®ℎ%B − ℎ%³± =í ��%«�þ_ ÐσT%«³ − O�Ñ + �_ℎ��(F(G%« − GN�_$ (64)

Balanço de energia aos tubos na câmara de combustão

Para o cálculo do balanço de energia aos tubos em volta da câmara de combustão é necessário a

distribuição de temperaturas de modo a calcular a condução de calor das paredes laterais da câmara

de combustão (e considerando que o tubo tem uma espessura aproximadamente desprezável).

Considerou-se então um modelo simples para definir a temperatura média da parede que é brasada

ao tubo.

Figura 3.4 – Exemplo do balanço às paredes laterais + tubos

A figura 3.4 apresenta um esquema de um tubo com placa com largura 2L que recebe o fluxo de

calor considerado uniforme q’’ e o transfere para o tubo pelo que a taxa de transferência de calor por

23

unidade de comprimento do tubo é dada por 2Lsepq’’. A temperatura média da placa é calculada com

base na análise da distribuição de temperaturas na direção normal aos tubos x. Considerando que o

calor transferido por radiação atua como uma fonte de calor interior, a distribuição de temperatura de

acordo com Incropera e DeWitt (2006) é dada por:

G(f$ = pqq2�~ ®����B − fB± + G_ (65)

Assim o balanço de energia transferida para um tubo (por comprimento de tubo) fica:

p�q = �2 P−�~ #G#f��þe\� = 2����pqq (66)

A distribuição de temperatura é importante para definir a temperatura média da parede Ts1 que é

dada por:

GN = G_ + pqq����B3�~ (67)

A temperatura da superfície do tubo pode ser aproximada pela temperatura média da água que

circula no seu interior T1=T(L$=TÒw= (Tw2+Tw3$ 2⁄ desprezando a resistência térmica oferecida pela

condução no tubo e a convecção na água. A temperatura média pode ser expressa então por:

GN�_ = ;_ − O_3�~ ����B + G5_ + G5B2 (68)

O balanço aos tubos é então dado pela seguinte expressão onde qenv,total representa as perdas da

envolvente que são explicadas na secção 3.6.

Aa 5=ℎ5« − ℎ5B$ = 4b�,��J����=;_ − O_$=f + �$�� ≅ �_=;_ − O_$ + �_ℎ��(F=G%« − GN�_$ − p���,����� (69)

3.4 Balanço ao permutador de calor primário

Nestes balanços é preciso ter em conta a transferência de calor por radiação da zona da câmara

de combustão pela superfície fictícia 2 e as trocas por convecção através de toda a superfície de

tubos e alhetas. Notar ainda que as propriedades da água são calculadas a Tm= (Tw1+Tw2) 2⁄ e que o

fator � pode ser obtido através do método �-NTU. Deste modo o balanço ao permutador de calor

primário pode ser feito considerando o lado da água e o lado dos gases de combustão. Note que a

transferência de calor para a água engloba neste caso trocas de calor radiativas e de convecção. No

entanto as trocas de radiação são efetuadas com a superfície das alhetas do permutador primário e

ignora-se os efeitos da radiação nos gases desta zona.

Aa 5=ℎ5B − ℎ5_$ = =;B − OB$ �B + c®G%³ − G5_±� (70)Aa %����®ℎ%³ − ℎ%¹± = c®G%³ − G5_±� (71)

A escolha da correlação a usar no cálculo do coeficiente global de convecção fica dependente da

configuração e geometria da caldeira de condensação. Assim, é de salientar que o permutador da

Bosch tem tubos alhetados mas nos espaços entre tubos existe um canal trapezoidal. Uma vez que

não existem correlações para uma configuração tão específica são escolhidas várias correlações, das

quais a melhor será posteriormente usada.

24

Para o coeficiente de convecção exterior considerou-se a correlação de Kim et al. (2008) para

tubos com alhetas do tipo herringbone, com os parâmetros geométricos indicados no anexo 1. Como

os tubos são elípticos aplicou-se a correção identificada na revisão bibliográfica (equação 10).

3.4.1 Método εεεε----NTUNTUNTUNTU----θθθθ

Uma vez determinado c é usado o método ε-NTU para determinar �, que representa a razão

entre a diferença média de temperatura entre os dois fluídos e a diferença entre as temperaturas de

entrada. Então, segundo o método ε-NTU:

bGc = c=Aa &�$(�� (72)

�� = =Aa &�$(��=Aa &�$(�� (73)

Assim, é possível deduzir �:

� = ü/ cü/®Aa &�±(��� = �bGc (74)

Assim, aplicam-se as expressões de Kays e London (1984) para escoamento de corrente cruzada

com uma fila de tubos. A equação 75 é utilizada quando (ma cp)min

encontra-se do lado dos gases de

combustão enquanto que a seguinte é utilizada quando (ma cp)min

existe para o lado do tubo, i.e. da

água.

� = 1�� Ð1 − 3²Ùô®_²�����±Ñ (75)

� = 1 − 3²®_²�����±Ùô (76)

3.5 Balanço ao condensador

Para o balanço à zona do condensador é necessário saber quantos tubos de condensadores

estão nessa zona para fazer os cálculos individualmente. É então necessário conhecer a temperatura

dos gases de combustão entre cada condensador. Como nem todo o vapor de água é condensado é

necessário um balanço de massa aos gases e ao vapor de água, uma vez que temos condensação,

de modo a determinar a humidade específica. Pela análise feita na indústria as caldeiras de

condensação raramente têm mais do que três condensadores pelo que é uma boa aproximação.

Assim sendo, o programa não será capaz de lidar com mais do que 3 condensadores sem que sejam

efetuadas ligeiras modificações no código e folha de cálculo. É necessário também, através de

cálculos, determinar qual o regime de escoamento no interior dos tubos e no caso de existirem

turbulators pode ser aplicada uma das correlações descritas para o permutador de calor primário.

Os números em parênteses são as posições em que são calculadas as massas e entalpias e

guardadas no programa para uso no condensador seguinte uma vez que entre condensadores existe

a necessidade de saber o que aconteceu no condensador anterior. A figura 3.5 elucida como se

processa o cálculo de maneira mais atrativa.

25

Figura 3.5 - Ilustração da zona de condensação

Os balanços apresentados de seguida são então a base para a criação das funções que regem a

zona de condensação.

Balanço de energia para um condensador

Os balanços seguintes refletem cada condensador p. Ao invés do método ε-NTU é utilizada a

temperatura média logarítmica nos balanços pois neste caso devido à condensação a capacidade

térmica dos gases varia e a diferença média logarítmica é uma aproximação para a diferença de

temperatura média como proposto por Herranz et al. (2000). Os índices i refletem as entradas dos

gases e j as saídas. Do lado da água k refere-se à saída do condensador e l à entrada. Note que Tw0

é a temperatura de entrada da água que é conhecida.

p� = Aa 5=ℎ50D − ℎ50�$ = c∆GeW (77) p� = Aa %����,�ℎ%� −Aa %����,9ℎ%9 −A�a ℎ� = c∆GeW (78)

∆GeW = ®G%� − G50D± − ®G%9 − G50�±�jЮG%� − G50D± ®G%9 − G50�±¾ Ñ (79)

Em que p=1,2,3, 6≤i≤9, 6≤j≤9, 0≤k≤2 e 0≤l≤2. Note também que para k=0 a temperatura será Tw0.

A entalpia de condensação à referência do vapor de água a 25°C é dada por:

ℎ� = −2442,3 + &n� PG���J + 298,152 \ =G���J − 298,15$ (80)

Onde a temperatura do condensado pode ser estimada pela seguinte expressão:

G���J = G50� + G%92 (81)

Para calcular a entalpia dos gases é necessário algum cuidado uma vez que só existe

efetivamente condensação quando TÒs≤Tsat. Em casos de carga parcial em que os gases arrefeçam

muito no permutador principal eventualmente podem-se atingir condições de saturação mas, admitiu-

se que não existe condensação sendo esta então contabilizada no primeiro condensador. Esta

aproximação não introduz erros significativos pois a diferença de temperatura entre os gases e a

água no permutador principal é muito elevada.

Então, para calcular a fração mássica d vapor de água presente após cada condensador utilizou-

se parte da expressão que calcula o coeficiente de convecção dos gases. Desta forma:

26

Aa � = M��ℎ���JΘ���J®GN% − GN�±ℎ�%q (82)

Para o modelo 1 a abordagem utilizada é aproximada, mas com uma ligeira alteração uma vez

que é assumido fluido puro neste modelo. Desta forma a expressão fica da seguinte forma:

Aa � = M��ℎ��(=G��� − GN�$ℎ�% (83)

Repare que o cálculo do caudal de condensado é sobrestimado uma vez que as frações mássicas

e/ou molares são avaliadas à entrada do condensador apesar do uso de temperaturas médias de

modo a reduzir o excesso induzido. Em suma, pode-se calcular o caudal de vapor na saída.

Aa �����@9 = Aa �����@� +Aa � (84)

Assim, é possível saber o caudal de vapor de água à saída do condensador, e desta forma é

possível saber o caudal de condensado uma vez que sabe-se o primeiro caudal de vapor de água à

entrada do primeiro condensador.

Posto isto e apesar das aproximações que são propostas, na presença de gases não

condensáveis existe sempre troca de calor sensível e calor latente. Serão então considerados dois

modelos diferentes que estão presentes na secção 2. O primeiro modelo considera apenas vapor de

água saturado nos gases de combustão (i.e. calor latente) enquanto o segundo considera a

existência de gases não condensáveis (i.e. calor sensível e calor latente). Em qualquer um dos

modelos não é usada nenhuma malha para refinar os resultados. Consideram-se apenas a entrada e

saída do condensador.

3.6 Perdas na envolvente da caldeira

Este modelo caracteriza as perdas por convecção natural e por radiação de modo simples e

considerando apenas a zona da câmara de combustão e ignorando a zona de condensação, uma vez

que a geometria da mesma é demasiado complexa para obter resultados fidedignos, e até porque é

na zona da câmara de combustão que se encontram as maiores perdas. Assim sendo para a

radiação temos a seguinte expressão segundo Kenna et al. (2007).

p���,��J = �������¢=GN���³ − G�³$ (85)

(86)

Quanto à convecção propõe-se que sejam calculados os números de Nusselt do envelope da

caldeira de condensação, que terá como hipótese um valor aproximadamente igual à emissividade do

lado interior da parede da caldeira de condensação. As propriedades do ar húmido serão retiradas

sempre à temperatura do filme Tf, ®TÒenv+T∞± 2⁄ . A temperatura da envolvente é considerada como a

temperatura média da parede.

A convecção natural pode ser calculada segundo Churchill et al. (1975). Considera-se apenas a

convecção natural em regime laminar e regime turbulento em que o número de Rayleigh crítico será

de 109 mas mesmo em regime turbulento o valor máximo estipulado para o número de Rayleigh será

de 1013.

27

�e = ��=GN��� − G�$�«��G� (87)

bgNNNNe = �0,825 + 0,3870 �e_ ¸⁄1 + =0,492 *l⁄ $Á _¸⁄ º/B·�B (88)

Deste modo o balanço de energia a considerar é o seguinte.

p���,����� = ℎN��� ���=GN��� − G�$ + �������¢=GN���³ − G�³$ (89)

Apesar destes cálculos é de salientar que os cálculos estão sobrestimados uma vez que as

perdas por radiação são desprezáveis nas paredes laterais dado que existe um involucro à volta da

câmara de combustão.

3.7 Propriedades dos fluídos

3.7.1 Propriedades dos gases

Os gases de combustão desempenham um papel determinante na transferência de calor. À parte

da radiação emitida pelo queimador (assumido como uma placa) é preciso reparar que as

transferências de calor a partir do permutador de calor (sem condensação) são dominadas pela ação

convectiva dos gases. Assim é necessário determinar as propriedades radiativas e convectivas dos

gases de combustão.

A entalpia dos gases de combustão pode ser calculada usando o modelo dos gases ideais e com

recurso a polinómios criados por Kee et al. (1987). Estes coeficientes derivam dos dados recolhidos

Burcat et al. (2005) pelo que foram usados coeficientes dos polinómios que garantem um intervalo de

temperaturas maior e uma precisão melhor. Os coeficientes encontram-se então no anexo 2. Tem-se

os seguintes polinómios para a entalpia dos gases que a expressão está em base molar pelo que é

preciso multiplicar pela massa molar da mistura.

ℎN0 0G = �_ + �B2 G + �«3 GB + �³4 G« + �¹5 G³ + �G (90)

A capacidade calorifica específica pode ser obtida da mesma fonte.

&�Ò 0 = �_ + �BG + �«GB + �³G« + �¹G³ (91)

Os limites de validade destes polinómios vão do 200 K até 1000 K sendo que nos 1000 K aos

6000K os coeficientes mudam.

É necessária ainda a viscosidade dinâmica e condutibilidade térmica para finalizar as

propriedades necessárias para definir o modelo. A viscosidade dinâmica pode ser determinada

recorrendo à equação 93. Note-se que Mk representa a massa molar de cada um dos componentes

dos gases de combustão. Já a condutibilidade térmica pode ser obtida com recurso à equação 94 em

que Z=1+0,402�H2O0,511®1-XH2O±.

28

�%�� =∑ �D Ö�D ª`D ×�Dþ_

∑ �D ª`D �Dþ_ (92)

�%�� = 1,085 Ëí=�D�D$�

Dþ_Ì � (93)

O coeficiente de difusão pode ser tirado com recurso à expressão presente em Incropera e De

Witt (2006). Esta expressão é simplista mas os valores não variam muito das mais atuais correlações.

Então DH2O-Ar=0,26×10-4 m2/s mas no entanto este valor deve ser corrigido para um coeficiente de

difusão do vapor de água interagindo com os gases de combustão e não com ar.

)ÚUݲR���� = )ÚUݲ£� P G298\« B⁄ Í�%����=G$���=G$ Î (94)

Quanto à concentração molar pode ser retirada usando igualmente expressões de gás ideal,

Cgases= ∑ Xiρgases/Mgases = n/ 0Gni . A massa específica é então retirada de igual modo.

�%���� =ín� � 0G��þ_ (95)

As propriedades para o cálculo das perdas com a envolvente da caldeira (i.e. cálculo do

coeficiente de convecção natural) podem ser dadas pelas expressões presentes no anexo 2

(Tsilingiris, 2007).

3.7.2 Propriedades da água

As propriedades da água há muito que estão determinadas e convencionadas segundo a

formulação IAPWS de Wagner et al. (2002). Assim sendo o erro induzido pelas seguintes expressões

é mínimo e está estipulado. Então para a água é necessário saber a massa específica,

condutibilidade térmica, tensão superficial, viscosidade (dinâmica) e claro, as entalpias.

A entalpia da água bem como a massa específica pode ser retirada a partir da expressão da

pressão de vapor da água (pressão parcial) e da sua derivada de onde ϑ=1- T Tc⁄ com Tc=647,096 K,

pc=22,064 MPa e ρc=322±3 Kg/m-3.No caso da entalpia é necessária mais uma expressão, a

diferença calorífica em que α0=1 KJ/Kg-1 e Λ= T Tc⁄ . Os coeficientes das expressões estão no anexo

2.

vj Pn�n� \ = G�G =�_¦ + �B¦_,¹ + �«¦« + �³¦«,¹ + �¹¦³ + �¸¦·,¹ (96)

#n�#G = −n��%(�G Óvj Pn�n�\ + �_ + 1,5�B¦0,¹ + 3�«¦B + 3,5�³¦B,¹ + 4�¹¦« + 7,5�¸¦¸,¹Ô (97)

��0 = #� + #_�²_Á + #B� + #«�³,¹ + #³�¹ + #¹�¹³,¹ (98)

29

Sabendo estas três expressões é possível retirar as expressões que definem a massa específica

e a entalpia da água enquanto líquido saturado ou vapor saturado. É possível saber a entalpia de

evaporação de forma trivial. Também é possível saber a temperatura de saturação aplicando o

método de Newton-Raphson, apresentado por Kelly (2003), com cálculo da derivada por diferenças

finitas progressivas de primeira ordem à expressão da pressão de saturação.

���� = 1 + �_¦_ «⁄ + �B¦B «⁄ + �«¦¹ «⁄ + �³¦_¸ «⁄ + �¹¦³« «⁄ + �¸¦__0 «⁄ (99)

vj P����\ = &_¦B ¸⁄ + &B¦³ ¸⁄ + &«¦º ¸⁄ + &³¦_º ¸⁄ + &¹¦«· ¸⁄ + &¸¦·_ ¸⁄ (100)

ℎ� = ��0 + 10« G�� P#n�#G \ (101)

ℎ� = ��0 + 10« G�� P#n�#G \ (102)

ℎ�% = −ℎ� + ℎ� (103)

No caso de tensão superficial da água é possível obter a seguinte expressão presente no

documento da IAPWS (1994) em que B=235,8 mN/m, b=-0,625 e µ=1,256. A expressão tem uma

incerteza de 0,30 e uma diferença máxima em relação ao valor experimental de 6% perto do ponto

crítico.

¢ = �¦�=1 + �¦$ (104)

Quanto à condutibilidade térmica esta pode ser calculada usando a expressão mais recente da

IAPWS de 2011. Neste trabalho não foi considerado o fator melhoramento crítico ÐλN2®TÒ,ρÒ±Ñ quer para

a condutibilidade térmica quer para a viscosidade dinâmica, sendo que este valor é 1 uma vez que

este valor só apresenta valores significativos perto do ponto crítico. Repare que kN= k kc⁄ com

kc=1×10-3Wm-1K-1, TÒ=θ e ρÒ=ρ ρc⁄ .

�N = �N0(GN) + �N_(GN, �) + �NB(GN, �) (105)

O primeiro fator representa a condutibilidade térmica no gás diluído pelo que não depende da

massa específica.

�N0(G) = ªGN∑ �DGND³Dþ0

(106)

O segundo fator representa a contribuição devido à massa específica ser finita.

�N_(GN, �) = 3fn ��í�P1GN − 1\�í ��9¹9þ0 =� − 1$9�³

�þ0 � (107)

A capacidade calorifica específica para a água líquida saturada ou vapor saturado pode ser

retirada usando a equação de Shomate e dados do National Institute of Standards and Technology

(NIST) para os coeficientes. Para o líquido saturado note que o erro obtido ao levar o limite inferior de

298 K para 273,16 K é mínimo pelo que pode muito bem ser levado até 273,16 K.

30

&� = + �~ + �~B + )~« + 2/~B (108)

Finalmente, a viscosidade dinâmica para líquido saturado pode ser calculada usando a seguinte

expressão de IPAWS (2008). Repare que µÒ=µ/µc em que µc=1×10-6 Pa.s.

�� = �0=GN$ + �_=GN, �$ + �B=GN, �$ (109)

�0=G$ = 100ªGN∑ DGND«Dþ0

(110)

�_(GN, �) = 3fn ��í�P1GN − 1\�í�9¸9þ0 =� − 1$9�¹

�þ0 � (111)

Os coeficientes de Lk, Lij, Hk e Hij estão no anexo 2. Assim sendo, validade destas expressões é

ampla e está descrita na tabela seguinte. Note que Tm=p$ é a temperatura de fusão dependente da

pressão e que à medida que a pressão aumenta o limite máximo da temperatura vai baixando até

atingir 348 K na gama 785≤p≤1000 MPa.

Tabela 3.1 - Limite de validade das expressões para a água

Propriedade(s) Limites de validade

pσ,ρl, ρv, hl, hv,σ,hfg 273,16 ≤ T ≤ 647,096 K

k,µ 273,16 ≤ T ≤ 1173,15 K para 0 ≤ p ≤ 611,657 Pa, Tm(p) ≤ T ≤ 1173,15 K para

611,657 Pa ≤ p ≤ 100 MPa

cpl 273,16 ≤ T ≤ 500 K

cpv 500 ≤ T ≤ 6000 K

3.8 Metodologia

3.8.1 Esquemática de resolução do modelo

De modo a instruir o utilizador do programa a melhor lidar com o mesmo serão apresentadas

duas figuras que mostram as entradas a efetuar em cada folha de cálculo de modo ao programa

funcionar na perfeição. Também constam nas figuras as saídas uma vez pedido os resultados na

folha “Solve” ao clicar no botão “Results”.

31

Figura 3.6 - Representação das entradas e saídas das folhas de cálculo representadas a negrito

Figura 3.7 - Representação das entradas e saídas das folhas de cálculo representadas a negrito

3.8.2 Método de resolução de sistemas de equações não lineares

Com a modelação feita e todos os parâmetros e tendo em conta as propriedades dos fluidos em

questão é necessário então resolver o sistema de equações não lineares. O método usado é um

método quasi-Newton conhecido pelo método de Broyden (1965).

Este método é uma generalização do método da secante, ilustrado na figura 3.9, que recorre a

uma relação de recorrência para resolver uma equação. Este método precisa de duas estimativas

iniciais e não tem convergência quadrática mas sim supra-linear o que leva este método a ser mais

lento a convergir do que o método de Newton. No entanto não é preciso conhecer a derivada da

equação uma vez que esta relação de recorrência pode ser pensada como diferenças finitas.

32

Figura 3.8 - Ilustração do método de Broyden

Assim sendo o método de Broyden tem uma iteração inicial a partir da qual temos apenas uma

estimativa inicial que calcula f��_ que pode ser entendido como a segunda estimativa inicial. De

modo a calcular a derivada são utilizadas diferenças finitas para a primeira iteração como para as

restantes. Logo, é possível ver que o método é muito sensível à primeira iteração, dado que esta

pode levar a segunda estimativa pobre, levando a que o sistema não convirja. Note-se que o cálculo

da matriz jacobiana é feito utilizando diferenças finitas mas no entanto para calcular xn+1 é necessário

efetuar o inverso da matriz jacobiana, tal como no método de Newton-Raphson. Este problema pode

ser ignorado no método da secante.

Então, o método de Broyden não melhorado funcionaria deste modo e mantém as características

do método da secante como ter convergência supra-linear para além não necessitar da expressão

exata da jacobiana. O único senão, em relação à secante, é o fato de ser necessário calcular a matriz

inversa da jacobiana, algo que o no método de Newton também é necessário. No entanto o método

usado na realidade é o método de Broyden melhorado que utiliza a fórmula de Sherman-Morrison,

presente em Kelly (2003), que se trata de um caso particular da fórmula de Woodbury e que consiste

em atualizar a matriz inversa usando um vetor perturbação que é adicionado à matriz evitando assim

novo cálculo da matriz inversa que é um processo exigente a nível computacional.

Assim sendo, o método evita dois aspetos importantes. O cálculo da jacobiana exata e o cálculo

da matriz inversa a cada iteração. Os problemas relativos à primeira iteração são então ainda mais

importantes uma vez que todo o método fica a depender da primeira iteração e pode não conseguir

convergir. Aliás, até no método da secante não garantido que ‖F‖ tenha uma direção descendente

(i.e. convirja) pelo que é necessário garantir esta situação com boas estimativas e uma matriz

jacobiana aproximada com uma precisão suficiente.

Assim, a matriz jacobiana é calculada com recurso a diferenças finitas progressivas presentes no

artigo de Bryan et al. (2004), que tornam o método inexato e mais propenso a erros ainda que as

33

diferenças finitas progressivas usadas sejam de terceira ordem de modo a tentar reduzir o erro ao

máximo relativamente aos valores exatos. Paradoxalmente o uso de diferenças finitas de terceira

ordem podem fazer o erro aumentar mais do que usando diferenças finitas de primeira ordem uma

vez que, se a aproximação inicial for pobre (má) levará a que o erro maior. Também é de referir que

um sistema, do qual as suas funções e restrições estejam mal desenhadas ao nível computacional,

fará qualquer tentativa de retirar uma matriz jacobiana próxima do valor exato infrutífera.

De forma a manter a convergência em praticamente qualquer situação são utilizadas duas

condições que podem ser ativadas de maneira independente e que só devem ser utilizadas caso o

método melhorado falhe.

A primeira condição é a condição de descida suficiente é proposta por Kelly (2003) e permite um

certo amortecimento da convergência (damping) o que em teoria diminui o número de iterações a

realizar até à convergência desde que a convergência esteja a efetuar-se na direção certa mas que

para casos em que não aconteça a condição irá aumentar o número de iterações podendo mesmo

levar à falha em convergir.

A outra condição é a condição de razão de convergência garantida proposta por Bader et al.

(2005) é uma medida que deve ser aplicada a sistemas pobres e que garante uma resposta muito

satisfatória em termos de convergência e de robustez uma vez que sempre a razão de convergência

tem um valor insatisfatório o método reinicia com os valores da última iteração. Este valor deve ir de

0,1 a 1 segundo o autor. Isto leva a que a precisão aumente muito e que o número de iterações locais

(após o reinício) diminua mas não as iterações globais que aumentam dependendo do valor dado à

razão de convergência, mas tal não implica que o sistema convirja. A desvantagem é que a matriz

jacobiana tem que ser calculada mais que uma vez mas ainda assim menos vezes do que no método

de Newton. Desta forma, para os casos em que existe convergência com ou sem a utilização da

condição, é necessário verificar se o uso da condição leva ou não a um maior tempo de convergência

dado que o mesmo exige o cálculo de mais inversas.

Ainda a título de precaução foi adicionada uma condição que impede a temperatura da água na

fase líquida de exceder valores não previstos nas funções que calculam as propriedades da água.

Esta condição apenas pode ser aplicada a conjuntos de temperaturas consecutivos (células

consecutivas de temperaturas da água). Quanto à convergência esta é garantida pela norma

euclidiana de F.

Os problemas que o método de Broyden (melhorado e com condições impostas) pode apresentar

são os seguintes:

• O jacobiano obtido no método de Broyden no final não tem necessariamente que ser

o verdadeiro jacobiano segundo Biegler (2000)

• O valor do incremento hDF para o jacobiano é muito importante e um valor mal

introduzido pode levar à não convergência.

34

• Não existe garantia de que o método convirja satisfatoriamente usando qualquer uma

das condições se o sistema for mal condicionado.

• Se a estimativa inicial se encontrar muito longe da solução o comportamento pode

ser errático. O aumento do número de iterações poderá apenas piorar a situação.

• O uso da condição de descida suficiente ou da razão de convergência pode aumentar

em muito o tempo de convergência em sistemas de equações bem condicionados e

que facilmente convergem com o método de Broyden melhorado.

• Das duas condições a última que deve ser usada é a condição de descida suficiente

uma vez que é uma condição simplista de um método muito mais complexo e porque

o sistema pode instabilizar dado que a falta de saltos significativos entre a estimativa

inicial e a solução aproximada pode significar que o sistema nunca converge ou

converge muito mais lentamente. Deixa-se então esta condição como última a ser

usada quando tudo o resto falhou.

São assim apresentadas de seguida as expressões que determinam o método usado remetendo-

se paras referências os teoremas que as suportam. Estas expressões são apresentadas de forma

generalizada, ignorando as especificidades do problema a ser resolvido. Como o programa foi

criado com base nessa premissa é então fácil adaptar o método para resolver outros problemas

inclusive problemas puramente matemáticos.

Admitindo que o sistema a resolver é definido por:

7=�$ = 0 (112)

E a matriz jacobiana é dada por:

�0 = 8q=f0$ = −8=f + 2ℎKL$ + 8=f + ℎKL$ − 38=f$ − 8=f − ℎKL$6ℎKL + �=ℎKL«$ (113)

O valor de ℎKL tem que ser dado pelo utilizador não deve ser muito pequeno nem muito grande,

sendo que deve situar-se num valor escalonado relativamente à escala da estimativa inicial. Posto

isto temos então a expressão para a qual é tirado o primeiro vetor o que leva às seguintes

expressões da qual dB,n é a direção de broyden e λB,n é o passo que pode ser introduzido pela

condição de descida suficiente (senão é igual 1).

f_ = f0 − �0²_8=f0$ (114) f��_ = f� − á1,���²_8=f�$ (115) #1,� = −��²_8=f�$ (116)

A matriz inversa é então atualizada da seguinte forma usando a fórmula de Sherman-Morrison.

���_²_ = ��²_ + ®!1,� − ��²_�1,�±!1,�� ��²_!1,�� ��²_�1,� (117)

Em que sB,n=xn+1-xn é a diferença entre dois cálculos consecutivos do vetor x, sB,nt a sua transposta e

yB,n=F=xn+1$-F=xn$ .

35

Já a razão de convergência é simplesmente a razão das normas euclidianas de duas iterações

consecutivas. χn= ‖Fn+1 ‖ ‖Fn ‖⁄ . A condição de descida suficiente pode então ser descrita através da

equação seguinte em que A é retirado iterativamente quando a equação verifica-se e assim obtém-se

λB=2-m.

�8=f� + 2²(#1,�� < =1 − 10²³2²($‖8=f�$‖ (118)

A implementação do método no ExcelTM foi realizada com sucesso e o modo de construção do

sistema de equações tentou ser o mais universal possível como já foi dito acima. A descrição de uso

do método está presente no anexo 3.

36

4 Resultados

O modelo criado foi feito com base em dados da caldeira de condensação fornecidos pela Bosch.

A Bosch também forneceu dados de testes efetuados à caldeira de condensação em questão pelo

que o modelo criado para esta caldeira pode ser testado para verificar a qualidade dos dados. Os

dados disponíveis permitem uma validação geral do modelo da caldeira de condensação mas, não

permitem verificar em detalhe as evoluções em cada uma das zonas de condensação separadas.

Uma vez que a principal finalidade do modelo é testar a influência do ar húmido, i.e. humidade

relativa, na eficiência da caldeira de condensação as comparações são efetuadas para o

desempenho global.

Para a zona da câmara de combustão foi necessário avaliar a emissividade da parede de modo a

permitir que o calor de acordo com o modelo seja o observado a partir da temperatura da água. O

modelo é aproximado pois considera apenas uma zona na radiação, tal como as propriedades das

paredes (incluindo as fictícias) são aproximadas. No entanto após a determinação da emissividade

pretende-se que o modelo reproduza variações nas condições de operação. No caso do permutador

de calor primário, o coeficiente global de transferência de calor é influenciado por ambos os

coeficientes de convecção interior e exterior, sendo comparadas diferentes correlações para ambos.

Desta forma procede-se à validação dos dados fornecidos (anexo 4). As configurações possíveis

na análise da temperatura do gás estão na tabela 4.1. A configuração mais completa é a configuração

3 pelo que será a mais usada ao longo dos resultados. Os dados estão então separados por ensaios

que estão no anexo 4. O ensaio mais usado será o ensaio #1 (configuração 3) pelas razões já

indicadas.

Tabela 4.1- Várias configurações na análise do gás

Configurações

1 2 3 4

Permutador de calor

primário com 111

alhetas + 2 Tubos

alhetados de Al

Permutador de calor

primário com 118

alhetas + 2 Tubos

alhetados de Al

Permutador de calor

primário com 118

alhetas + 2 Tubos

alhetados de Al com

turbulators

Permutador de calor

primário com 118

alhetas + 3 Tubos

alhetados de Al com

turbulators

Desta forma é necessário quantificar a rendimento da caldeira de condensação sendo que esta é

estabelecida pela expressão seguinte.

� = üa�a (119)

Em que Qa representa a soma das parcelas úteis de cada zona da caldeira de condensação e é

um valor fixo enquanto Ca pode ser visto como PCI*ma fuel ou PCS*ma fuel, dependendo para que região

37

do mundo os testes são efetuados, sendo por isso variável consoante seja usado PCI ou PCS. A

conversão pode então ser feita da seguinte maneira.

�ðÙ� = �ðÙã P*�k*�^\ (120)

Tem-se ainda que existem perdas no gás ao longo das zonas de transição entre os permutadores

que não estão contabilizadas e que devem e podem fazer parte do erro obtido na verificação das

temperaturas bem como perdas do combustível devido a uma queima deficiente.

No entanto, antes da apresentação dos resultados da caldeira em si, é necessário testar o

método numérico para verificar se a convergência ocorre para um número de iterações aceitável,

bem como para uma tolerância razoável.

4.1 Análise do método numérico

As figuras 4.1 e 4.2 pretendem atestar estas situações utilizando o ensaio #3 (em anexo) como o

sistema a resolver. Deste modo verifica-se que a condição de razão de convergência leva a que seja

necessárias menos iterações ainda que leve mais tempo a resolver o sistema. Embora não presente

nas figuras a razão de convergência também permite melhorar a tolerância mínima. No entanto como

a figura 4.1 ilustra, é necessário escolher uma razão de convergência adequada. O valor pode variar

mas no entanto a tendência é quanto menor melhor até um mínimo de 0,1 segundo Bader et al.

(2005). A figura 4.2 também mostra que a condição de descida suficiente varia muito pouco a cada

iteração o que leva a que a não ser que a estimativa inicial esteja próxima da solução da equação

então esta condição torna o método inviável daí a sugestão deixada acima.

Em suma, o método de Broyden melhorado funciona, em geral, de maneira aceitável e sem

grandes custos computacionais. Utilizando a condição de razão de convergência é possível reduzir o

número de iterações ainda que aumente o tempo para atingir a solução pelo que esta condição deve

ser utilizada quando a diferença entre o número de iterações do método melhorado e a condição de

razão de convergência for mais que um certo valor a decidir pelo utilizador. A condição de descida

suficiente é muito demorada e leva o sistema a não convergir pelo que não é aconselhada a não ser

em última instância em sistemas mal condicionados.

38

Figura 4.1 - Número de iterações utilizando várias combinações do método de Broyden disponíveis

Figura 4.2 - Comparação do número de iterações necessário para tolerância de 108

4.2 Análise da Combustão

Uma maneira simples de confirmar que a combustão está a ser bem calculada é calcular a

variação da temperatura adiabática de chama com o excesso de ar uma vez que esta variação é

conhecida, ou seja, o gradiente par λ>1 a e λ<1 é diferente. Posto isto para metano e utilizando os

dados do teste a 51,3kW a 0% de humidade relativa obtém-se a figura 4.3.

0

5

10

15

20

25

30

1,00E-08

me

ro d

e i

tera

çõe

s

Tolerância

Broyden melhorado

Broyden melhorado com Razão

Conv. 0,9

Broyden melhorado com Razão

Conv. 0,5

Broyden melhorado com Razão

Conv. 0,2

Broyden melhorado com Razão

Conv. 0,1

-10

-8

-6

-4

-2

0

2

4

6

1 6 11 16

Log

(||

F|

|)

Número de iterações

Broyden melhorado

Broyden melhorado

com Razão Conv.

0,2

Broyden melhorado

com Desc. Suf.

39

Figura 4.3 - Gráfico exemplificativo da evolução da temperatura adiabática com o excesso de ar para

HR=0%

Desta maneira é possível verificar que o máximo ocorre para um excesso de ar unitário.

Considerando dissociação este máximo ocorreria para um excesso de ar ligeiramente inferior a 1 uma

vez que a capacidade calorifica dos produtos de combustão decresce mais rapidamente para λ<1 até

ao máximo, sendo que após este máximo este decréscimo é mais lento. A figura 4.4 ilustra o porquê

do gradiente não ser igual para um excesso de ar maior ou menor que 1 ainda que sem contabilizar a

dissociação. Não esquecer que a percentagem de N2 está dividida por 10.

Figura 4.4 - Gráfico da evolução das frações molares com o excesso de ar para HR=0%

1500

1650

1800

1950

2100

2250

2400

0,7 0,8 0,9 1 1,1 1,2 1,3 1,4 1,5 1,6

T a

dia

tica

(K

)

Coeficiente de excesso de ar

0,0%

2,5%

5,0%

7,5%

10,0%

12,5%

15,0%

17,5%

20,0%

22,5%

25,0%

0,7 0,8 0,9 1 1,1 1,2 1,3 1,4 1,5 1,6

Fra

cçã

o m

ola

r (%

)

Excesso de ar

%CO

%CO2

%H2O

%N2/10

%O2

40

É preciso ainda verificar se o comportamento da temperatura adiabática de chama é consistente

ao dar-se um aumento da humidade relativa. Segundo a figura 4.5 essa evolução corresponde à

evolução esperada, uma vez que um aumento e/ou existência de vapor de água nos reagentes levará

a que a temperatura máxima atingida seja cada vez menor, pois parte da energia deixa de estar

disponível para aumentar a temperatura dos produtos de combustão.

Figura 4.5 - Temperatura adiabática de chama vs. Humidade relativa

Ainda considerando a validação da combustão, e utilizando os dados do ensaio #1 com uma

humidade relativa de 0%, é interessante verificar a influência dos diferentes combustíveis nas

parcelas de transmissão de calor para água por convecção e/ou por radiação na câmara de

combustão e permutador primário. Segundo Miranda et. al. (2003) um aumento de emissividade na

parede leva a um aumento da parcela de radiação pelo que, é do interesse então ver o que acontece

com a mudança de combustível uma vez que, a parte convectiva é influenciada pelo coeficiente

global de convecção do permutador de calor primário e assim é trivial assumir que um aumento deste

coeficiente faz aumentar a parcela de transmissão de calor por convecção. Posto isto, a figura 4.6

indica que a utilização de um combustível mais “pesado” (i.e. propano) leva a que a parcela de

radiação diminua e que a influência do caudal de água é desprezável. Tal é esperado uma vez que,

apesar de a percentagem de gases participantes na zona da câmara de combustão aumentar, a

temperatura na câmara de combustão é inferior quando é usado propano. Este resultado pode ser

verificado recorrendo aos resultados de Miranda et. al. (2003) que indica um maior valor de CO e de

NOx para o butano e metano do que para o propano. Os valores obtidos estão de acordo com as

percentagens obtidas por Kenna et al. (2007).

2260

2280

2300

2320

2340

0% 20% 40% 60% 80% 100%

Te

mp

era

tura

(K

)

HR (%)

41

Figura 4.6 - Comparação da influência dos combustíveis nos mecanismos de transmissão de calor

Note-se que para este caso a variação de temperatura na água e nos gases de combustão não

respeita a metodologia para deduzir os métodos θ-ε-NTU, uma vez que existe uma componente de

radiação no balanço do permutador primário.

4.3 Análise das Perdas

As perdas foram calculadas utilizando o ensaio #1 e o resultado encontra-se na Figura 4.7. Como

é fácil inferir, a influência da humidade relativa nas perdas por convecção natural e radiação é

diminuta e desprezável em situações de potências elevadas como as de aquecimento de águas

sanitárias ou aquecimento central. No entanto estas perdas são com o exterior, ou seja, não estão

contabilizadas as perdas devido à ineficiência da caldeira às diferentes cargas (i.e. potências). A

percentagem máxima de perdas, em relação ao input, foi de 0,9%.

Desta forma é fácil reparar, ainda assim, que as perdas diminuem com o aumento da humidade

relativa, algo justificado por uma diminuição do coeficiente de convecção natural com o aumento da

humidade relativa e uma diminuição também da temperatura da parede dado que a temperatura de

“entrada” dos gases de combustão baixa com o aumento da humidade.

20%

22%

24%

26%

28%

30%

32%

34%

36%

38%

40%60%

62%

64%

66%

68%

70%

72%

74%

76%

78%

80%

0 5 10 15 20 25 30 35

Ra

dia

ção

(%

)

Co

nv

ecç

ão

(%

)

Caudal mássico de água (l/min) para ΔT=41,37°C e ε=0,8

Propano

Metano

42

Figura 4.7 - Perdas pela envolvente da zona de combustão vs. Humidade relativa

4.4 Análise do lado dos gases de combustão

De modo a validar o modelo é essencial validar as propriedades e correlações usadas. Como já

foi referido, apenas é possível efetuar uma análise qualitativa pois os dados relativos às temperaturas

intermédias são médios e para potências desconhecidas. Segundo a Bosch a variação de

temperatura, da zona de combustão é de aproximadamente 100°C e para o permutador de calor

primário esta variação é de aproximadamente 1000°C, ambos em relação à temperatura dos gases à

entrada da câmara de combustão, para potências elevadas. Desta forma, e admitindo um desvio de

10°C (10%) na figura 4.8 e um desvio de 100°C (10%) para a figura 4.9, foram calculados os valores

da emissividade da parede da zona de combustão e os valores de temperatura obtidos de acordo

com cada uma das três correlações disponíveis para o permutador de calor primário. Foram utilizados

os dados do ensaio #1 (anexo 4).

302

304

306

308

310

312

314

0% 20% 40% 60% 80% 100%

Pe

rda

s e

nv

elo

pe

(W

)

HR (%)

43

Figura 4.8 - Variação da temperatura na zona de combustão vs. Emissividade da parede em

situações extremas de humidade relativa

Figura 4.9 - Variação de temperatura no permutador de calor primário vs. Correlação usada

De acordo com a figura 4.8 é verifica-se que o valor qualitativo de 100°C não é cumprido nos

testes disponíveis a para altas potências. Sendo assim foi fixado um valor de 0,7 para a emissividade

das paredes que, sendo realista pois este valor é utilizado muitas vezes por defeito, também é o mais

próximo do desvio em relação aos 100°C. Também de acordo com a mesma figura é possível

verificar a influência do dióxido de carbono na variação de temperatura. Um aumento da percentagem

de gases participantes, i.e. xH2O, leva a um aumento das trocas de calor aumentando a variação de

70

80

90

100

110

120

130

0 0,1 0,2 0,3 0,4 0,5 0,6 0,7 0,8 0,9 1

Dif

ere

nça

de

te

mp

era

tura

(°C

)

Emissividade da parede

ensaio #1,

HR=0%

ensaio #1,

HR=100%

ensaio #2,

HR=0%

850

900

950

1000

1050

1100

1150

1 2 3

Te

mp

era

tura

(°C

)

Correlação

HR=0% Pot.

49,59 kW

HR=100% Pot.

49,59 kW

Kim 1997 Kim et al. 2008 Youn e Kim 2007

- 10%

+10%

- 10%

+ 10%

44

temperatura. Já quanto à influência da adição de trocas de calor por convecção na câmara de

combustão tem-se que, a temperatura após a câmara de combustão irá baixar mais pelo que as

trocas de calor no permutador primário serão ligeiramente menores.

Quanto à figura 4.9, esta demonstra que qualquer uma das correlações é válida para o modelo e

ensaio em questão assumindo um desvio de 10% em relação aos 1000°C qualitativos. Desta forma,

será utilizada a correlação de Kim et al. (2008) como já foi referido na apresentação do modelo. Esta

figura tem por objetivo validar tal escolha uma vez que à partida na secção 2 existem mais

parâmetros que estão dentro dos limites de validade.

Assim sendo, utilizando os dados do ensaio #1 (anexo 4), a tabela 4.2 indica a variação das

temperaturas intermédias e de saída dos gases e rendimento em função da humidade relativa, para a

potência nominal de 49,59 kW e uma fração molar de CO2 de aproximadamente 7,65%, em que se

conhece o valor experimental da temperatura de saída dos gases mas não a humidade relativa. Para

todos os ensaios foi utilizado o modelo de condensação que considera os gases não condensáveis

uma vez que este modelo é o único que caracteriza de maneira realista o problema em questão.

Tabela 4.2 - Resultados dos ensaios #1, #2 e #5

Ensaio Simulação /

Experimental

Temperatura intermédia

dos gases de

combustão (°C)

Temperatura de

saída dos gases de

combustão (°C)

Rendimento

relativo a

PCI

#1

(49,59 kW)

Experimental -------- 52,40 102,16%

Simulação, HR=0% 196,00 85,70 101,21%

Simulação,

HR=100% 198,72 93,24 101,43%

#2

(58,28 kW)

Experimental 250,00 50,00 103,73%

Simulação, HR=0% 220,48 96,17 100,08%

Simulação,

HR=100% 223,41 103,04 100,16%

#5

(58,28 kW)

Experimental 255,00 52,00 101,53%

Simulação, HR=0% 219,20 154,77 95,52%

Simulação,

HR=100% 220,07 162,44 95,33%

A tabela 4.2 permite então verificar que para o ensaio #1 a temperatura dos gases de saída

aumenta com o aumento da humidade relativa, assim como a temperatura intermédia (i.e.

imediatamente após o permutador primário). De salientar que o desvio máximo para a temperatura de

saída dos gases nestes ensaios é de 110°C devido aos desvios dos cálculos nas duas zonas.

Pela tabela 4.2 ainda é possível verificar que apesar de existir condensação em todos os ensaios,

o ensaio #5 apresenta uma diminuição de rendimento. Tal fica a dever-se ao fato de ao não existirem

45

turbulators nos tubos alhetados (i.e. condensadores) o coeficiente de convecção no interior baixar o

que aumenta a temperatura dos tubos, aumentando a temperatura da interface (referenciada no

modelo dos gases não condensáveis) para o qual a driving force (ou driving forces se forem

considerados os fatores de correção) diminuem pelo que as trocas por calor latente diminuem. Desta

forma existe condensação mas não de maneira significativa para aumentar o rendimento. Posto isto

pode-se inferir que parte do vapor de água age como massa de inerte, diminuindo o rendimento

apesar de a temperatura de saída dos gases aumentar (o que é contraditório).

Quanto à evolução do rendimento o modelo subestima os valores experimentais e também é de

notar que o aumento de rendimento com a humidade relativa é mais significativo para o ensaio #1.

Tal é justificado pelo fato de este aumento ser mais acentuado à medida que a potência na caldeira

desce (i.e. temperaturas no condensador). A fundamentação para esta situação é a mesma que a

explicada para o ensaio #5 pois a driving force é muito sensível à temperatura da interface entre o

filme de condensado e a camada de gases não condensáveis.

Já tabela 4.3 indica a temperatura dos gases de combustão nas diferentes posições da caldeira.

É sugerido que a figura 3.1 seja consultada para uma melhor compreensão. Note que apesar de

algumas temperaturas não revelarem serem as verdadeiras aquando de um ensaio experimental a

diferença obtida nessa mesma zona é a correta. Tal deve-se à maneira como o modelo foi

desenvolvido.

Tabela 4.3 - Temperatura dos gases de combustão em diferentes posições para o ensaio #1

(H.R.=50%)

Diferentes posições Temperatura (°C) G%B 1197,54 G%³ 1084,29 GN�B 31,34 GN�_ 125,87 G%¹ 197,37 G%· 122,60 G%º 90,12

Para o ensaio #2 (configuração 3) já existem valores medidos para a temperatura intermédia.

Como a configuração deste ensaio é a mesma do ensaio #1 é possível retirar mais informação dos

resultados apresentados pelo modelo. A evolução das temperaturas e rendimentos é a mesma do

ensaio #1. Já os desvios máximos obtidos (no ensaio #2) para as temperaturas intermédias e de

saída dos gases são, respetivamente, de 30°C e de 53°C. Apesar de uma grande discrepância na

temperatura de saída dos gases de combustão a diferença da potência não será tão significativa. O

desvio máximo em relação ao rendimento experimental é de 6,2%. Note que, como descrito no anexo

46

4, qualquer para parâmetro omisso para o ensaio #2, utilizam-se os mesmos parâmetros que no

ensaio #1.

Quanto às correlações é possível verificar que a correlação utilizada para o permutador primário

prevê razoavelmente os resultados experimentais. Já quando à correlação para os permutadores

seguintes é de salientar que existe um desvio mais acentuado em relação aos resultados

experimentais umas que, devido aos mecanismos complexos (i.e. trocas de calor por convecção e

difusão de massa) devido à presença de gases não condensáveis é mais difícil obter bons resultados.

Uma vez que o ensaio #5 mantém o número de permutadores de calor e que apenas não

apresenta turbulators na região de condensação em relação aos outros, permite verificar que a

diminuição do rendimento (média) é de, aproximadamente 4,7% relativamente ao ensaio com

turbulators à mesma potência (ensaio #2). Esta variação é consistente embora maior do que a

variação apresentada pela Bosch em que o rendimento desce 2,2%.

Para o caso dos ensaios #3 e #4 apresentam-se as figuras 4.10 e 4.11 em alternativa às tabelas.

Os desvios nas temperaturas intermédias são praticamente iguais uma vez que a única alteração

entre o ensaio #3 e o #4 é o facto de o último ter 3 condensadores. Assim sendo o desvio é à volta de

8,5°C. O desvio máximo para a temperatura de saída dos gases varia entre os 9,1°C e os 9,9°C

enquanto o desvio máximo do rendimento para o ensaio #3 é de 7,3% e para o ensaio #4 é de 4,6%.

As variações entre uma humidade relativa nula e uma humidade relativa de 100% foram de 7,3%

(ensaio #3) a 8,4% (ensaio #4). Assim o aumento no rendimento (média) devido ao uso de mais um

condensador foi de 1,3% para uma humidade relativa de 100%. Relativamente aos dados da Bosch

esse aumento foi de 4,4%.

Figura 4.10 - Gráfico da temperatura de saída do permutador primário e de saída dos gases de

combustão para os ensaios #3 (linhas vermelhas) e #4 (linhas amarelas) e a tracejado os valores

medidos

0

5

10

15

20

25

30

35

40

45

50

55

60

65

0% 20% 40% 60% 80% 100%

Te

mp

era

tura

(°C

)

HR (%)

47

Figura 4.11 -Gráfico com o rendimento da caldeira de condensação vs. HR para os ensaios #3 (linha

vermelha) e #4 (linha amarela) e a tracejado os rendimentos medidos

Na análise da evolução da temperatura de saída dos gases é necessário analisar também o

rendimento. Isto deve-se ao facto de a temperatura de saída dos gases ser inferior ao valor medido e

o rendimento também ser inferior para uma gama de humidades relativas no ensaio #4. Uma vez que

não é conhecida a fração de vapor de água na exaustão não é possível adiantar mais informações

mas, no entanto, o fato de as linhas de rendimento ultrapassarem a referência indica que a fração de

calor latente efetivamente trocada nos condensadores é mais elevada do que o esperado. A título

indicativo então é mostrada a tabela 4.4 com as percentagens de calor latente relativamente ao calor

total trocado nos condensadores. Também é possível verificar que com o aumento da humidade

relativa a contribuição das trocas de calor latente aumentam, o que é expetável.

Tabela 4.4 - Comparação calor total trocado e calor latente para ensaio #3

HR Contribuição do calor latente

0% 43,1%

40% 44,3%

80% 45,6%

100% 46,3%

Analisando ainda a partir do lado dos gases de combustão, as diferenças obtidas mas desta vez

ao nível da diferença de potências então é possível deparar-se com os seguintes dados presentes na

tabela 4.5 que ilustram a entalpia de saída dos gases de combustão à referência de 25°C. Estes

dados refletem que embora os desvios na temperatura de saída dos gases pareçam assinaláveis,

não são tão acentuadas ao analisar-se em termos de potência, pelas razões já referidas. Esta análise

é feita a uma humidade relativa de 50% obtendo-se assim, aproximadamente, uma média da

potência.

106%

108%

110%

112%

114%

116%

118%

0% 20% 40% 60% 80% 100%

Re

nd

ime

nto

(%

)

HR (%)

48

Tabela 4.5 – Entalpias de saída dos gases de combustão e diferença obtida em unidades de potência

para diferentes ensaios

Ensaio Valor obtido (kJ/Kg) Valor experimental (kJ/Kg) Diferença (kW)

1 (49,59kW) 68,53 29,83 0,88

2 (58,28kW) 80,02 30,01 1,23

3 (5,82kW) -6,74 1,03 0,45

4 (5,82kW) -12,79 -2,07 0,16

5 (58,28kW) 144,40 30,01 2,86

4.5 Análise do lado da água

Para o lado da água, devido ao elevado calor específico da mesma a variação será menor mas

não menos importante. Aliás, o interesse último de uma caldeira de condensação é a obtenção

transferência de calor eficiente entre os gases de combustão e a água.

A tabela 4.6 mostra as variações obtidas entre uma humidade relativa de 0% e de 100% na

temperatura de saída da água da caldeira de condensação para o primeiro documento em que existe

o valor da temperatura de saída experimental e também para as diferentes configurações das quais o

valor experimental é tirado indiretamente. Como era esperado, as variações na temperatura de saída

da água são mínimas devido ao elevado calor específico da água relativamente aos gases de

combustão, pelo que para a potência mínima o efeito da humidade relativa é negligenciável.

Tabela 4.6 - Temperatura de saída da água da caldeira de condensação para diferentes ensaios em

situações extremas de humidade

Ensaio Humidade relativa (%) Valor experimental (°C) Diferença obtida (°C)

0% 100%

1 (49,59kW) 53,38 53,47 53,72 0,25 – 0,34

2 (58,28kW) 60,03 60,07 60,73 0,70 – 0,66

3 (5,82kW) 18,0 18,4 17,45 0,55 – 0,95

4 (5,82kW) 17,46 17,74 17,45 0,01 – 0,29

5 (58,28kW) 57,86 57,77 60,73 2,96 – 2,87

4.6 Outras caldeiras

O objetivo da dissertação é a modelação de um modelo corrente de uma caldeira de

condensação que sirva de indicativo para a indústria e não como modelo universal que deva ser

utilizado em qualquer situação. No entanto, certas partes do ficheiro ExcelTM podem ser utilizadas de

modo geral uma vez que existem fatores que não devem mudar de modelo para modelo. Em suma o

modelo criado deve ser usado com precaução e devem ser feitas simulações exaustivas de modo a

determinar as melhores correlações e adaptações para uma determinada caldeira e/ou grupo de

caldeiras. Neste caso de generalização são utilizadas as mesmas correlações que as do modelo

testado anteriormente de modo a manter a coerência do trabalho e consequentemente da análise.

49

Posto isto, a Bosch enviou documentos de mais 3 tipos de caldeiras cujos dados encontram-se

no anexo 4. A primeira caldeira usa gás propano (ensaio #6) e as outras duas caldeiras usam gás

metano (ensaios #7 e 8). Estes três modelos não contêm condensadores pelo que a tentativa de

generalizar o modelo centrar-se-á exclusivamente na câmara de combustão e no permutador

primário. Os tubos alhetados têm uma geometria semelhante ao da caldeira de condensação usada

para a modelação mas não são exatamente iguais, uma vez que as alhetas diferem ligeiramente. As

dimensões dos três equipamentos também são diferentes entre si.

As figuras 4.12 e 4.13 apresentam uma comparação entre os rendimentos calculados (para

humidades relativas extremas) para as caldeiras nas quais não existem condensadores e para a

caldeira de condensação, na qual se baseou o modelo, em relação aos valores experimentais.

Calculando o desvio-padrão definido por:

)* = Ê∑ hf���J f���¾ − 1iBb (121)

Considerando ar seco para combustão observa-se um desvio padrão dos rendimentos de 3,63%,

e à humidade máxima foi de 4,46%.

Figura 4.12 - Comparação de rendimentos experimentais com rendimentos calculados para +-10%

para HR=0%

60%

70%

80%

90%

100%

110%

120%

60% 70% 80% 90% 100% 110% 120%

ηp

red

(%

)

η exp (%)

+10%

-10%

50

Figura 4.13 - Comparação de rendimentos experimentais com rendimentos calculados para +-10%

para HR=100%

Estes resultados mostram uma boa concordância entre os rendimentos obtidos

experimentalmente e os rendimentos obtidos através do modelo criado para o efeito. Assim sendo, as

figuras indicam que o modelo prevê razoavelmente os rendimentos de caldeiras sem condensação e

com condensação.

De acordo com a tabela seguinte as diferenças obtidas à humidade mínima (0%) para a

temperatura de saída dos gases de combustão são maiores que na caldeira sobre a qual o modelo foi

baseado. Já a temperatura da água estimada teve uma diferença máxima de 1°C no ensaio #7.

Tabela 4.7 - Diferença obtida para os gases de combustão na saída para os diferentes ensaios

Ensaio Valor obtida (°C) Valor experimental (°C) Diferença (°C) Diferença (kW)

6 (22,2 kW) 276,1 202,0 74,1 0,99

7 (27,0 kW) 159,9 195,0 35,1 0,72

8 (32,8 kW) 276,77 223,0 53,8 0,96

60%

70%

80%

90%

100%

110%

120%

60% 70% 80% 90% 100% 110% 120%

ηp

red

(%

)

η exp (%)

-10%

+10%

51

5 Conclusões

A partir dos resultados obtidos é possível reter alguns aspetos importantes relativos à influência

da humidade relativa nas caldeiras de condensação. O primeiro prende-se com o facto de que para

potências elevadas, quando a influência da humidade relativa aumenta, o rendimento pode diminuir

quando a preponderância do calor latente não é significativa, devido a uma fraca driving force. Assim

sendo, a massa adicional de vapor de água devido ao aumento da humidade relativa é considerada

parcialmente como sendo massa inerte, absorvendo energia e fazendo diminuir o rendimento.

No entanto, para potências baixas o rendimento aumenta com o aumento da humidade relativa,

variando entre 7,3% e 8,4% para as humidades relativas extremas para os dois ensaios a baixa

potência. Para estas potências o calor latente já é preponderante pelo que também é possível

verificar que as temperaturas irão baixar menos do que o esperado se não existisse nenhum calor

latente, uma vez que o calor sensível não é dominante. Já o incremento verificado devido ao

condensador adicional no ensaio #4 foi de 1,3%. Note-se que estes valores dependem do número de

condensadores e das condições de entrada.

Posto isto, aumento na temperatura de saída da água relaciona-se pouco com o aumento da

humidade relativa a não ser quando existe condensação devido à sua elevada capacidade calorífica.

Nesse caso a influência é assinalável em termos proporcionais.

Constata-se então que os maiores rendimentos estão para as potências mais baixas uma vez que

existe mais condensação a essas potências. Isto significa que a caldeira de condensação deve

controlar a potência em função das condições de entrada/funcionamento de modo a obter o melhor

rendimento possível pois a potências elevadas os gases de combustão não atingem temperatura de

saturação do vapor de água, perdendo-se o calor latente. As conclusões de Kuck (1996) devem

portanto de ser abordadas em benefício do rendimento dadas as condições de entrada.

A verificação dos modelos para cada um dos componentes não pode ser efetuada de forma

precisa. Com base nas correlações da literatura e no modelo de uma zona da câmara de combustão

obtiveram-se desvios de 10% em relação aos valores que são típicos de acordo com o fabricante. É

conveniente notar que diferenças da primeira zona da câmara de combustão influenciam as

seguintes. Ainda assim o modelo consegue desvios-padrão para o rendimento entre 3,63% e 4,46%

consoante a humidade relativa para as mais diversas geometrias e sem a necessidade de ensaios

laboratoriais a priori ao contrário dos modelos referenciados na secção 1.2. Quanto às perdas (em

termos de potência) estas são mínimas quando comparadas com a potência nominal.

Verificou-se também que o modelo de condensação pura (modelo 1) não foi validado uma vez

que não corresponde à realidade do problema. Já o modelo que tem em conta a presença de gases

não condensáveis (modelo 2) obtém bons resultados que podem ser melhorados posteriormente

através de um estudo detalhado aos condensadores.

Deixa-se então ao cuidado de trabalhos futuros o melhoramento do modelo através de ensaios

experimentais exaustivos nomeadamente na zona de condensação e ao nível das perdas uma vez

52

que a criação de correlações de modo a prever perdas consoante o tipo de modelo (ventilação

assistida, recuperação de calor nas laterais) era essencial para melhor prever o rendimento obtido

nas simulações. O estudo das zonas, segmentando-as em malha e testando em laboratório as

temperaturas intermédias a várias potências e diferentes excessos de ar também iria ser uma mais-

valia uma vez que é possível quantificar os erros em zonas mais pequenas e verificar de uma melhor

maneira os resultados obtidos.

No ponto da generalização são necessários ainda mais dados de caldeiras (com ou sem

condensação) de modo a obter valores dos desvios e erros mais fidedignos. Assim sendo, a

colaboração com mais empresas pode ajudar a evoluir o modelo matemático que, apesar de

simplista, consegue dar resultados aproximados de modo a caracterizar completamente uma caldeira

de condensação a um custo muito menor que os modelos CFD.

Espera-se portanto que este trabalho tenha seguimento com a continuação do desenvolvimento

do modelo matemático que, apesar de simples é prático e pode ser usado para analisar a influência

de parâmetros geométricos ou de condições de operação no desempenho da caldeira. Este trabalho

deve ser sempre comparado com resultados de ensaios mas constituí uma metodologia que pode

minimizar o número de ensaios e de variantes a testar.

53

6 Referências Bibliográficas

1. Adamek, T., Bestimmung der kondensationgröβen auf feingewellten oberflächen zur

auslegung optimaler wand-profile, Wärme-und Stoffübertragung 15, 255-270 (1981)

2. Adamek, T. A.; Webb, R. L., Prediction of film condensation on horizontal integral-fin tubes,

Int. J. of Heat and Mass Transfer 33, 1721-1735 (1990)

3. Anthony, Stephan, An instantaneous condensing gas-fired water heater: Modeling and

performance, PhD Dissertation – Purdue University (1986)

4. Azevedo, T., Complementos de radiação – Uso de superfícies fictícias, IST (2000)

5. Bader, B.; Pawlowski, Roger P. and Kolda, Tamara G., Robust Large-scale parallel nonlinear

solvers for simulations, Sandia Report (2005)

6. Bella, B.; Cavallini, A.; Longo, G. A.; Rossetto, L., Pure vapor condensation of refrigerants 11

and 113 on horizontal integral finned tube at hight vapor velocity, J. of Enhanced Heat Transfer 1,

77-86 (1994)

7. Biegler, L. T., Systems of Nonlinear Equations – Part II, Carnegie Mellon University Document

(2000)

8. Bir, R. B.; Stewart, W. E.; Lightfoot, E. D., Transport Phenomena. Wiley, New York (1960)

9. Brauer, H., Compact heat exchangers, Chem. Process Eng. (August), 451-460 (1964)

10. Briggs, D. E.; Young, E. H., Convection heat transfer and pressure drop of air flowing across

triangular pitch banks of finned tubes, Chem. Eng. Prog. Symp. Ser. 41, Vol 59, 1-10 (1963)

11. Broyden, C. G., A class of methods for solving nonlinear simultaneous equations, Math.

Comput. 19, 577-593 (1965)

12. Bryan, K. and Shibberu, Y., Finite Difference Approximation of Derivatives, Mathematical

Modelling Classes (2004)

13. Burcat, A.; Ruscic, B., Third Millennium Ideal Gas and Condensed Phase Thermochemical

Database for Combustion with updates from Active Thermochemical Tables, Argonne National

Laboratory (2005)

14. Churchill, S. W.; Chu, H. H. S., Int. Journal of Heat and Mass Transfer, 18, 1323 (1975)

15. Coelho, P.; Costa, M., Combustão, 1ª Edição (2007)

16. Collier, J. G.; Thome, J. R., Convective Boiling and Condensation, 3rd Ed., Oxford Univ. Press

(1994)

17. Eiamsa-ard, S.; Promvonge, P., Influence of Double-sided Delta wing Tape insert with

Alternate-axes on Flow and Heat Transfer Characteristics in a Heat Exchanger Tube, Chinese, J.

of Chem. Eng. 19(3), 410-423 (2011)

18. Eiamsa-ard, S.; Kongkaitpaiboon, V.; Promvonge, P., Thermal Performance Assessment of

Turbulent Tube Flow Through Wire Coil Turbulators, Heat Transfer Eng. 32:11-12, 957-967 (2011)

19. Fitzgerald, Claire L.; Briggs, A.; Rose, John W.; Wang, Hua Sheng, Effect of vapour velocity

on condensate retention between fins during condensation on low-finned tubes, International

Journal of Heat and Mass Transfer 55 (2012) 1412-1418

20. Gogonin, J. L.; McQuiston, F. C., An experimental investigation of air dehumidfication in

parallel plate exchanger, ASHRAE Transactions 29, 146-149 (1973)

54

21. Goldschmidt, V. W., Harris, D. K., Measurements of the overall heat transfer from combustion

gases confined within elliptical tube heat exchanger, Exp. Therm. Fluid Sci. 26, 33-37 (2002)

22. Goldstein, L.; Sparrow, E. M., Experiments on the transfer characteristics of a corrugated fin

and tube heat exchanger configuration, J. Heat Transfer 98, 26-34 (1976)

23. Gnielski, V., New equations for heat and mass transfer in turbulent pipe and channel flow, Int.

Chem. Eng. 16 (2), 359-368 (1976)

24. Gnielski, V., A new calculation procedure for the heat transfer in the transition region between

laminar and turbulent pipe flow, Forschung im Ingenieurwesen 6, 240-248 (1995)

25. Hanby, V. I., Modelling the performance of condensing boilers, Journal of the Energy Inst. 80,

229-231 (2007)

26. Herranz, L. E., Desarrollo y validación de modelos avanzados de transmisión de calor en los

condensadores de refrigeración passiva de la contención SBWR, Ph.D. Thesis, PolyThecnical

Uni. Of Madrid, (1996)

27. Herranz, L. E.; Anderson, M. H.; Corradini, M. L., A diffusion layer model for steam

condensation within the AP600 containment, Nucl. Eng. Design 183, 133-150 (1998)

28. Herranz, L. E.; Muñoz-Cobo, J. L.; Palomo, M. J., Modeling condensation heat transfer on a

horiztontal finned tube in the presence of noncondensable gases, Nuclear Engineering and

Design 201, 273-288 (2000)

29. Hubert, D.; Walter, H., Forced Convection Heat Transfer in the Transition Region Between

Laminar and Turbulent Flow for a vertical circular tube, Latest Trends on Theoretical and Applied

Mechanics – Fluid Mech. And Heat & Mass Transfer (2010)

30. IAPWS, Release on IAPWS Formulation 2008 for Viscosity of Ordinary Water Substance

(2008)

31. IAPWS, IAPWS Release on Surface Tension of Ordinary Water Substance (1994)

32. IAPWS, Release on the IAPWS Formulation 2011 for the Thermal Conductivity of Ordinary

Water Substance (2011)

33. Ibrahim, Talaat A. and Gomma, A., Thermal performance criteria of elliptic tube bundle in

crossflow, International Journal of Thermal Sciences 48 (2009) 2148-2158

34. Incropera; DeWitt; Bergman; Lavine, Fundamentals of Heat and Mass Transfer, 6th Edition

(2006)

35. Inman, Phillip, The Guardian Article, (2005)

http://www.guardian.co.uk/money/2005/apr/02/consumerissues.jobsandmoney

36. Kays, W. M., London, A. L., Compact Heat Exchangers 3rd Ed., McGraw-Hill NY (1984)

37. Kee, R. J.; Rupley F. M.; Miller, J. A., CHEMKIN-II: A Fortran Chemical Kinetics Package for

the Analysis of Gas-Phase Chemical Kinetics, Sandia National Laboratories (1989)

38. Kelly, C. T., Solving Nonlinear Equations with Newton’s Method, Fundamentals of Algorithms

(2003)

39. Kenna, R.; Elbur, M.; Li, W.; Holsteijn, R., Preparatory Study on Eco-Design of Water Heaters,

Report prepared for European Commission (2007)

55

40. Kim, N. H.; Yun, J. H.; Webb, R. L., Heat transfer and friction correlations for wavy plate fin-

and-tube heat exchangers, Journal Heat Transfer 119 (1997), 560-567

41. Kim, Y.; Kim, Y., Heat transfer characteristics of flat plate finned-tube heat exchangers with

large fin pitch, Int. Journal of Refrigeration 28, 851-858 (2005)

42. Kim, Nae-Hyun; Ham, Jung-Ho; Cho, Jin-Pyo, Experimental investigation on airside

performance of fin-and-tube heat exchangers having herringbone wave fins and proposal of a new

heat transfer and pressure drop correlation, Journal of Mechanical Science and Technology 22

(2008) 545-555

43. Khan, Mesbah G.; Fartaj, A.; Ting, David S.-K., An experimental characterization of cross-flow

cooling of air via an in-line elliptical tube array, Int. Journal of Heat and Fluid Flow 25, 636-648

(2004)

44. Kuck, J., Efficiency of vapour-equipped condensing boilers, Applied Thermal Eng. 16, 233-244

(1996)

45. Makaire, D.; Ngendakumana, P., Thermal performances of condensing boilers, 32nd TLM –

IEA Energy Cons. And Emissions Reduc. In Combustion – Japan (Nara) (2010)

46. Meng, Ji-An; Liang, X.-G.; Chen, Z.-J.; Li, Z.-X., Experimental study on convective heat

transfer in alternating elliptical axis tubes, Experimental Thermal and Fluid Science 29, 457-465

(2005)

47. Miranda, M.; Barbosa, R.; Pinho, C., Thermal Behavior of a Pressurized Water Heater,

Proceedings of COBEM, 10-14 (2003)

48. Morisot, O., Modèle de batterie froide à eau glacée à la maûtrise des consommantion

d’énergie en conception de bâtiments climatisés et en conduite d’installation, Thèse de doctorat –

Ecole des Mines de Paris (2000)

49. Nakai, S.; Okazaki, T., Heat transfer from horizontal circular wire at small Reynolds and

Grashof numbers-I. Pure convection, Int. J. Heat and Mass Transfer 18, 387-396 (1975)

50. Nishiyama, R. J.; Ota, T.; Matsuno, T., Heat transfer and flow around elliptic cylinders in

tandem arrangement, JSME Int. J., Ser. II 31 (3), 410-419 (1988)

51. Peterson, P. F., Theoretical basis for the Uchida correlation for condensation in reactor

containments, Nucl. Eng. Design 162, 301-306 (1996)

52. Promvonge, P, Thermal augmentation in circular tube with twisted tape and wire coil

turbulators, Energy Conv. And Management, 49, 2949-2955 (2008)

53. Promvonge, P.; Eiamsa-ard, S., Heat transfer augmentation in a circular tube using V-nozzle

turbulator inserts and snail entry, Exp. Thermal and Fluid Science 31, 332-340 (2007)

54. Warnsholz, R., High Efficiency Condensing Boilers, Superintendents meeting – R. M. Cotton

Company (2009)

55. Radiative Heat Transfer, Michael F. Modest, Academic Press, 2nd Ed. (2003)

56. Schimdt, Th. E., Heat transfer calculations for extended surfaces, Refrigeration Eng. 351

(1949)

57. Sieder, E. N.; Tate, G. E., Heat transfer and pressure drop of liquids in tubes, Ind. Eng. Chem.

28, 1429-1436 (1936)

56

58. Shah, R. K.; London, A. L., Laminar Forced Convection in Ducts, Academic Press NY (1978)

59. Shi, X.; Che, D.; Agnew, B.;Gao J., An investigation of performance of compact heat

exchanger for latent recovery from exhaust flue gases, Int. Journal of Heat and Mass Transfer 54,

606-615 (2011)

60. Tao, Y. B.; He, Y. L., Wu, Z. G.; Tao, W. Q., Three-dimensional numerical study and field

synergy principle analysis of wavy fin heat exchangers with elliptic tubes, International Journal of

Heat and Fluid Flow 28 (2007) 1531-1544

61. The Warren Report, EIBI July/August (2007)

62. Tsilingiris, P. T., Thermophysical and transport properties of humid air at a temperature range

between 0°C and 100°C, Exergy Conversion and Management, 49, 1098-110 (2007)

63. Wagner, W. and Pruβ A., The IAPWS Formulation 1995 for the Thermodynamic Properties of

Ordinary Water Substance for General and Scientific Use (2002)

64. Webb, R. L.; Rudy, T. M., An analytical model to predict condensate retention on horizontal

integral-fin tubes, Proceedings ASME-JSME Thermal Eng. Joint Conference 1, 373-378 (1983)

65. Webb, R. L.; Rudy, T. M.; Kedzierski, M. A., Prediction of condensation coefficient on

horizontal integral-fin tubes, ASME J. of Heat Transfer 107, 369-376 (1985)

66. Williams, C., Evaluate functions and Formulas fun: How to make Excel’s Evaluate method

twice as fast (2011)

http://fastexcel.wordpress.com/2011/11/02/evaluate-functions-and-formulas-fun-how-to-make-

excels-evaluate-method-twice-as-fast/

67. Xie, G.; Wang, Q., Sunden, B., Parametric study on multiple correlations on air side heat

transfer and friction characteristics of fin-and-tube heat exchangers with large-diameter tube rows,

Applied Thermal Engineering 29 (2009) 1-16

68. Youn, B. and Kim, N. H., An experimental investigation airside performance of fin-and-tube

heat exchangers having sinusoidal wave fins, Heat Mass Transfer 43 (2007) 1249-1262

57

Anexos

Anexo 1. Outras correlações

Correlações utilizadas para comparar os resultados obtidos nos ensaios com as simulações

efetuadas na zona do permutador primário.

Kim et al. (1997)

�W�« = ¿0,987 − 0,01bn�l� 3Kw ≥ 1000, b = 1,21,35 − 0,162bn�l� 3Kw < 1000, b = 1,2 (122)

�« = 0,394 3Kw²0,«¹· P*e*e\²0,B·B P8�)�\²0,B0¹ P��*K\²0,¹¹º P*K8� \²0,_«« (123)

Youn e Kim (2007)

�W�« = ¿0,911 − 0,0061bn�l� 3Kw > 1000, b = 1,21,59 − 0,137bn�l� 3Kw ≤ 1000, b = 1,2 (124)

�« = 0,214 3Kw²0,«¹B P*e*e\0,«¸_ P8�)�\²0,¸³0 P��*K\²_,0B P*K8� \²0,¸¹¹ (125)

Tabela A.1 - Limite de validade das correlações para escoamento exterior

3Jw dc* Pt Pl

* Pf NL Φ =°$ Xf*

Kim (1997) 500-6000 9,53-12,7 25,4-31,3 22-27,5 2,08-4,22 1-4 5,7-34,7 3,67-5,5 3Jw *d *e⁄ * 8�/)� ��/*K *K/8� NL Φ =°$

Youn (2007) 700-4000 1,15-1,33 0,12-0,33 1,44-10,0 0,23-1,74 1-4 5,7-39,8 * - Indicativo de que esta gama não é cumprida para o modelo em causa.

Correlação de Leckner para cálculo da emissividade dos gases de combustão na câmara de

combustão.

Seja a emissividade dada por:

�% = �ÚB0 + �ÙÝB − ∆�=n5�, n��$ (126)

f� = j�j = n�n (127)

� = �( = 3,6���(F ��(F (128)

Em que i pode corresponder a BÞ ou �ÞB e 0 corresponde à dita espécie a G� = 1000], n0 = 1 ��l. Outras expressões que mantêm-se para qualquer espécie são =n��) = 1 ��l &A, ~ = G/G0

e � = ��U !��U !��" U. Note que n� corresponde à pressão parcial da espécie.

�� = �0�=n��, 1��l, G%«$ P ��0\� ®n��, n, G%«± (129)

58

�0 = exp �í í &9� PG%«G0 \9 Plog_0 n��=n��$� \� W

9þ0

î

�þ0� (130)

P ��0\ = 1 − =� − 1$=1 − *&$� + � − 1 + *& exp Ë−& Ívu�_0 =n��)( =n��$� ÎBÌ (131)

∆� = P �10,7 + 101� − 0,0089�_0,³\�log_0hnÚUÝ! + nÙÝUi �

=n��$0 �B,·¸ (132)

Os dados relativos à correlação de Leckner encontram-se na seguinte tabela.

Tabela A.2 - Tabela com dados da correlação de Leckner

'() CO2 *,+ 2,2 2,3 ,--…,/-

,-0…,/0

-2,2118 -1,1987 0,035596 0,85667 0,93048 -0,14391

-0,10838 -0,17156 0,045915

-3,9893 2,7669 -2,1081 0,39163 1,2710 -1,1090 1,0195 -0,21897

-0,23678 0,19731-0,19544 0,044644

PE Pn + 2,56n�√~ \ n0¾ =n + 0,28n�$ n0⁄

=234$5=234$6 13,2~B 70,054~B , ~ < 0,70,225~B , ~ > 0,7

3 ¿ 2,144, ~ < 0,751,888 − 2,053 vu�_0~ , ~ > 0,75 =1 + 0,1$ t_,³¹⁄

b 1,10 ~_,³⁄ 0,23

c 0,5 1,47

59

Anexo 2. Propriedades

Propriedades da água

Coeficientes da pressão de vapor:

a1=-7,85951783, a2=1,84408259,a3=-11,7866497,a4=22,6807411,a5=-15,9618719,a6=1,80122502

Coeficientes da massa específica (líquido saturado):

b1=1,99274064, b2=1,09965342, b3=-0,510839303, b4=-1,75493479,

b5=-45,5170352,b6=-6,74694450×105

Coeficientes da massa específica (vapor saturado):

c1=-2,03150240,c2=-2,68302940,c3=-5,38626492,c4=-17,2991605,c5=-44,7586581,c6=-63,9201063

Coeficientes da propriedade calorífica especial:

dα=-1135,905627715,d1=-5,65134998×10-8,d2=2690,66631,d3=127,287297,d4=-135,003439,d5=0,981825814

Coeficientes da condutibilidade térmica 8 para os dois primeiros fatores da expressão:

Tabela A.3 - Coeficientes de �D

k Hk

0 1,67752

1 2,20462

2 0,6366564

3 -0,241605

Tabela A.4 - Coeficientes de �9 j

0 1 2 3 4 5 6 i

0 5,20094×10-1 2,22531×10-1 -2,81378×10-1 1,61913×10-1 -3,25372×10-2 0 0

1 8,50895×10-2 9,99115×10-1 -9,06851×10-1 2,57399×10-1 0 0 0

2 -1,08374 1,88797 -7,72479×10-1 0 0 0 0

3 -2,89555×10-1 1,26613 -4,89837×10-1 0 6,98452×10²B 0 -4,35673×10-34 -2,7203370 4,57586331 -2,57040×10-1 0 0 8,72102×10-3 0

5 0 1,20573×10-1 0 0 0 0 -5,93264×10²

Coeficientes da viscosidade dinâmica para os dois primeiros fatores da expressão:

Tabela A.5 - Coeficientes de �D

k Lk

0 2,443221×10-3

1 1,323095×10-2

2 6,770357×10-3

60

3 -3,454586×10-3

4 4,096266×10-4

Tabela A.6 - Coeficientes de ��9 j

0 1 2 3 4 5 i

0 1,60397357 -0,646013523 0,111443906 0,102997357 -0,054123634 0,00609859258

1 2,33771842 -2,78843778 1,53616167 -0,463045512 0,0832827019 -0,00719201245

2 2,19650529 -4,54580785 3,55777244 -1,40944978 0,275418278 -0,0205938816

3 -1,21051378 1,60812989 -0,621178141 0,0716373224 0 0

4 -2,7203370 4,57586331 -3,18369245 1,1168348 -0,19268305 0,012913842

Coeficientes da capacidade térmica específica cpv:

Tabela A.7 - Coeficientes da equação de Shomate

Temperatura 500-1700 1700-6000

A 30,09200/18 41,96426/18

B 6,832514/18 8,622053/18

C 6,793435/18 -1,499780/18

D -2,534480/18 0,098119/18

E 0,082139/18 -11,15764/18

Coeficientes da capacidade térmica específica cpl:

Tabela A.8 - Coeficientes da equação de Shomate

Temperatura 273,16-500

A -203,6060/18

B 1523,290/18

C -3196,413/18

D 2474,455/18

E 3,855326/18

61

Propriedades dos gases de combustão

Coeficientes para o cálculo das propriedades dos gases de combustão:

Tabela A.9 - Coeficientes dos polinómios dos gases de combustão para uma dada gama de temperaturas

CO a1 a2 a3 a4 a5 a6 a7

1000-6000 3,0484859 1,3517281E-3 -4,8579405E-7 7,8853644E-11 -4,6980746E-15 -1,4266117E+4 6,0170977

200-1000 3,5795335 -6,1035369E-4 1,0168143E-6 9,0700586E-10 -9,0442449E-13 -1,4344086E+4 3,5084093

CO2

1000-6000 4,6365111 0,0027414569 -0,00000099589759 1,6038666E-10 -9,1619857E-15 -49024,904 -1,9348955

200-1000 2,356813 0,0089841299 -0,0000071220632 2,4573008E-09 -1,4288548E-13 -48371,971 0,99009035E+01

H2O

1000-6000 2,6770389 0,0029731816 -0,00000077376889 9,4433514E-11 -4,2689991E-15 -29885,894 6,88255

200-1000 4,1986352 -0,0020364017 0,0000065203416 -5,4879269E-09 1,771968E-12 -30293,726 -0,84900901

N2

1000-6000 2,95257637 0,0013969004 -0,000000492631603 7,86010195E-11 -4,60755204E-15 923,948688 5,87188762

200-1000 3,53100528 -0,000123660988 -0,000000502999433 2,43530612E-09 1,40881235E-12 -1046,97628 2,96747038

O2

1000-6000 3,66096083 0,000656365523 -0,000000141149485 2,05797658E-11 -1,29913248E-15 1215,97725 3,41536184

200-1000 3,78245636 -0,00299673415 0,000009847302 -9,68129508E-09 3,24372836E-12 1063,94356 , 3,65767573

Tabela A.10 - Coeficientes dos polinómios dos combustíveis para uma dada gama de temperaturas

CH4 a1 a2 a3 a4 a5 a6 a7

1000-6000 1,911786 0,0096026796 -0,00000338387841 5,3879724E-10 -3,19306807E-14 10099,2136 8,48241861

62

200-1000 5,14825732 -0,013700241 0,0000493749414 -4,91952339E-08 1,70097299E-11 10245,3222 -4,63322726

C2H6

1000-6000 4,04666411 0,0153538802 -0,00000547039485 8,77826544E-10 -5,23167531E-14 12447,3499 -0,968698313

200-1000 4,29142572 -0,00550154901 0,0000599438458 -7,08466469E-08 2,68685836E-11 -11522,2056 2,66678994

C3H8

1000-6000 6,6691976 0,0206108751 -0,00000736512349 1,18434262E-09 -7,0691463E-14 -16275,4066 -13,1943379

200-1000 4,21093013 0,00170886504 0,0000706530164 -9,20060565E-08 3,64618453E-11 -14381,0883 5,61004451

CH3OH

1000-6000 3,52726795 0,0103178783 -0,00000362892944 5,77448016E-10 3,42182632E-14 26002,8834 5,16758693

200-1000 5,65851051 -0,0162983419 0,0000691938156 -7,58372926E-08 2,8042755E-11 25611,9736 -0,897330508

C2H2

1000-6000 4,65878489 0,00488396667 -0,00000160828888 2,46974544E-10 -1,38605959E-14 25759,4042 -3,99838194

200-1000 0,808679682 0,0233615762 -0,0000355172234 2,80152958E-08 2,80152958E-08 26428,9808 13,9396761

C2H4

1000-6000 3,99182724 0,0104833908 -0,00000371721342 5,94628366E-10 -3,53630386E-14 4268,65851 -0,269081762

200-1000 3,95920063 -0,00757051373 0,0000570989993 -6,91588352E-08 2,6988419E-11 5089,77598 4,09730213

C8H18

1000-6000 20,9430708 0,0441691018 -0,0000153261633 2,30544803E-09 -1,29765727E-13 35575,5088 -81,0637726

200-1000 12,524548 -0,0101018826 0,00022199261 -0,000000284863722 1,12410138E-10 -29843,4398 -19,7109989

C10H21

1000-6000 23,3759939 0,0575362038 -0,000021096802 3,46309821E-09 -2,0943403E-13 -18765,8614 -86,8825727

200-1000 14,8510661 -0,00385330874 0,000243710502 -0,000000317547576 1,25908377E-10 -12786,1722 -23,5344919

63

Propriedades do ar húmido

Utilizando o trabalho de Tsilingiris (2007) é então possível retirar as propriedades necessárias

para o cálculo das perdas, considerando a existência de vapor de água no ar. As seguintes equações

representam a proposta de cálculo das propriedades segundo o autor. Note que o subscrito �

representa o ar seco, 9 o vapor de água e A o ar húmido.

�( = =3,484 − 1,317��$ n1000G (133)

�( = �� + ��=0,7887�� − ��$1 − 0,2113�� (134)

&�( = &�� + 0,622 n�n − n� (135)

�( = �� + ��=0,8536�� − ��$1 − 0,1464�� (136)

A partir destas quatro equações a difusividade térmica e o número de Prandtl são facilmente

calculados.

�( = �(�(&�( (137)

*l( = �( &�(�( (138)

De seguida são apresentadas as correlações para o ar seco e o vapor de água para a massa

específica, viscosidade dinâmica, capacidade calorífica e condutibilidade térmica. A gama de

validade destas expressões está para o intervalo de temperaturas correspondente a 0℃ <T<100℃. A

figura 3.6 Mostra a evolução da massa específica com o aumento da temperatura e o efeito do

aumento da humidade relativa. Desta forma é fácil inferir que a variação torna-se maior para

temperaturas maiores que 60°C. Para as outras propriedades esta interpretação pode ser

considerada válida também.

64

Figura A.1 - Exemplo da variação da massa específica com temperatura e humidade relativa -

Tsilingiris (2007)

Correlações ar seco:

Viscosidade dinâmica (Ns/m2 x 10

6)

�� = 0 + _G + BGB + «G« + ³G³ (139)

• 0 = −9,8601 × 10²_ • _ = 9,080125 × 10²B • B = −1,17635575 × 10²³

• « = 1,2349703 × 10²· • ³ = −5,791299 × 10²__

Condutibilidade térmica (W/mK)

�� = 0 + _G + BGB + «G« + ³G³ + ¹G¹ (140)

• 0 = −2,276501 × 10²« • _ = 1,2598485 × 10²³ • B = −1,4815235 × 10²·

• « = 1,73550646 × 10²_0 • ³ = −1,066657 × 10²_« • ¹ = 2,47663035 × 10²_·

Capacidade calorífica (kJ/KgK)

&�� = 0 + _G + BGB + «G« + ³G³ (141)

• 0 = 0,103409 × 10

• _ = −0,284887 × 10²« • B = 0,7816818 × 10²¸

• « = −0,4970786 × 10²Á • ³ = 0,1077024 × 10²_B

Correlações vapor de água:

Viscosidade dinâmica (Ns/m2 x 10

6) (T em °C)

65

�� = 0 + _G (142)

• 0 = 8,058131868 × 10_ • _ = 4,00549451 × 10²_

Condutibilidade térmica (W/mK x 103) (T em °C)

�� = 0 + _G + BGB (143)

• 0 = 1,761758242 × 10_ • _ = 5,558941059 × 10²B • B = 1,663336663 × 10²³

Capacidade calorífica (kJ/KgK) (T em °C)

&�� = 0 + _G + BGB (144)

• 0 = 1,86910989

• _ = −2,578421578 × 10²³ • B = 1,941058941 × 10²¹

Anexo 3. Instruções

Assim sendo, a folha “Solve” contém o sistema de equações que está dividido em coeficiente,

variável (ou função) e parte constante da equação. No caso de a função ser personalizada ou, UDF

(User-Defined Function), é necessário que a função seja escrita como se fosse para introduzir no

VBATM (sem “;” só “,”) uma vez que é aplicada a função Evaluate do VBA

TM. Este facto pode ser um

fator que aumenta a complexidade do uso do método mas no entanto é possível escrever a função

normalmente numa célula e depois apenas usar a função de substituir do ExcelTM para deixar a

função só com vírgulas e sem o símbolo de igual (“=”) e pronta a substituir certos valores que são

incógnitas por “X_?” em que “?” representa um número inteiro maior ou igual que 1. Já os coeficientes

respetivos às incógnitas/funções bem como a parte constante da equação podem ser introduzidos

normalmente.

Também para garantir que o VBA não vai procurar e introduzir valores da folha de cálculo errada

é necessário estipular a que folha cada equação pertence. No caso de existirem funções de

diferentes folhas numa equação a função alheia à folha estipulada terá que ver as suas referências às

células modificadas ao acrescentar “Sheet!” às referências das células.

Desta maneira é possível diminuir o tempo de acesso à função Evaluate em metade segundo

Williams (2011). Para além destes inputs é possível decidir o incremento ℎKL para as diferenças

finitas, a tolerância, que está ligada à norma euclidiana de 7=�$, o número de iterações máximas (que

no caso do uso da condição de razão de convergência corresponde às iterações globais), bem como

se o método usa qualquer uma das condições acima descritas. Durante a implementação verificou-se

que a utilização do tipo double aumentou a precisão dos resultados embora à custa de alocação de

mais memória.

66

Como o método foi desenhado para resolver sistemas de equações o mesmo não está preparado

para resolver uma só equação (N=1). Note ainda que a introduzir valores constantes nas UDF é

necessário que os valores sejam introduzidos seguindo as normas americanas, ou seja, em vez de

vírgula usar o ponto para definir as casas decimais. Este passo é absolutamente crucial.

A título de exemplo foi acrescentada a figura A.1 de modo a demonstrar o programa e

nomeadamente o sistema de equações e a maneira como está estruturado para ser resolvido.

Figura A.2 - Exemplo da folha "Solve" em configuração 0

Apresentam-se assim as equações que constam no sistema de equações não lineares. É de

salientar que por cada permutador com condensação (i.e. condensador) adicionado, existe mais duas

equações iguais mas relativas a temperaturas diferentes.

äåååååæåååååç

Aa %����®ℎ%_ − ℎ%B± = �0=O0 − ;0$Aa %����®ℎ%B − ℎ%³± =í ��%«

�þ_ ÐσT%«³ − O�Ñ + �_ℎ��(F=G%« − GN�_$GN�_ = ;_ − O_3�~ ����B + G5_ + G5B2Aa 5=ℎ5« − ℎ5B$ = �_=;_ − O_$ + �_ℎ��(F=G%« − GN�_$ − p���,�����Aa 5=ℎ5B − ℎ5_$ = =;B − OB$ �B + c®G%³ − G5_±�Aa %����®ℎ%³ − ℎ%¹± = c®G%³ − G5_±�Aa 5=ℎ50D − ℎ50�$ = c∆GeWAa %����,�ℎ%� −Aa %����,9ℎ%9 −A�a ℎ� = c∆GeW

(145)

67

Anexo 4. Ensaios Caldeiras

Dados da caldeira relevantes para correlações. Note que não serão disponibilizados todos os dados

por razões de confidencialidade.

• ���F�� = 320 4/A]

• �£� = 180 4/A]

• *d/*e = 0,897

• 8�/)� = 0,140

• ��/*J = 4,000

• *J/8� = 0,769

Dados das condições de teste. Repare que a partir do ensaio #1 apenas apareceram as

modificações/alterações feitas de ensaio para ensaio, i.e. apenas apareceram nos ensaios

alterações feitas com base no ensaio anterior.

Ensaio #1 (Configuração 3)

Tabela A.11 - Correlações usadas no ensaio #1

Correlação

Permutador de calor primário: Lado água Delta Wing (2011)

Permutador de calor primário: Lado gases Youn (2007)

Condensador: Lado água Wire coil (2011)

Condensador: Lado gases Briggs e Young (1963)

• G�(F = 20,94℃

• G5« = 53,72℃

• G%º = 52,4℃

• fÙÝU = 7,65%

• Aa 5 = 15,99 v/Aj

• n5 = 2,08 ��l • � = 93,07%

• ��������� = 93,07%

• ü��(���� = 49,59 �4

• üú��� = 46,15 �4

• ΔG = 41,37℃

• ΔGNNNN = 33,04℃

• Tipo de gás: Metano (G20)

Ensaio #2 (Configuração 3)

Alterações feitas em relação ao ensaio #1.

• ü��(���� = 58,28 �4

• fÙÝU = 8,45%

• G%º = 50℃

• G%¹ = 250℃

Ensaio #3 (Configuração 3)

Alterações feitas em relação ao ensaio #2.

• ü��(���� = 5,82 �4

• fÙÝU = 1,5%

• G%º = 26℃

• G%¹ = 60℃

Tabela A.12 - Correlações alteradas para o ensaio #3

Correlação

Condensador: Lado gases Briggs e Young (1963) e Modelo1/Modelo 2

Ensaio #4 (Configuração 4)

68

Alterações feitas em relação ao ensaio #3.

• fÙÝU = 1,5%

• G%º = 23℃

• G%¹ = 60℃

Tabela A.13 - Correlações alteradas para o ensaio #4

Correlação

Condensador: Lado gases Briggs e Young (1963) e Modelo 2

Ensaio #4.1 (Configuração 4)

Alterações feitas em relação ao ensaio #3.

• fÙÝU = 1,6%

Ensaio #5 (Configuração 2)

Alterações feitas em relação ao ensaio #4.

• ü��(���� = 58,28 �4

• fÙÝU = 8,45%

• G%º = 52℃

• G%¹ = 255℃

Tabela A.14 - Correlações alteradas para o ensaio #5

Correlação

Condensador: Lado água Gnielinski (1995)

Ensaio #6

A partir do ensaio #6 serão testadas as caldeiras sem condensação pelo que os dados que não

estiverem presentes a partir do ensaio #6 podem ser assumidos como iguais aos do ensaio #1.

Qualquer dado alterado desde o ensaio #6 é também indicado. As novas dimensões não serão

apresentadas no entanto.

• G�(F = 25℃

• G5« = 60℃

• G50 = 20℃

• fÙÝU = 8,50%

• Aa 5 = 0,115 ]�/!

• Número de alhetas: 63

• � = 86,9%

• ü��(���� = 22,2�4

• ΔG = 40℃

• G%º = 202℃

• Tipo de gás: Propano (G30)

Ensaio #7

Alterações feitas em relação ao ensaio #6.

• fÙÝU = 5,65%

• Aa 5 = 0,136 ]�/!

• � = 86,9%

• Número de alhetas: 73

• ü��(���� = 27 �4

• G%º = 195℃

• Tipo de gás: Metano (G20)

Ensaio #8

Alterações feitas em relação ao ensaio #7.

• fÙÝU = 6,95% • Aa 5 = 0,170 ]�/!

69

• � = 86,7%

• Número de alhetas: 93

• ü��(���� = 32,8 �4

• G%º = 223℃

• Tipo de gás: Metano (G20)