Upload
trinhhuong
View
227
Download
0
Embed Size (px)
Citation preview
UNIVERSIDADE TECNOLÓGICA FEDERAL DO PARANÁ
CAMPUS CURITIBA DEPARTAMENTO DE PESQUISA E PÓS GRADUAÇÃO
PROGRAMA DE PÓS GRADUAÇÃO EM ENGENHARIA MECÂNICA E DE MATERIAIS – PPGEM
HENDY TISSERANT RODRIGUES
SIMULAÇÃO NUMÉRICA DO ESCOAMENTO BIFÁSICO GÁS-LÍQUIDO NO PADRÃO DE
GOLFADAS UTILIZANDO UM MODELO LAGRANGEANO DE SEGUIMENTO DE PISTÕES
CURITIBA OUTUBRO - 2009
HENDY TISSERANT RODRIGUES
SIMULAÇÃO NUMÉRICA DO ESCOAMENTO BIFÁSICO GÁS-LÍQUIDO NO PADRÃO DE
GOLFADAS UTILIZANDO UM MODELO LAGRANGEANO DE SEGUIMENTO DE PISTÕES
Dissertação apresentada como requisito parcial
à obtenção do título de Mestre em Engenharia,
do Programa de Pós-Graduação em
Engenharia Mecânica e de Materiais, Área de
Concentração em Engenharia Térmica, do
Departamento de Pesquisa e Pós-Graduação,
do Campus de Curitiba, da UTFPR.
Orientador: Prof. Rigoberto E. M. Morales, Dr.
Co-orientador: Prof. Eugênio S. Rosa, PhD.
CURITIBA OUTUBRO - 2009
TERMO DE APROVAÇÃO
HENDY TISSERANT RODRIGUES
SIMULAÇÃO NUMÉRICA DO ESCOAMENTO BIFÁSICO GÁS-LÍQUIDO NO PADRÃO DE
GOLFADAS UTILIZANDO UM MODELO LAGRANGEANO DE SEGUIMENTO DE PISTÕES
Esta Dissertação foi julgada para a obtenção do título de Mestre em Engenharia, e
aprovada em sua versão final pelo Programa de Pós-Graduação em Engenharia
Mecânica e de Materiais, área de concentração em Engenharia Térmica.
_________________________________
Prof. Giuseppe Pintaúde, Dr.Eng.
Coordenador de curso
Banca Examinadora
______________________________ ______________________________
Prof. Rigoberto E. M. Morales, Dr. Eng. Prof. Eugênio Spanó Rosa, PhD
PPGEM/UTFPR DE/FEM/UNICAMP
______________________________ ______________________________
José Roberto Fagundes Netto, PhD Prof. Cezar Otaviano R. Negrão, PhD
PETROBRAS/CENPES PPGEM/UTFPR
Curitiba, 3 de Fevereiro de 2010
iii
AGRADECIMENTOS
Agradeço ao LACIT/PPGEM/UTFPR e ao 2PFG/FEM/UNICAMP pela
possibilidade de realização desse trabalho. Agradeço à Petrobras e ANP pelo
suporte técnico e financeiro para o desenvolvimento do tema.
Agradeço aos professores Rigoberto E. M. Morales, Eugênio S. Rosa e
Ricardo A. Mazza pelo apoio, incentivo, troca de informações e amizade ao longo
dos últimos anos.
Dedico esse trabalho à minha esposa Joyce, pelo apoio e carinho
incondicional, e aos meus pais, pelo cuidado e preocupação em todos os momentos.
Agradeço a Deus por tudo que tem proporcionado em minha vida.
iv
RODRIGUES, Hendy Tisserant, Simulação numérica do escoamento bifásico gás-
líquido no padrão de golfadas utilizando um modelo lagrangeano de seguimento de
pistões, 2009, Dissertação (Mestrado em Engenharia) - Programa de Pós-graduação
em Engenharia Mecânica e de Materiais, Universidade Tecnológica Federal do
Paraná, Curitiba, 180p.
RESUMO
No escoamento bifásico de líquido e gás em tubulações, tem-se com grande
freqüência a ocorrência do padrão de golfadas. Esse padrão de escoamento tem
como característica a sucessão intermitente no tempo e no espaço de duas
estruturas distintas: bolha alongada e pistão de líquido. O presente trabalho
apresenta a modelagem matemática lagrangeana unidimensional do escoamento
em golfadas, onde as equações de conservação da massa e quantidade de
movimento são aplicadas a cada bolha e pistão, tornando-os elementos
computacionais definidos por volumes de controle discretos que evoluem ao longo
da tubulação. O modelo leva em conta efeitos desprezados em trabalhos anteriores,
como a variação da fração de líquido no pistão e a intermitência intrínseca do
escoamento, obtida através da introdução aleatória de bolhas e pistões na entrada
da tubulação. As equações diferenciais obtidas na modelagem são discretizadas
através do método de diferenças finitas. O sistema de equações algébricas
resultantes da discretização é resolvido utilizando o algoritmo TDMA. São calculados
parâmetros característicos do escoamento em golfadas, como os comprimentos e
velocidades das bolhas e pistões e a queda de pressão. Essas variáveis são
monitoradas através de valores médios ou distribuições em determinados pontos da
tubulação e seguindo-se uma célula ao longo de sua passagem pela tubulação. Os
resultados numéricos obtidos são comparados a dados experimentais, fornecidos
pelo 2PFG/FEM/UNICAMP, para o escoamento de ar-água, ar-glicerina e N2-óleo
em um duto horizontal, inclinado e vertical.
Palavras-chave: Escoamento intermitente em golfadas, modelo de seguimento de
pistões, escoamento bifásico.
v
RODRIGUES, Hendy Tisserant, Simulação numérica do escoamento bifásico gás-
líquido no padrão de golfadas utilizando um modelo lagrangeano de seguimento de
pistões, 2009, Dissertação (Mestrado em Engenharia) - Programa de Pós-graduação
em Engenharia Mecânica e de Materiais, Universidade Tecnológica Federal do
Paraná, Curitiba, 180p.
ABSTRACT
In two-phase gas-liquid flow through pipelines the slug flow pattern frequently
occurs. This flow pattern is characterized by an intermittent, in space and time,
succession of two distinct structures: elongated bubble and liquid slug. This work
presents the one-dimensional Lagrangean mathematical model for slug flow. The
mass and momentum conservation equations are applied for each bubble and slug,
which become computational elements defined by discrete control volumes that
evolve through the pipe. The model takes into account some effects neglected in
previous works such as the slug liquid fraction variation and the inherent
intermittence, obtained by the random introduction of bubbles and slugs at the pipe
inlet. The differential equations obtained in the mathematical model are discretized
using the finite difference method and the resulting linear system is solved with the
TDMA algorithm. Typical parameters of slug flow are calculated, such as the bubbles
and slugs lengths and velocities and pressure drop. These variables are monitored
through its mean values or distributions in determined locations along the pipe, or by
the following of one bubble passage through the pipe. Numerical results are
compared against experimental results from 2PFG/FEM/UNICAMP for air-water, air-
glycerin and nitrogen-oil flows in horizontal, inclined and vertical pipes.
Keywords: Slug flow, slug tracking model, two-phase flow.
vi
SUMÁRIO
AGRADECIMENTOS..............................................................................................................................III
RESUMO................................................................................................................................................ IV
ABSTRACT............................................................................................................................................. V
SUMÁRIO............................................................................................................................................... VI
LISTA DE FIGURAS.............................................................................................................................. IX
LISTA DE TABELAS ............................................................................................................................ XII
LISTA DE SÍMBOLOS......................................................................................................................... XIII
1 INTRODUÇÃO ............................................................................................................................... 1
1.1 OBJETIVOS .............................................................................................................................. 5 1.2 JUSTIFICATIVAS........................................................................................................................ 6 1.3 ESTRUTURA DO TRABALHO........................................................................................................ 6
2 REVISÃO BIBLIOGRÁFICA.......................................................................................................... 8
2.1 ESTUDOS ANTERIORES SOBRE ESCOAMENTO EM GOLFADAS ....................................................... 8 2.1.1 Modelos estacionários..................................................................................................... 13 2.1.2 Modelos transientes......................................................................................................... 15
2.2 DEFINIÇÃO DE PARÂMETROS IMPORTANTES ASSOCIADOS AO ESCOAMENTO BIFÁSICO EM GOLFADAS
19 2.3 EQUAÇÕES DE FECHAMENTO .................................................................................................. 22
2.3.1 Velocidade de translação das bolhas de gás.................................................................. 23 2.3.2 Fração de líquido no pistão ............................................................................................. 24 2.3.3 Freqüência da célula unitária .......................................................................................... 26
3 MODELAGEM MATEMÁTICA..................................................................................................... 28
3.1 EQUAÇÕES DE CONSERVAÇÃO APLICADAS A UMA CÉLULA GENÉRICA ......................................... 28 3.1.1 Conservação da massa na célula ................................................................................... 30 3.1.2 Conservação da quantidade de movimento no pistão .................................................... 40 3.1.3 Acoplamentos entre o pistão e as bolhas adjacentes ..................................................... 43
3.2 EQUAÇÕES AUXILIARES........................................................................................................... 49 3.2.1 Deslocamento da frente da bolha.................................................................................... 50 3.2.2 Deslocamento da traseira da bolha................................................................................. 51 3.2.3 Definição do atrito............................................................................................................ 51 3.2.4 Velocidade do filme de líquido......................................................................................... 53
3.3 MODELAGEM DAS SINGULARIDADES......................................................................................... 54
vii
3.3.1 Coalescência de bolhas .................................................................................................. 54
4 MÉTODO DE SOLUÇÃO NUMÉRICA ........................................................................................ 57
4.1 DISCRETIZAÇÃO DAS EQUAÇÕES ACOPLADAS........................................................................... 57 4.2 DISCRETIZAÇÃO DAS EQUAÇÕES AUXILIARES ........................................................................... 64 4.3 CONDIÇÕES INICIAIS E DE CONTORNO...................................................................................... 65 4.4 PROCESSO DE ENTRADA DE BOLHAS E PISTÕES NO DOMÍNIO DE CÁLCULO ................................. 67
4.4.1 Cálculo das variáveis na entrada de uma nova célula .................................................... 69 4.5 PROCESSO DE SAÍDA DE BOLHAS E PISTÕES DO DOMÍNIO DE CÁLCULO....................................... 79 4.6 PROCESSO DE INÍCIO DA SIMULAÇÃO ....................................................................................... 82 4.7 ALGORITMO DE SIMULAÇÃO..................................................................................................... 83 4.8 SONDAS VIRTUAIS .................................................................................................................. 86 4.9 SIMULAÇÕES PRELIMINARES ................................................................................................... 87
4.9.1 Análise do passo de tempo ............................................................................................. 87 4.9.2 Análise das constantes do fator de esteira ..................................................................... 90 4.9.3 Análise da condição de contorno aleatória versus condição periódica........................... 94
5 RESULTADOS PARA O ESCOAMENTO HORIZONTAL .......................................................... 98
5.1 CONFIGURAÇÕES EXPERIMENTAIS........................................................................................... 98 5.2 COMPARAÇÃO ENTRE AS SIMULAÇÕES E OS RESULTADOS EXPERIMENTAIS ................................ 99
5.2.1 Resultados para o escoamento de ar e água ................................................................. 99 5.2.2 Resultados para o escoamento de ar e glicerina .......................................................... 106 5.2.3 Resultados para o escoamento de N2 e óleo SAE 20-50 ............................................. 111
6 RESULTADOS PARA O ESCOAMENTO VERTICAL E INCLINADO ..................................... 116
6.1 CONFIGURAÇÕES EXPERIMENTAIS......................................................................................... 116 6.2 COMPARAÇÃO ENTRE AS SIMULAÇÕES E RESULTADOS EXPERIMENTAIS ................................... 117
6.2.1 Resultados para o escoamento de ar e água na vertical .............................................. 117 6.2.2 Resultados para o escoamento de N2 e óleo SAE 20-50 na vertical ............................ 123 6.2.3 Resultados para o escoamento de ar e água inclinado ................................................ 128
7 CONCLUSÕES .......................................................................................................................... 133
REFERÊNCIAS BIBLIOGRÁFICAS................................................................................................... 136
APÊNDICE A – RESULTADOS PARA DIVERSAS CONDIÇÕES DE ESCOAMENTO................... 142
AR E ÁGUA NA HORIZONTAL................................................................................................................ 142 AR E GLICERINA NA HORIZONTAL ........................................................................................................ 149 NITROGÊNIO E ÓLEO SAE 20-50 NA HORIZONTAL ............................................................................... 156 AR E ÁGUA NA VERTICAL .................................................................................................................... 160 NITROGÊNIO E ÓLEO SAE 20-50 NA VERTICAL .................................................................................... 164 AR E ÁGUA INCLINADO ....................................................................................................................... 166
viii
APÊNDICE B – RESUMO DE CORRELAÇÕES PARA CÁLCULO DAS VARIÁVEIS DE FECHAMENTO ................................................................................................................................... 170
APÊNDICE C – MODELOS DE DOIS FLUIDOS E DRIFT FLUX...................................................... 176
MODELO DE DOIS FLUIDOS (TWO-FLUID MODEL) .................................................................................. 176 MODELO DE DRIFT FLUX..................................................................................................................... 178
ix
LISTA DE FIGURAS
FIGURA 1.1 – PADRÕES DE ESCOAMENTO BIFÁSICO LÍQUIDO-GÁS EM TUBULAÇÕES HORIZONTAIS (LÍQUIDO EM
BRANCO E GÁS EM CINZA): A) BOLHAS DISPERSAS, B) GOLFADAS E C) ESTRATIFICADO .................................. 2 FIGURA 1.2 – PADRÕES DE ESCOAMENTO BIFÁSICO LÍQUIDO-GÁS EM TUBULAÇÕES VERTICAIS (LÍQUIDO EM
BRANCO E GÁS EM CINZA): A) BOLHAS DISPERSAS, B) GOLFADAS, C) AGITADO E D) ANULAR ......................... 2 FIGURA 1.3 – REPRESENTAÇÃO DAS PRINCIPAIS VARIÁVEIS DO ESCOAMENTO EM GOLFADAS.................................. 3 FIGURA 2.1 – REPRESENTAÇÃO DO ESCOAMENTO EM GOLFADAS A) HORIZONTAL OU INCLINADO E B) VERTICAL.... 9 FIGURA 2.2 – REPRESENTAÇÃO DAS LINHAS DE CORRENTE NO LÍQUIDO EM UM REFERENCIAL A) SE MOVENDO COM
A BOLHA E B) ESTACIONÁRIO, E C) REPRESENTAÇÃO DOS PERFIS DE VELOCIDADE NO LÍQUIDO .................... 10 FIGURA 2.3 – LINHAS DE CORRENTE PARA O ESCOAMENTO VERTICAL ................................................................... 12 FIGURA 2.4 – REPRESENTAÇÃO DA QUEDA DE PRESSÃO AO LONGO DO ESCOAMENTO EM GOLFADAS .................... 13 FIGURA 2.5 – VOLUME DE CONTROLE AO LONGO DE UM TRECHO DE TUBULAÇÃO ................................................. 22 FIGURA 3.1 – REPRESENTAÇÃO DA J-ÉSIMA CÉLULA GENÉRICA DO ESCOAMENTO EM GOLFADAS .......................... 29 FIGURA 3.2 – VOLUME DE CONTROLE DEFINIDO PELAS FRONTEIRAS DO PISTÃO .................................................... 31 FIGURA 3.3 – VOLUME DE CONTROLE DEFINIDO PELAS FRONTEIRAS DO FILME E BOLHA ALONGADA..................... 35 FIGURA 3.4 – FORÇAS QUE ATUAM NO VOLUME DE CONTROLE DEFINIDO PELO PISTÃO.......................................... 43 FIGURA 3.5 – REPRESENTAÇÃO DO VOLUME DE CONTROLE NO ACOPLAMENTO DA TRASEIRA DO PISTÃO COM A
FRENTE DA BOLHA ........................................................................................................................................ 44 FIGURA 3.6 – REPRESENTAÇÃO DO VOLUME DE CONTROLE NO ACOPLAMENTO ENTRE A FRENTE DO PISTÃO E A
TRASEIRA DA BOLHA..................................................................................................................................... 45 FIGURA 3.7 – REPRESENTAÇÃO DO VOLUME DE CONTROLE NO ACOPLAMENTO ENTRE A FRENTE DO PISTÃO E A
FRENTE DA BOLHA ........................................................................................................................................ 47 FIGURA 3.8 – REPRESENTAÇÃO DO VOLUME DE CONTROLE NA INTERFACE ENTRE O PISTÃO E A BOLHA QUE O
SEGUE ........................................................................................................................................................... 53 FIGURA 3.9 – COALESCÊNCIA DE BOLHAS .............................................................................................................. 56 FIGURA 4.1 – REPRESENTAÇÃO ESQUEMÁTICA DO ESCOAMENTO AO LONGO DO TUBO EM DIFERENTES INSTANTES
DE TEMPO E A DEFINIÇÃO DAS CONDIÇÕES DE CONTORNO E INICIAL............................................................. 66 FIGURA 4.2 – REPRESENTAÇÃO ESQUEMÁTICA DO PROCESSO DE ENTRADA DE PISTÃO DE LÍQUIDO E BOLHAS DE
GÁS NO DOMÍNIO COMPUTACIONAL .............................................................................................................. 68 FIGURA 4.3 – REPRESENTAÇÃO DA GEOMETRIA UTILIZADA NA INTEGRAÇÃO DO PERFIL DA BOLHA ...................... 72 FIGURA 4.4 – FLUXOGRAMA PARA CÁLCULO DOS PARÂMETROS DA NOVA CÉLULA QUE SERÁ INSERIDA NA
TUBULAÇÃO.................................................................................................................................................. 74 FIGURA 4.5 – REPRESENTAÇÃO DE DISTRIBUIÇÕES EXPERIMENTAIS ...................................................................... 75 FIGURA 4.6 – REPRESENTAÇÃO ESQUEMÁTICA DO PROCESSO DE SAÍDA DE BOLHAS E PISTÕES DO DOMÍNIO
COMPUTACIONAL.......................................................................................................................................... 82 FIGURA 4.7 – REPRESENTAÇÃO ESQUEMÁTICA DA CONDIÇÃO INICIAL................................................................... 83 FIGURA 4.8 – ALGORITMO GERAL DE SIMULAÇÃO ................................................................................................. 85
x
FIGURA 4.9 – VARIAÇÃO DO TEMPO DE SIMULAÇÃO E DO RESÍDUO RELATIVO EM FUNÇÃO DO PASSO DE TEMPO
UTILIZADO. ................................................................................................................................................... 90 FIGURA 4.10 – EVOLUÇÃO DOS COMPRIMENTOS DE BOLHAS E PISTÕES, TAXA DE COALESCÊNCIAS E VELOCIDADE
DA BOLHA AO LONGO DO TUBO PARA DIVERSOS VALORES DAS CONSTANTES DE ESTEIRA NO ESCOAMENTO DE
AR E ÁGUA .................................................................................................................................................... 92 FIGURA 4.11 – EVOLUÇÃO DOS COMPRIMENTOS DE BOLHAS E PISTÕES, TAXA DE COALESCÊNCIAS E VELOCIDADE
DA BOLHA AO LONGO DO TUBO PARA DIVERSOS VALORES DAS CONSTANTES DE ESTEIRA NO ESCOAMENTO DE
AR E GLICERINA ............................................................................................................................................ 93 FIGURA 4.12 – RESULTADOS PARA UMA CÉLULA SEGUIDA AO LONGO DE SUA PASSAGEM PELO TUBO COM AS
CONDIÇÕES DE CONTORNO PERIÓDICA E INTERMITENTE ............................................................................... 96 FIGURA 4.13 – TAXA DE COALESCÊNCIA AO LONGO DO TUBO PARA AS CONDIÇÕES DE CONTORNO PERIÓDICA E
INTERMITENTE .............................................................................................................................................. 96 FIGURA 4.14 – FUNÇÕES DENSIDADE DE PROBABILIDADE PARA OS COMPRIMENTOS DE BOLHAS E PISTÕES E PARA A
VELOCIDADE DA BOLHA EXPERIMENTAIS E SIMULADOS UTILIZANDO-SE AS CONDIÇÕES DE ENTRADA
PERIÓDICA E ALEATÓRIA............................................................................................................................... 97 FIGURA 5.1 – REPRESENTAÇÃO DAS SEÇÕES DE TESTE EXPERIMENTAL E NUMÉRICA PARA O ESCOAMENTO DE AR E
ÁGUA NA HORIZONTAL ............................................................................................................................... 100 FIGURA 5.2 – RESULTADOS MÉDIOS PARA O ESCOAMENTO HORIZONTAL DE AR E ÁGUA NA CONDIÇÃO A@W#2 102 FIGURA 5.3 – RESULTADOS MÉDIOS PARA O ESCOAMENTO HORIZONTAL DE AR E ÁGUA NA CONDIÇÃO A@W#4 102 FIGURA 5.4 – RESULTADO DAS FUNÇÕES DENSIDADE DE PROBABILIDADE PARA O ESCOAMENTO HORIZONTAL DE
AR E ÁGUA NA CONDIÇÃO A@W#2 ............................................................................................................ 104 FIGURA 5.5 - RESULTADO DAS FUNÇÕES DENSIDADE DE PROBABILIDADE PARA O ESCOAMENTO HORIZONTAL DE AR
E ÁGUA NA CONDIÇÃO A@W#4 ................................................................................................................. 105 FIGURA 5.6 – REPRESENTAÇÃO DAS SEÇÕES DE TESTE EXPERIMENTAL E NUMÉRICA PARA O ESCOAMENTO DE AR E
GLICERINA NA HORIZONTAL ....................................................................................................................... 106 FIGURA 5.7 – RESULTADOS MÉDIOS PARA O ESCOAMENTO HORIZONTAL DE AR E SOLUÇÃO DE GLICERINA NA
CONDIÇÃO A@G#2 .................................................................................................................................... 107 FIGURA 5.8 – RESULTADOS MÉDIOS PARA O ESCOAMENTO HORIZONTAL DE AR E SOLUÇÃO DE GLICERINA NA
CONDIÇÃO A@G#4 .................................................................................................................................... 107 FIGURA 5.9 – RESULTADO DAS FUNÇÕES DENSIDADE DE PROBABILIDADE PARA O ESCOAMENTO HORIZONTAL DE
AR E SOLUÇÃO DE GLICERINA NA CONDIÇÃO A@G#2 ................................................................................ 109 FIGURA 5.10 – RESULTADO DAS FUNÇÕES DENSIDADE DE PROBABILIDADE PARA O ESCOAMENTO HORIZONTAL DE
AR E SOLUÇÃO DE GLICERINA NA CONDIÇÃO A@G#4 ................................................................................ 110 FIGURA 5.11 – REPRESENTAÇÃO DAS SEÇÕES DE TESTE EXPERIMENTAL E NUMÉRICA PARA O ESCOAMENTO DE
NITROGÊNIO E ÓLEO SAE 20-50 NA HORIZONTAL....................................................................................... 111 FIGURA 5.12 - RESULTADOS MÉDIOS PARA O ESCOAMENTO HORIZONTAL DE NITROGÊNIO E ÓLEO SAE 20-50
N@O#1...................................................................................................................................................... 112 FIGURA 5.13 – RESULTADOS MÉDIOS PARA O ESCOAMENTO HORIZONTAL DE NITROGÊNIO E ÓLEO SAE 20-50
N@O#4...................................................................................................................................................... 113
xi
FIGURA 5.14 – RESULTADO DAS FUNÇÕES DENSIDADE DE PROBABILIDADE PARA O ESCOAMENTO HORIZONTAL DE
NITROGÊNIO E ÓLEO SAE 20-50 NA CONDIÇÃO N@O#1 ............................................................................ 114 FIGURA 5.15 – RESULTADO DAS FUNÇÕES DENSIDADE DE PROBABILIDADE PARA O ESCOAMENTO HORIZONTAL DE
NITROGÊNIO E ÓLEO SAE 20-50 NA CONDIÇÃO N@O#4 ............................................................................ 115 FIGURA 6.1 – REPRESENTAÇÃO DAS SEÇÕES DE TESTE EXPERIMENTAL E NUMÉRICA PARA O ESCOAMENTO DE AR E
ÁGUA NA VERTICAL .................................................................................................................................... 117 FIGURA 6.2 – RESULTADOS MÉDIOS PARA O ESCOAMENTO VERTICAL DE AR E ÁGUA NA CONDIÇÃO A@WV#2.. 119 FIGURA 6.3 – RESULTADOS MÉDIOS PARA O ESCOAMENTO VERTICAL DE AR E ÁGUA NA CONDIÇÃO A@WV#8.. 119 FIGURA 6.4 – RESULTADO DAS FUNÇÕES DENSIDADE DE PROBABILIDADE PARA O ESCOAMENTO VERTICAL DE AR E
ÁGUA NA CONDIÇÃO A@WV#2 ................................................................................................................. 121 FIGURA 6.5 – RESULTADO DAS FUNÇÕES DENSIDADE DE PROBABILIDADE PARA O ESCOAMENTO VERTICAL DE AR E
ÁGUA NA CONDIÇÃO A@WV#8 ................................................................................................................. 122 FIGURA 6.6 – REPRESENTAÇÃO DAS SEÇÕES DE TESTE EXPERIMENTAL E NUMÉRICA PARA O ESCOAMENTO DE
NITROGÊNIO E ÓLEO SAE 20-50 NA VERTICAL ........................................................................................... 123 FIGURA 6.7 – RESULTADOS MÉDIOS PARA O ESCOAMENTO VERTICAL DE NITROGÊNIO E ÓLEO SAE 20-50 NA
CONDIÇÃO N@OV#2 ................................................................................................................................. 124 FIGURA 6.8 – RESULTADOS MÉDIOS PARA O ESCOAMENTO VERTICAL DE NITROGÊNIO E ÓLEO SAE 20-50 NA
CONDIÇÃO N@OV#6 ................................................................................................................................. 125 FIGURA 6.9 – RESULTADO DAS FUNÇÕES DENSIDADE DE PROBABILIDADE PARA O ESCOAMENTO VERTICAL DE
NITROGÊNIO E ÓLEO SAE 20-50 NA CONDIÇÃO N@OV#2 ......................................................................... 126 FIGURA 6.10 – RESULTADO DAS FUNÇÕES DENSIDADE DE PROBABILIDADE PARA O ESCOAMENTO VERTICAL DE
NITROGÊNIO E ÓLEO SAE 20-50 NA CONDIÇÃO N@OV#6 ......................................................................... 127 FIGURA 6.11 – REPRESENTAÇÃO DAS SEÇÕES DE TESTE EXPERIMENTAL E NUMÉRICA PARA O ESCOAMENTO DE AR E
ÁGUA INCLINADO........................................................................................................................................ 128 FIGURA 6.12 RESULTADOS MÉDIOS PARA O ESCOAMENTO INCLINADO DE AR E ÁGUA NA CONDIÇÃO A@WI#2 .. 129 FIGURA 6.13 – RESULTADOS MÉDIOS PARA O ESCOAMENTO INCLINADO DE AR E ÁGUA NA CONDIÇÃO A@WI#5 129 FIGURA 6.14 – RESULTADO DAS FUNÇÕES DENSIDADE DE PROBABILIDADE PARA O ESCOAMENTO INCLINADO DE AR
E ÁGUA NA CONDIÇÃO A@WI#2 ................................................................................................................ 131 FIGURA 6.15 – RESULTADO DAS FUNÇÕES DENSIDADE DE PROBABILIDADE PARA O ESCOAMENTO INCLINADO DE AR
E ÁGUA NA CONDIÇÃO A@WI#5 ................................................................................................................ 132
xii
LISTA DE TABELAS
TABELA 4.1 – DEFINIÇÃO DOS PARÂMETROS UTILIZADOS NAS SIMULAÇÕES PRELIMINARES.................................. 87 TABELA 5.1 – DEFINIÇÃO DAS CONFIGURAÇÕES EXPERIMENTAIS .......................................................................... 98 TABELA 5.2 – DEFINIÇÃO DAS CONDIÇÕES DE CONTORNO PARA O ESCOAMENTO DE AR E ÁGUA ........................... 99 TABELA 5.3 – DEFINIÇÃO DAS CONDIÇÕES DE CONTORNO PARA O ESCOAMENTO DE AR E GLICERINA.................... 99 TABELA 5.4 – DEFINIÇÃO DAS CONDIÇÕES DE CONTORNO PARA O ESCOAMENTO DE NITROGÊNIO E ÓLEO ............. 99 TABELA 6.1 – DEFINIÇÃO DAS CONFIGURAÇÕES EXPERIMENTAIS ........................................................................ 116 TABELA 6.2 – CONDIÇÕES DE CONTORNO PARA AR E ÁGUA NA VERTICAL ........................................................... 116 TABELA 6.3 – CONDIÇÕES DE CONTORNO PARA NITROGÊNIO E ÓLEO SAE 20-50 ................................................ 116 TABELA 6.4 – CONDIÇÕES DE CONTORNO PARA AR E ÁGUA INCLINADO............................................................... 117
xiii
LISTA DE SÍMBOLOS
A Área
C Coeficiente de atrito
C0 Constante da velocidade da bolha
D Diâmetro
F Força
f Freqüência
Fr Número de Froude
g Aceleração da gravidade
h Constante de esteira
j Numeração das células
J Velocidade superficial
L Comprimento
M Massa
m Vazão mássica
n Número de bolhas dentro do tubo
P Pressão
Q Vazão volumétrica
R Fração volumétrica da fase
Re Número de Reynolds
S Perímetro
t Tempo
U Velocidade absoluta
x Coordenada da frente do pistão
y Coordenada da frente da bolha
β Fator de intermitência
Δ Variação
ρ Massa específica
μ Viscosidade dinâmica
σ Tensão superficial
xiv
τ Tensão cisalhante
θ Ângulo de inclinação do tubo com a horizontal
Sub-índices
B Região da bolha alongada
D Deslizamento
G Fase gasosa
L Fase líquida
S Região do pistão de líquido (slug)
T Translação da bolha alongada
U Unidade
Capítulo 1 Introdução
1
1 INTRODUÇÃO
Escoamentos multifásicos ocorrem com grande freqüência tanto na natureza
como em aplicações industriais. Exemplos de sua ocorrência na natureza estão no
transporte de sedimentos em rios e correntes marinhas, e movimentação de massas
térmicas na atmosfera terrestre. Em aplicações industriais, escoamentos multifásicos
podem ocorrer em geradores de vapor, condensadores e no transporte de misturas
em tubulações.
Um caso particular de escoamentos multifásicos são os escoamentos
bifásicos de líquido e gás, nos quais o escoamento pode se arranjar
geometricamente em diversos padrões, a depender das condições do escoamento,
como vazões das fases, configurações geométricas e propriedades dos fluidos. Nas
Figuras 1.1 e 1.2 são representados os padrões de ocorrência mais comum nos
escoamentos em dutos horizontais e verticais. Na Figura 1.1 são mostrados os
padrões: bolhas dispersas, golfadas (slug flow) e estratificado, enquanto na Figura
1.2 são mostrados os padrões: bolhas dispersas, golfadas (slug flow), agitado (churn
flow) e anular.
O escoamento estratificado ocorre em dutos horizontais sob determinadas
condições de vazão de líquido e gás. Aumentando-se a vazão de líquido o
escoamento se torna instável até transacionar para o escoamento em golfadas, que
permanece para uma grande faixa de vazões de gás. Com o aumento ainda maior
da vazão de líquido, a turbulência associada quebra as bolhas alongadas formando
bolhas dispersas.
No caso do escoamento em dutos verticais, para baixas vazões de gás ocorre
o escoamento em bolhas dispersas. Aumentando-se a vazão de gás as bolhas
começam a coalescer, formando as bolhas alongadas (escoamento em golfadas).
Em altas vazões de gás o escoamento passa por uma transição (também chamado
de escoamento agitado) até atingir o padrão anular.
Capítulo 1 Introdução
2
a)
b)
c)
Figura 1.1 – Padrões de escoamento bifásico líquido-gás em tubulações horizontais
(líquido em branco e gás em cinza): a) Bolhas dispersas, b) golfadas e c)
estratificado
a) b) c) d) Figura 1.2 – Padrões de escoamento bifásico líquido-gás em tubulações verticais
(líquido em branco e gás em cinza): a) Bolhas dispersas, b) golfadas, c) agitado e d)
anular
Dentre os padrões de escoamento apresentados, o escoamento em golfadas
ocorre em uma grande faixa de vazões de líquido e gás. Dessa forma, em diversas
atividades industriais, o correto tratamento desse padrão de escoamento é de suma
importância.
Capítulo 1 Introdução
3
O escoamento em golfadas é caracterizado pela sucessão de duas regiões
distintas: o pistão de líquido e a bolha alongada. O pistão de líquido é composto por
uma região com grande quantidade de líquido, e pode conter pequenas bolhas de
gás dispersas. Porém, mesmo com a presença de bolhas dispersas, o pistão ainda
se caracteriza como uma barreira entre as duas bolhas de gás adjacentes. Para o
escoamento horizontal, a região da bolha alongada é considerada um escoamento
estratificado, com líquido na região inferior e gás na região superior. No caso do
escoamento vertical, na região da bolha o escoamento é anular com líquido no
exterior e gás no interior. Tanto no escoamento horizontal como no vertical, a região
de líquido que escoa ao lado da bolha é chamada de filme de líquido, e é
considerada livre de gás disperso.
A Figura 1.3 apresenta uma célula unitária do escoamento em golfadas,
contendo uma bolha com comprimento BL e um pistão com comprimento SL . Para a
caracterização desse tipo de escoamento, é importante também quantificar as
velocidades do líquido no pistão, LSU , do líquido no filme, LBU , e de translação da
bolha, TU , além da pressão do gás no interior da bolha, GBP .
Figura 1.3 – Representação das principais variáveis do escoamento em golfadas
O escoamento em golfadas ocorre de maneira intermitente, pois a repetição
das estruturas bolha alongada e pistão de líquido não é periódica no tempo e nem
no espaço. Ou seja, todas as bolhas e pistões que passam por um referencial fixo
apresentam comprimentos e velocidades diferentes. Acredita-se que a intermitência
Capítulo 1 Introdução
4
é independente da forma de entrada das duas fases na tubulação. Ou seja, mesmo
que as fases sejam injetadas na tubulação de forma periódica ou contínua, o
escoamento se torna intermitente ao longo da tubulação.
Os primeiros estudos sobre o escoamento em golfadas foram realizados por
Wallis (1969), que definiu o conceito de célula unitária. Logo depois, Dukler e
Hubbard (1975) e Fernandes et al. (1983) realizaram estudos experimentais e
desenvolveram modelos simplificados para o escoamento em golfadas na horizontal
e vertical, respectivamente. Mais tarde, Taitel e Barnea (1990) propuseram um
modelo genérico para qualquer ângulo da tubulação. Esses trabalhos são
conhecidos como os modelos de estado estacionário, pois desconsideram a
intermitência do escoamento, e calculam os parâmetros como se as bolhas e pistões
se repetissem periodicamente no tempo e no espaço.
Mais recentemente, com o desenvolvimento computacional, modelos de
simulação numérica surgiram com o intuito de obter resultados mais realistas,
levando em conta a intermitência do escoamento. Uma classe de modelos surgiu
baseada nos modelos de dois fluidos e drift flux, que utilizam uma malha euleriana
ao longo da tubulação e a solução numérica das equações diferenciais de
conservação da massa e quantidade de movimento para cada fase. Esses modelos
apresentam um custo computacional muito alto e apresentam instabilidades na
solução das equações de acordo com o arranjo das fases.
Outra classe de modelos que surgiu recentemente são os modelos
lagrangeanos de seguimento de pistões (slug tracking). Nessa classe de modelos,
são obtidas equações para volumes de controle definidos geometricamente pelos
pistões e bolhas, e esses volumes de controle são seguidos ao longo da simulação.
Como os volumes de controle têm a mesma ordem de grandeza das bolhas e
pistões, a simulação apresenta um menor custo computacional. Existem diversos
modelos de seguimento de pistões na literatura, iniciando-se por Barnea e Taitel
(1993), depois Zheng et al. (1994), Nydal e Banerjee (1995), Taitel e Barnea (2000),
Grenier (1997), Franklin e Rosa (2004), Ujang et al. (2006) e Rodrigues (2006).
Dentre os diversos modelos de seguimento de pistões propostos na literatura,
são utilizadas diferentes simplificações na modelagem. Efeitos importantes, como a
consideração do pistão de líquido aerado e a variação da fração de líquido no pistão,
ou a característica intermitente do escoamento, não são tratados com muita
Capítulo 1 Introdução
5
relevância. Além disso, não são encontrados modelos para simulação do
escoamento inclinado ou vertical, e sua comparação com resultados experimentais.
1.1 Objetivos
O objetivo do presente trabalho é o desenvolvimento de um modelo
unidimensional utilizando o método de seguimento de pistões (slug tracking) com a
finalidade de simular o escoamento bifásico em golfadas em dutos horizontais,
inclinados e verticais.
No desenvolvimento do modelo foram considerados os efeitos de inércia do
pistão de líquido, troca de quantidade de movimento entre os pistões e as bolhas,
efeito do filme de líquido na variação da quantidade de movimento em uma célula, e
efeito da variação da fração de líquido no pistão.
O modelo apresentado é genérico para qualquer ângulo de inclinação da
tubulação. As equações diferenciais obtidas na modelagem matemática são
discretizadas através do método de diferenças finitas com o esquema semi-implícito
de Crank-Nicholson, gerando um sistema linear que é resolvido numericamente a
cada passo de tempo. O modelo foi implementado em um programa computacional
na linguagem FORTRAN.
A intermitência do escoamento é introduzida no modelo através da condição
de entrada das bolhas e pistões na tubulação. Sempre que uma nova célula vai
entrar na tubulação, os parâmetros dessa célula são sorteados de listas de valores
independentes para cada variável. Dessa forma, cada nova célula apresenta valores
diferentes para cada variável.
Como dados de saída das simulações, tem-se o monitoramento das variáveis
importantes em algumas posições ao longo da tubulação, ou o acompanhamento de
uma célula ao longo de sua passagem pela tubulação. Isso permite a avaliação dos
valores médios ou as distribuições das variáveis. Os resultados obtidos são
comparados com dados experimentais para escoamentos horizontais de ar e água,
ar e glicerina e nitrogênio e óleo, verticais de ar e água e nitrogênio e óleo, e
inclinados de ar e água.
Capítulo 1 Introdução
6
1.2 Justificativas
O escoamento no padrão de golfadas ocorre em diversos tipos de aplicações
industriais. Uma dessas aplicações em que o escoamento em golfadas tem grande
ocorrência é a indústria petrolífera. Esse padrão de escoamento é característico nas
linhas de produção, que transportam o petróleo desde o poço até a unidade
estacionária de produção (UEP), onde as fases são separadas e processadas.
A variação intrínseca de vazão e pressão do escoamento em golfadas
interfere na planta de processamento, pois ora tem-se a chegada de apenas líquido,
e ora tem-se a chegada de grande quantidade de gás, o que afeta os parâmetros no
dimensionamento da planta. Além disso, a freqüência de alternância das regiões do
pistão de líquido e da bolha alongada em determinada região da linha de produção
causa vibrações. A freqüência dessas vibrações deve ser monitorada, e deve-se
atentar para o risco de ressonância e possíveis danos. Devido a esses fatores,
torna-se importante a modelagem do escoamento em golfadas e o entendimento dos
fenômenos que ocorrem, de forma a prever situações e se evitar danos ou perda de
produção.
Do ponto de vista acadêmico, este trabalho faz parte de um esforço realizado
pelo LACIT/UTFPR e 2PFG/FEM/UNICAMP, financiado pela Petrobras, no sentido
de estudar e caracterizar o escoamento bifásico de líquido e gás. Nos últimos anos,
tem-se realizado estudos experimentais e numéricos com essa finalidade. O modelo
aqui apresentado deve ser implementado em um programa computacional dando
origem a um simulador de escoamento bifásico, com interface amigável ao usuário,
de forma que possa ser utilizado por engenheiros de campo em análises de
elevação e escoamento.
1.3 Estrutura do trabalho
O presente trabalho está dividido em sete capítulos. O primeiro capítulo
apresenta a introdução ao tema, bem como os objetivos e justificativas do trabalho.
No segundo capítulo é apresentada uma revisão bibliográfica sobre o escoamento
em golfadas. São apresentados alguns conceitos desse padrão de escoamento,
Capítulo 1 Introdução
7
modelos propostos para a simulação do escoamento em golfadas e a revisão das
equações de fechamento para a velocidade de translação das bolhas, fração de
líquido no pistão e freqüência da célula unitária. O apêndice C apresenta os modelos
numéricos de dois fluidos e drift flux, e o apêndice B apresenta um resumo das
correlações existentes na literatura para o cálculo das variáveis de fechamento.
O terceiro capítulo apresenta a modelagem matemática unidimensional do
escoamento em golfadas. São aplicadas as equações de conservação da massa e
quantidade de movimento a uma célula genérica do escoamento, e, ao final do
desenvolvimento, são obtidas duas equações diferenciais para a velocidade do
pistão e para a pressão do gás na bolha. Além disso, são apresentadas equações
auxiliares, utilizadas para descrever o movimento das bolhas e pistões a partir das
velocidades e pressões.
No capítulo 4 é apresentado o método de solução numérica, com a
discretização das equações diferenciais e implementação em um sistema linear. É
apresentado também o algoritmo de simulação e a definição das condições iniciais,
de contorno, e de entrada e saída de bolhas e pistões na simulação. São
apresentadas também simulações preliminares com o objetivo de se definir
parâmetros para a realização das simulações, como o passo de tempo, a definição
das constantes do fator de esteira e a comparação entre as condições de contorno
periódica e aleatória.
Os capítulos 5 e 6 apresentam a comparação entre resultados experimentais
e numéricos para diversas condições de escoamento. O capítulo 5 apresenta os
resultados para escoamento horizontal de ar e água, ar e glicerina e nitrogênio e
óleo SAE 20-50, enquanto o capítulo 6 apresenta resultados para escoamento
vertical de ar e água e nitrogênio e óleo SAE 20-50, e resultados para o escoamento
inclinado de ar e água. Nesses capítulos são apresentados apenas alguns
resultados do banco de simulações, e os demais resultados são apresentados no
apêndice A.
O capítulo 7 apresenta as conclusões do trabalho e recomendações para
trabalhos futuros.
Capítulo 2 Revisão bibliográfica 8
2 REVISÃO BIBLIOGRÁFICA
Nesse capítulo é apresentada a revisão da bibliografia utilizada na
modelagem do escoamento em golfadas para a implementação em um modelo de
seguimento de pistões (slug tracking). Inicialmente, são apresentados estudos que
evidenciam alguns conceitos e características importantes do escoamento em
golfadas. Em seguida, são revisados os modelos para se calcular parâmetros do
escoamento em golfadas, como os modelos de estado estacionário e os modelos de
simulação numérica transientes. Na seqüência são desenvolvidos parâmetros e
conceitos importantes para escoamentos bifásicos em geral. Por fim, são
apresentadas as equações geralmente utilizadas no cálculo das variáveis de
fechamento dos modelos de escoamento em golfadas: velocidade de translação das
bolhas de gás, fração de líquido no pistão e freqüência da célula unitária.
2.1 Estudos anteriores sobre escoamento em golfadas
Nessa seção serão apresentadas algumas características básicas do
escoamento no padrão de golfadas, que resultam de estudos experimentais e
numéricos apresentados na literatura.
A principal característica desse padrão de escoamento é a repetição
intermitente, no espaço e no tempo, de duas estruturas distintas: bolha alongada e
pistão de líquido. A bolha alongada escoa ao lado de um filme de líquido em uma
geometria anular (na vertical) ou estratificada (na horizontal). O pistão de líquido
possui pequena concentração de gás na forma de pequenas bolhas dispersas.
No escoamento vertical a bolha alongada possui simetria em relação ao eixo
central do tubo, enquanto no escoamento horizontal a bolha alongada é deslocada
para cima devido ao efeito da gravidade, e encosta-se à parede do tubo. No
escoamento inclinado as características são semelhantes às do escoamento
horizontal até determinado ângulo de inclinação crítico em que a bolha alongada
deixa de tocar a parede do tubo. Esse ângulo crítico é desconhecido, pois depende
de diversos parâmetros do escoamento (vazões, geometria, propriedades dos
Capítulo 2 Revisão bibliográfica 9
fluidos, condições de operação). Além disso, no escoamento vertical o pistão de
líquido é mais aerado do que no escoamento horizontal.
Imagens fotográficas do escoamento em golfadas são encontradas nos
trabalhos de Fagundes Netto (1999) para o escoamento horizontal e Talvy et al.
(2000) para o escoamento vertical. De forma ilustrativa, as imagens apresentadas
pelos autores são representadas nas Figuras 2.1a e 2.1b. Nota-se que a frente da
bolha alongada tem formato cilíndrico, pois a bolha tende a minimizar o arrasto de
forma. Na traseira da bolha percebe-se um formato mais complexo, com
desprendimento de pequenas bolhas e uma esteira.
a)
Bolha alongada
Pistão de líquido
b)
Figura 2.1 – Representação do escoamento em golfadas a) horizontal ou inclinado e
b) vertical
O estudo do movimento do líquido ao redor da bolha alongada é importante,
pois auxilia a entender como o escoamento em golfadas ocorre. Diversos autores
estudaram o movimento do líquido ao redor da bolha alongada, tanto
experimentalmente, com a medição das velocidades, como numericamente, através
da simulação do escoamento ao redor da bolha.
Mao e Dukler (1989) realizaram um estudo experimental e mediram diversos
parâmetros do escoamento vertical, como as velocidades e comprimentos de bolhas
Capítulo 2 Revisão bibliográfica 10
e pistões, frações de vazio e as tensões cisalhantes na parede. Bugg et al. (1998)
implementaram um modelo numérico bidimensional baseado em diferenças finitas
para determinar as linhas de corrente ao redor da bolha. Taha e Cui (2006) também
obtiveram as linhas de corrente ao redor de uma bolha alongada numericamente,
porém, utilizaram o programa comercial Fluent para uma análise tridimensional.
Polonsky et al. (1999), van Hout et al. (2002) e Nogueira et al. (2006) utilizaram a
técnica de PIV (particle image velocimetry) para a visualização do escoamento ao
redor da bolha em seus experimentos.
Alguns resultados apresentados por esses autores são esquematizados nas
Figuras 2.2 e 2.3. As Figuras 2.2a e 2.3a apresentam as linhas de corrente a partir
de um referencial que se move junto com a bolha alongada. O escoamento é suave
ao redor do nariz da bolha devido ao formato cilíndrico. No entanto, na traseira da
bolha forma-se uma esteira onde ocorrem recirculações que geram o
desprendimento de pequenas bolhas.
xy
xya)
b)
c)
Figura 2.2 – Representação das linhas de corrente no líquido em um referencial a)
se movendo com a bolha e b) estacionário, e c) representação dos perfis de
velocidade no líquido
As Figuras 2.2b e 2.3b apresentam as linhas de corrente a partir de um
referencial estacionário, e percebe-se que a bolha alongada empurra o líquido que
Capítulo 2 Revisão bibliográfica 11
está à sua frente. No escoamento vertical, o líquido que passa pelo filme é acelerado
para trás e, ao chegar ao pistão seguinte, se desacelera e muda de direção,
causando a recirculação. Além disso, no escoamento horizontal, a velocidade do
líquido no filme é sempre positiva, enquanto no escoamento vertical a velocidade do
filme pode ser negativa devido ao efeito gravitacional.
Nas Figuras 2.2c e 2.3c estão apresentados os perfis de velocidade do líquido
no filme e no pistão a partir de um referencial estacionário. O líquido em contato com
a parede é estacionário, mas na região do filme, as camadas superiores são
arrastadas pela interface com a bolha alongada de forma semelhante ao
escoamento de Couette. No caso do escoamento vertical, mesmo a bolha se
movendo para cima, na fronteira do gás com o filme de líquido a velocidade é
negativa. No pistão, o perfil de velocidades necessita de uma distância para se
recuperar da aceleração abrupta e atingir o perfil característico de regime
permanente. É importante ressaltar que os desenhos apresentados são apenas
ilustrativos.
Capítulo 2 Revisão bibliográfica 12
x
x
y
y
a) b) c)
Figura 2.3 – Linhas de corrente para o escoamento vertical
Outro parâmetro importante no escoamento em golfadas é o comportamento
da pressão ao longo da célula unitária. Dukler e Hubbard (1975) mediram a variação
temporal da pressão entre determinada posição no tubo e a saída no escoamento
horizontal. Pinto e Campos (1996) utilizaram um medidor de pressão diferencial em
uma posição do tubo no escoamento vertical.
Os resultados obtidos por Dukler e Hubbard (1975) e Pinto e Campos (1996)
são resumidos através da Figura 2.4, que apresenta as bolhas 1 e 2, com pressões
1P e 2P , separadas por um pistão. A pressão ao longo da bolha alongada é
constante, e, no pistão de líquido, há uma variação de pressão, SPΔ , devido ao atrito
nas paredes e força gravitacional (no escoamento vertical ou inclinado). No entanto,
logo no início do pistão existe uma variação de pressão mais acentuada. Dukler e
Hubbard (1975) chamaram essa variação de MixPΔ , que surge devido à mistura do
líquido que vem do filme e é incorporado ao pistão seguinte. O termo MixPΔ está
relacionado à expansão súbita de um jato de líquido, que é um fenômeno complexo
Capítulo 2 Revisão bibliográfica 13
devido à geometria em que a expansão ocorre. Mais tarde, Taitel e Barnea (1990b)
concluíram que a variação de pressão MixPΔ também pode ser calculada a partir das
forças gravitacionais e de atrito que atuam no filme de líquido.
1
P1
Pre
ssão
P2
x
ΔPS
LB LS
LMix
ΔPMix
P3
2 3
Figura 2.4 – Representação da queda de pressão ao longo do escoamento em
golfadas
Os estudos apresentados anteriormente tiveram como motivação o
conhecimento de parâmetros específicos do escoamento em golfadas, como queda
de pressão ou perfis de velocidades. Esses parâmetros são importantes, pois
servem de base para a correta modelagem do fenômeno. Outros trabalhos foram
realizados com o intuito de modelar o escoamento para o cálculo de variáveis
importantes para utilização em projetos de engenharia. A seguir são apresentados
alguns modelos com essa característica.
2.1.1 Modelos estacionários
Os primeiros modelos propostos para calcular parâmetros de interesse do
escoamento em golfadas são os modelos chamados de estado estacionário ou
modelos de célula unitária. Essa definição é devido à consideração de que o
escoamento é estacionário e periódico, ou seja, uma única célula (bolha e pistão) se
repete tanto no tempo como no espaço. Com essa simplificação, todos os cálculos
Capítulo 2 Revisão bibliográfica 14
são realizados para uma única célula e extrapolados para todo o comprimento do
tubo.
Um dos primeiros estudos foi realizado por Wallis (1969), que definiu o
conceito de célula unitária. O autor utilizou correlações existentes na época para
calcular a velocidade da bolha e propôs uma forma simples de se calcular a queda
de pressão gravitacional e por atrito.
Em seguida, Dukler e Hubbard (1975) realizaram um estudo experimental do
escoamento em golfadas na horizontal. Os autores apresentam características do
escoamento que são aceitas até hoje, como a queda de pressão devido à
aceleração do líquido que passa do filme para o pistão que o segue. A partir do
modelo proposto por eles pode-se calcular os parâmetros de interesse, como
comprimentos de bolha e pistão e queda de pressão na célula. No entanto, esse
modelo despreza a intermitência (aleatoriedade) do escoamento e necessita de
algumas equações adicionais de fechamento, como a freqüência da célula e a
fração de líquido no pistão.
Fernandes et al. (1983) realizaram um estudo do escoamento em golfadas na
vertical. O trabalho segue a linha proposta por Dukler e Hubbard (1975), porém inclui
algumas características importantes ao escoamento vertical. Os autores também
propõem uma rede de equações para a solução do escoamento, porém são
necessárias equações de fechamento.
Mais recentemente, Taitel e Barnea (1990a) apresentam um modelo
estacionário mais genérico para o escoamento em golfadas. Os autores apresentam
a modelagem genérica para o escoamento horizontal, vertical e inclinado. Além
disso, os autores propõem três modelos para levar em consideração o formato da
bolha alongada (também chamado de modelo de bolha), o que não havia sido feito
pelos autores anteriores. No entanto, mesmo sendo o trabalho mais completo sobre
a modelagem em estado estacionário, ainda são necessárias equações de
fechamento.
Os modelos de estado estacionário possuem grandes restrições e
simplificações. Porém, a sua resolução é muito simples (apenas equações
algébricas) e resulta no cálculo de valores médios “razoáveis” de propriedades
físicas do escoamento em golfadas, o que explica sua grande aceitação e utilização
na indústria.
Capítulo 2 Revisão bibliográfica 15
2.1.2 Modelos transientes
Com o desenvolvimento dos computadores, surgiram novas metodologias
para o cálculo dos parâmetros no escoamento bifásico durante processos
transientes. As três principais metodologias são as de dois fluidos, drift flux e
seguimento de pistões (slug tracking). Uma revisão dos modelos de dois fluidos e
drift flux é apresentada no apêndice C e, a seguir, é apresentada a revisão dos
modelos de seguimento de pistões.
Nos modelos de seguimento de pistões, cada pistão e bolha do escoamento
são considerados objetos distintos, que são propagados ao longo da tubulação.
Dessa forma, cada volume de controle computacional tem comprimento da ordem
dos comprimentos de bolhas e pistões, o que diminui consideravelmente o número
de equações a serem resolvidas e diminui erros na solução numérica.
Um dos primeiros trabalhos desenvolvidos com a metodologia de seguimento
de pistões foi o de Barnea e Taitel (1993). Os autores apresentam um modelo bem
simplificado, em que o líquido e o gás são considerados como incompressíveis, os
pistões não possuem bolhas de gás dispersas, e a fração de líquido no filme abaixo
da bolha é constante. Devido a essas considerações, todos os pistões possuem a
mesma velocidade, e a traseira da bolha tem a mesma velocidade de translação que
a frente da bolha. A velocidade de translação da frente da bolha é calculada, a cada
passo de tempo e para cada célula, em função da velocidade do pistão à sua frente
e uma correção devido ao comprimento do pistão à frente da bolha devido à esteira
(quanto menor o pistão, maior é a velocidade da bolha). Dessa forma, pistões
pequenos tendem a desaparecer, devido às coalescências.
Como condição de entrada das células na tubulação, os autores testaram a
entrada de uma população de pistões seguindo uma distribuição uniforme ou uma
distribuição normal, porém, obtiveram resultados semelhantes para as duas
condições. Os comprimentos das bolhas na entrada são calculados de forma
proporcional aos comprimentos dos pistões. Além disso, no modelo apresentado,
nenhuma consideração foi feita com relação à queda de pressão ou à variação de
quantidade de movimento nos pistões. Os autores testaram o modelo em uma
tubulação vertical e para o avanço das frentes das bolhas foi utilizado o método de
discretização por diferenças finitas explícito no tempo.
Capítulo 2 Revisão bibliográfica 16
Em seguida, Zheng et al. (1994) aprimoraram o modelo de Barnea e Taitel
(1993) para a simulação do escoamento em terrenos acidentados (hilly terrain), com
variação discreta na inclinação em determinados locais da tubulação. As inovações
desse modelo são a consideração dos pistões aerados, com a fração de líquido no
pistão calculada em função da velocidade do pistão (que continuava sendo igual
para todos os pistões) e a adição de modelos para a geração e dissipação de
pistões durante a passagem em uma mudança brusca de inclinação.
Logo depois, Nydal e Banerjee (1995) apresentaram um modelo considerando
o pistão como não-aerado, a fração de líquido no filme variável e o gás
compressível. Os autores utilizam quatro equações para calcular os parâmetros ao
longo do tempo, através de uma discretização explícita no tempo: equação de
conservação da quantidade de movimento no pistão de líquido para o cálculo da
velocidade do pistão, equação de conservação da massa na bolha para o cálculo da
pressão do interior da bolha, equação de conservação da massa de líquido no filme
para calcular a altura do filme, e equação de conservação da quantidade de
movimento no filme para calcular a velocidade do filme. As velocidades das frentes
das bolhas são calculadas em função da velocidade do pistão. Os autores aplicam o
modelo ao escoamento em terrenos acidentados, e introduzem uma população de
pistões na entrada do tubo distribuídos uniformemente.
Mais tarde, Taitel e Barnea (1998) propuseram um novo modelo de
seguimento de pistões, considerando diversos efeitos que haviam sido
negligenciados. Nesse modelo, o gás é considerado como compressível e ideal, a
fração de líquido no filme é variável e considerada igual à altura de filme de
equilíbrio, a fração de líquido no pistão é calculada em função da velocidade do
pistão (que agora é variável devido à expansão do gás) e a queda de pressão é
calculada. Os valores das velocidades dos pistões e das pressões nas bolhas são
calculados através de um sistema linear de equações que surgem devido à
aplicação do balanço de quantidade de movimento na célula, desprezando-se os
termos de variação e fluxo de quantidade de movimento (ou seja, a somatória de
forças na célula é nula), e do “balanço de volume” entre dois pistões adjacentes, que
representa a variação da velocidade superficial de um pistão até outro devido à
expansão do gás.
Capítulo 2 Revisão bibliográfica 17
Os autores apresentam resultados para esse modelo apenas no escoamento
horizontal, e com entrada de todos os pistões com mesmo comprimento. Pouco
tempo depois, Taitel e Barnea (2000) estenderam o seu modelo para o escoamento
em terrenos acidentados, incorporando os mecanismos de geração e
desaparecimento de pistões apresentado em Zheng et al. (1994).
Nydal et al. (2003) apresentaram um modelo semelhante ao apresentado por
Nydal e Banerjee (1995), porém, incorporaram a aeração do pistão e um modelo
para simular processos de pigging (limpeza de tubulações utilizando-se uma esfera
com diâmetro semelhante ao da tubulação). Nesse caso, a esfera é considerada
como um pistão, com alguns atributos modificados. Além disso, os autores
apresentam comparações com experimentos em tubulações com mudança de
inclinação, como risers.
Mais recentemente, Al-Safran et al. (2004) utilizaram o modelo proposto em
Taitel e Barnea (2000) e compararam os resultados das simulações com dados
experimentais obtidos no TUFFP (Tulsa University Fluid Flow Projects) para dois
casos de escoamento em terrenos acidentados. Segundo os autores, a comparação
dos resultados é boa.
Grenier (1997) apresentou um modelo de seguimento de pistões, que foi
reapresentado por Franklin e Rosa (2004), e segue alguns conceitos do modelo de
Nydal e Banerjee (1995). O gás é considerado compressível e ideal, o pistão é não-
aerado e a fração de líquido no filme é constante. O modelo é baseado na aplicação
das equações de conservação da massa e quantidade de movimento na forma
integral a volumes de controle que seguem as bolhas e pistões. A equação de
conservação da quantidade de movimento é aplicada ao pistão de líquido, sem
desconsiderar os termos de variação e fluxo da quantidade de movimento no pistão,
resultando em uma equação para a variação temporal da velocidade do pistão.
Franklin e Rosa (2004) mostraram que o termo de fluxo de quantidade de
movimento pode ser escrito como a soma dos termos de queda de pressão devido à
aceleração do líquido que passa do filme para o pistão e devido à diferença de
pressão hidrostática entre o filme na traseira da bolha e o pistão. A equação de
conservação da massa é aplicada à fase de gás na célula resultando em uma
equação para a variação temporal da pressão no interior da bolha. A entrada das
Capítulo 2 Revisão bibliográfica 18
bolhas e pistões na tubulação é realizada de maneira periódica (todas as bolhas e
pistões possuem os mesmos comprimentos).
Mais recentemente, Ujang et al. (2006) propuseram seu modelo de
seguimento de pistões. Os autores consideraram o gás como incompressível e o
pistão é aerado, porém as velocidades do líquido e do gás no pistão são iguais. A
fração de líquido no filme é calculada através do balanço de quantidade de
movimento em regime estacionário e a transferência de gás da bolha para o pistão
seguinte é modelada através de uma equação experimental. As variáveis primárias
são a massa de gás e o comprimento de cada objeto (bolha ou pistão) e a
discretização é explícita no tempo. Os autores propõem uma forma de inserção de
bolhas a partir de uma distribuição não-correlacionada de Poisson, baseada em
dados experimentais. O modelo só é válido para o escoamento horizontal, devido às
correlações utilizadas para modelar a transferência de gás da bolha para o pistão
seguinte.
Finalmente, Rodrigues (2006) utilizou o modelo proposto por Franklin e Rosa
(2004), porém, o pistão de líquido foi considerado como aerado e a equação de
conservação da quantidade de movimento foi aplicada a toda a célula,
considerando-se os termos de variação de quantidade de movimento e fluxo de
quantidade de movimento no filme de líquido. Além disso, o autor propôs um modelo
para a inserção das células de forma intermitente (aleatória), com uma população de
bolhas e pistões calculados com base em uma distribuição normal para a velocidade
das bolhas na entrada.
No entanto, os resultados mostraram que os termos de inércia do filme de
líquido não são relevantes, e o modelo de intermitência utilizado introduzia a
intermitência somente em uma variável, a velocidade da bolha. Além disso, alguns
resultados para o escoamento vertical foram apresentados, porém, sem muitos
dados experimentais para comparação.
O presente trabalho é uma evolução do trabalho de Rodrigues (2006). Na
modelagem apresentada, as variáveis principais são as velocidades do líquido nos
pistões e as pressões no interior das bolhas, que são resultado da resolução de um
sistema de equações resolvido a cada passo de tempo. As frações de líquido no
pistão e no filme são consideradas variáveis no tempo e a entrada de bolhas e
pistões é realizada a partir de populações independentes para os comprimentos de
Capítulo 2 Revisão bibliográfica 19
bolhas, comprimentos de pistões e velocidades superficiais de cada célula. Além
disso, são apresentados resultados e comparações para escoamentos horizontais e
verticais de ar e água, ar e solução de glicerina, e nitrogênio e óleo SAE20-50.
2.2 Definição de parâmetros importantes associados ao escoamento bifásico
em golfadas
Em escoamentos bifásicos, a diferenciação entre as velocidades de cada
fase, velocidade da mistura e velocidade superficial de cada fase, e a relação
dessas velocidades com as vazões volumétricas de cada fase são importantes. Em
uma seção transversal genérica de um escoamento bifásico, existe uma área
ocupada pelo líquido, LA , e uma área ocupada pelo gás, GA . As velocidades locais
das fases de líquido e gás são calculadas, respectivamente, por:
e GLL G
L G
QQU UA A
= = , (2.1)
onde LQ e GQ são as vazões volumétricas de líquido e gás na seção transversal
analisada. Além disso, define-se a fração de uma fase como a razão entre a área
ocupada por essa fase e a área de seção transversal, e, dessa forma:
e GLL G
AAR RA A
= = , (2.2)
onde A é a área da seção transversal. Substituindo-se as definições da Equação
(2.2), as velocidades de cada fase são escritas como:
e GLL G
L G
QQU UR A R A
= = . (2.3)
Capítulo 2 Revisão bibliográfica 20
Uma variável mais comum de se utilizar é a velocidade superficial da fase.
Essa velocidade é definida como a velocidade que a fase teria se estivesse
escoando sozinha na tubulação, e é escrita como:
e GLL G
QQJ JA A
= = , (2.4)
onde LJ e GJ são as velocidades superficiais de líquido e gás.
A vazão volumétrica total na seção transversal é dada por:
L GQ Q Q= + , (2.5)
e, dessa forma, define-se a velocidade média da mistura (ou velocidade de mistura)
na seção transversal como:
L GQ QQJA A
+= = , (2.6)
e utilizando-se a Equação (2.4) obtém-se:
L GJ J J= + , (2.7)
ou seja, a velocidade da mistura é a somatória das velocidades superficiais de cada
fase. Além dessas definições mais importantes, é definida ainda a velocidade de
deslizamento de cada fase, como a diferença entre a velocidade da fase e a
velocidade da mistura, por exemplo, para o gás:
GD GU U J= − , (2.8)
e a velocidade relativa, como a diferença entre as velocidades das duas fases, por
exemplo, para o gás:
Capítulo 2 Revisão bibliográfica 21
GR G LU U U= − . (2.9)
Combinando-se as Equações (2.8) e (2.9) obtém-se uma relação entre as
velocidades de deslizamento e relativa em função da fração de líquido:
GD L GRU R U= . (2.10)
As velocidades superficiais de líquido e gás, LJ e GJ , são variáveis de grande
utilização em escoamentos bifásicos, pois são facilmente relacionadas à vazão
volumétrica das fases. No escoamento de líquido e gás, geralmente, o líquido é
tratado como incompressível, e o gás como compressível. Dessa forma, a vazão
volumétrica do líquido (e, conseqüentemente, a velocidade superficial) é constante
ao longo de todo o escoamento. Por outro lado, a vazão volumétrica do gás varia ao
longo do escoamento, de acordo com a pressão. A Figura 2.5 apresenta um volume
de controle genérico englobando uma longa seção de uma tubulação com
escoamento bifásico de líquido e gás. A equação de conservação da massa de gás
em regime permanente aplicada a esse volume de controle resulta em:
1 1 2 2G G G GQ Qρ ρ= , (2.11)
onde Gρ é a massa específica do gás e GQ é a vazão volumétrica do gás, e os sub-
índices 1 e 2 referem-se às seções da tubulação mostradas na Figura 2.5.
Utilizando-se as definições de velocidade superficial, Equação (2.4), obtém-se:
1 1 2 2G G G GJ Jρ ρ= . (2.12)
Além disso, é comum considerar-se o gás como ideal, e, dessa forma, a equação de
estado para escoamento isotérmico resulta em:
cteG
G
Pρ
= , (2.13)
Capítulo 2 Revisão bibliográfica 22
onde GP é a pressão do gás. Dessa forma, a Equação (2.12) é reescrita como:
1 1 2 2G G G GP J P J= , (2.14)
ou:
12 1
2
GG G
G
PJ JP
= , (2.15)
e a velocidade superficial do gás em qualquer seção do escoamento é calculada em
função da velocidade superficial em outra seção corrigida pela razão entre as
pressões. A Equação (2.15) mostra que, quando ocorre variação de pressão, a
velocidade superficial do gás também varia. No entanto, para pequenas variações
de pressão, a variação na velocidade superficial do gás também será pequena, e,
em alguns casos, desconsiderada.
ρG1 G1Q ρG2 G2Q1 2
Figura 2.5 – Volume de controle ao longo de um trecho de tubulação
2.3 Equações de fechamento
Tanto nos modelos estacionários quanto nos transientes, são necessárias
equações de fechamento. Nessa seção é apresentada uma revisão das equações
de fechamento para a velocidade das bolhas de gás, para a fração de líquido no
pistão e para a freqüência da célula unitária.
Capítulo 2 Revisão bibliográfica 23
2.3.1 Velocidade de translação das bolhas de gás
As bolhas de gás no escoamento em golfadas são as bolhas alongadas e as
pequenas bolhas que podem aparecer dispersas no interior do pistão de líquido. A
seguir a forma como as velocidades dessas bolhas são calculadas é apresentada.
A velocidade das bolhas dispersas no interior do pistão de líquido GSU pode
ser escrita em função da velocidade da mistura no pistão, J , como:
GS GDU J U= + , (2.16)
onde GDU é a velocidade de deslizamento das bolhas dispersas, ou em função da
velocidade do líquido no pistão, LSU , como:
GS LS GRU U U= + , (2.17)
onde GRU é a velocidade relativa das bolhas dispersas no pistão. Uma relação para
a velocidade relativa das bolhas dispersas no escoamento vertical foi proposta por
Zuber e Hench apud Fernandes et al. (1983) e tem a forma:
( ) 0,250,5
21,53 L GGR LS
L
gU R
σ ρ ρρ
−⎡ ⎤= ⎢ ⎥
⎣ ⎦, (2.18)
onde σ é a tensão superficial entre o líquido e o gás, g é a aceleração da gravidade
e Lρ e Gρ são as massas específicas do líquido e do gás. Uma relação semelhante
pode ser obtida para a velocidade de deslizamento das bolhas dispersas utilizando-
se a Equação (2.10).
A velocidade de translação da bolha alongada tem um comportamento
diferente das bolhas dispersas, devido ao seu formato (é comum utilizar-se o termo
“bolha de Taylor” para se referenciar as bolhas alongadas do escoamento em
golfadas na vertical). Por isso, a velocidade de translação da bolha alongada não é
descrita conforme a Equação (2.17), mas é calculada a partir da superposição de
três efeitos, de acordo com a seguinte equação:
Capítulo 2 Revisão bibliográfica 24
( )( )0 1T DU C J U h= + + , (2.19)
onde TU é a velocidade de translação da bolha alongada, J , a velocidade da
mistura no pistão à sua frente, 0C , a constante de influência do líquido em
movimento, DU , a velocidade de deslizamento e h , a constante que caracteriza a
influência da esteira da bolha precedente. A definição da Equação (2.19) é baseada
em três efeitos que ocorrem, a princípio, no escoamento vertical. No entanto,
equações semelhantes à (2.19) têm sido utilizadas para o escoamento horizontal.
Quando uma bolha está inserida em um líquido estacionário, a bolha se
movimenta devido às forças de empuxo. Esse efeito é caracterizado pela velocidade
de deslizamento, dada por DU . Além disso, se o líquido em que a bolha está inserida
estiver em movimento, a velocidade de ascensão da bolha é aumentada. A
constante 0C quantifica a relação entre a velocidade da mistura a frente da bolha e a
velocidade da bolha. Finalmente, em uma seqüência de bolhas, a esteira da bolha
precedente pode influenciar a velocidade da bolha seguinte. Esse efeito é
quantificado pela constante de esteira, h .
Diversos estudos têm sido realizados para a definição das constantes para
cálculo da velocidade da bolha. Alguns modelos estão apresentados na Tabela B.1
do apêndice B. De modo geral, considera-se que a velocidade da bolha alongada é
aproximadamente igual à velocidade máxima (na seção transversal) do pistão à sua
frente. Dessa forma, em escoamentos turbulentos, a constante 0C é
aproximadamente 1,2, enquanto em escoamentos laminares é próxima de 2.
Além disso, no escoamento horizontal, nota-se um número de Froude
(número adimensional que representa a razão entre as forças de inércia e
gravitacionais) crítico no qual a velocidade da bolha tem uma grande variação. Isso
se deve a movimentação vertical da frente da bolha, que tende a desencostar do
tubo quando o número de Froude é muito alto.
2.3.2 Fração de líquido no pistão
A fração de líquido no pistão, LSR , também conhecida por hold-up de líquido,
é um parâmetro importante no escoamento em golfadas. Seu conhecimento é
Capítulo 2 Revisão bibliográfica 25
necessário para o cálculo correto da fração de líquido na região do filme de líquido,
uma vez que todo o gás presente na célula unitária está dividido entre a bolha
alongada e as pequenas bolhas dispersas no pistão. Além disso, a fração de líquido
no pistão é um fator primordial para o cálculo da queda de pressão devido à força da
gravidade, no escoamento vertical ou inclinado.
Após algumas observações experimentais, nota-se que no escoamento
vertical ascendente o pistão de líquido contém muitas bolhas dispersas, ao contrário
do escoamento horizontal. Esse comportamento é esperado, pois no escoamento
vertical o filme de líquido ao redor da bolha alongada tem velocidade negativa, ou
seja, está escoando em direção oposta ao pistão de líquido que o segue. No
momento em que o líquido passa do filme ao pistão de líquido seguinte, uma zona
de grande agitação é formada. Essa agitação é responsável pela remoção de
pequenas porções de gás da traseira da bolha, formando bolhas dispersas no
pistão. O mesmo fenômeno não acontece no escoamento horizontal, pois o líquido
no filme escoa na mesma direção da bolha.
Um problema ainda em aberto é a determinação da fração de líquido no
pistão em diversas inclinações de escoamento e de seu comportamento em
singularidades do tubo, como mudanças de inclinação e curvas de raio curto.
As pequenas bolhas dispersas no pistão podem, ainda, apresentar uma
distribuição axial ou radial. No escoamento horizontal, devido à atuação da
gravidade, as pequenas bolhas tendem a escoar pela parte superior do tubo,
enquanto no escoamento vertical a maior concentração de bolhas está na parte
central. Além disso, devido à grande agitação formada na região próxima à bolha
que precede o pistão (principalmente no escoamento vertical), a concentração de
bolhas nessa região é maior. Isso acontece, pois uma parte das pequenas bolhas
que se desprende da bolha alongada acaba voltando para a mesma. A Tabela B.2
do apêndice B apresenta alguns modelos para se calcular a fração de líquido no
pistão.
Capítulo 2 Revisão bibliográfica 26
2.3.3 Freqüência da célula unitária
A freqüência de ocorrência das estruturas no escoamento bifásico líquido-gás
intermitente em golfadas é definida como o inverso do tempo de passagem de uma
célula, como:
1
U
ft
=Δ
, (2.20)
onde UtΔ é o tempo de passagem da célula unitária. Esse parâmetro é importante
para a modelagem do escoamento, pois se relaciona com a velocidade de
translação da bolha e comprimentos de pistão e líquido. De acordo com equações
dos modelos estacionários, a freqüência da célula unitária pode ser calculada como:
( )1, , TT T
U B S
UU Uf f fL L L
ββ −= = = , (2.21)
onde β é o fator de intermitência, definido como:
B
S B
LL L
β =+
. (2.22)
Portanto, conhecendo-se a freqüência do escoamento e algum outro
parâmetro pode-se, a partir das equações acima, obter uma estimativa dos
comprimentos médios característicos.
Além disso, a escolha do material das tubulações onde a mistura bifásica
escoa depende da freqüência do escoamento, pois é necessário evitar materiais
cuja freqüência natural de vibração seja próxima à freqüência de passagem das
células.
Existem diversos modelos para cálculo da freqüência em função de
parâmetros conhecidos do escoamento. Os mais antigos baseiam-se em resultados
de experimentos. A partir do ajuste de dados obtidos experimentalmente em função
de parâmetros conhecidos, foram propostas correlações para o cálculo da
Capítulo 2 Revisão bibliográfica 27
freqüência. As correlações podem ser simples, com o cálculo da freqüência em
função de LJ e GJ , ou mais elaboradas, considerando as configurações geométricas
e reológicas do escoamento.
Uma outra forma para o cálculo da freqüência em escoamentos horizontais
surge a partir da modelagem da formação de pistões a partir do escoamento
estratificado. Essa metodologia é mais consistente fisicamente, porém, apresenta
maior dificuldade de implementação. A Tabela B.3 do apêndice B apresenta alguns
modelos para o cálculo da freqüência da célula. É importante notar que os últimos 4
modelos apresentados na Tabela B.3 são modelos baseados na transição do
escoamento estratificado para o escoamento em golfadas.
Capítulo 3 Modelagem matemática 28
3 MODELAGEM MATEMÁTICA
Esse capítulo apresenta a modelagem matemática do escoamento em
golfadas visando à implementação em um modelo de seguimento de pistões. As
equações de conservação da massa e quantidade de movimento na forma integral e
unidimensional são aplicadas a volumes de controle móveis e deformáveis que
envolvem o pistão e a região da bolha de gás e filme de líquido. O desenvolvimento
dessas equações resulta em duas equações diferenciais em relação ao tempo
acopladas, para as variáveis LSU , velocidade média do líquido no pistão, e GBP , a
pressão na bolha, para a j-ésima célula (pistão seguido de uma bolha) genérica do
escoamento. Equações auxiliares necessárias para modelar o deslocamento da
frente e traseira da bolha ao longo do tubo em função de LSU e GBP , modelar as
forças de atrito na parede e calcular a velocidade do filme de líquido são
apresentadas. A modelagem das singularidades que ocorrem no escoamento em
golfadas também é apresentada.
As hipóteses utilizadas na modelagem são:
• Modelo unidimensional.
• Líquido incompressível.
• Gás ideal.
• Escoamento isotérmico e sem mudança de fase.
• Quantidade de movimento associada ao gás e à região do filme de líquido
desprezível.
• Pressão do gás não varia axialmente no interior de uma bolha alongada.
• Concentração de bolhas de gás em cada pistão não varia axialmente.
• Filme de líquido não aerado.
• Forças interfaciais entre o gás e o líquido desprezadas.
3.1 Equações de conservação aplicadas a uma célula genérica
Nessa seção as equações de conservação da massa e de quantidade de
movimento são aplicadas a volumes de controle que envolvem a j-ésima célula do
Capítulo 3 Modelagem matemática 29
escoamento, representada na Figura 3.1. O comprimento do pistão de líquido j é
SjL , o comprimento da bolha é BjL , a pressão do gás no interior da bolha j é GBjP ,
LSjU e GSjU são as velocidades médias do líquido e do gás, respectivamente, no
pistão j , LBjU é a velocidade média do líquido no filme e θ é o ângulo de inclinação
do tubo com a horizontal.
A equação de conservação da massa é aplicada às fases de líquido e de gás
em toda a célula, enquanto a equação de conservação da quantidade de movimento
é aplicada apenas à fase de líquido no interior do pistão. Para facilitar a modelagem,
primeiramente a equação de conservação da massa é aplicada às fases de líquido e
de gás no pistão, e, na seqüência, às fases de líquido e de gás na região da bolha
alongada e filme de líquido. Em seguida, a equação de conservação da quantidade
de movimento é aplicada ao pistão. Ao final da modelagem, duas equações
diferenciais em relação ao tempo são obtidas para as variáveis LSjU e GBjP .
j ésima célula
LSj
LBj
ULSj
ULBj
PGBj
UGSj
Figura 3.1 – Representação da j-ésima célula genérica do escoamento em golfadas
A equação de conservação da massa na forma integral aplicada a um volume
de controle é dada por:
0VC SC
dV V dAt
ρ ρ∂+ ⋅ =
∂ ∫ ∫ , (3.1)
Capítulo 3 Modelagem matemática 30
onde ρ é a massa específica, V é a velocidade, A é a área superficial do volume
de controle, V é o volume e t é o tempo. O primeiro termo da Equação (3.1) é a
variação temporal da massa no interior do volume de controle, e o segundo termo é
o fluxo líquido de massa que cruza a superfície de controle.
A equação de conservação da quantidade de movimento unidimensional na
forma integral em relação a um referencial inercial é:
VC SC
udV uV dA Ft
ρ ρ∂+ ⋅ =
∂ ∑∫ ∫ , (3.2)
onde u é a componente axial da velocidade e F∑ representa as forças externas
que atuam no volume de controle. A seguir as Equações (3.1) e (3.2) são aplicadas
ao escoamento em golfadas.
3.1.1 Conservação da massa na célula
A seguir, a equação de conservação da massa é aplicada para as fases de
líquido e de gás presente nas regiões do pistão e da bolha e filme. O objetivo é obter
uma equação que relacione a pressão do gás na bolha, GBP , com a velocidade do
líquido no pistão, LSU .
Líquido no pistão
A Figura 3.2 apresenta o volume de controle (tracejado) definido pelas
fronteiras do pistão j , ou seja, a parede do tubo na direção radial, e as coordenadas
jx e jy na direção axial. São apresentados, também, os fluxos de massa de líquido
e gás que cruzam as fronteiras do volume de controle nas coordenadas jx e jy .
Capítulo 3 Modelagem matemática 31
yj LSj
xj
mLxj
mGxj
mGyj
mLyj
.
.
.
.
Figura 3.2 – Volume de controle definido pelas fronteiras do pistão
A equação de conservação da massa para a fase de líquido no interior do
volume de controle apresentado na Figura 3.2 é:
( ) 0j jL LSj Sj L LSxj LSxj L LSyj LSyj
dx dyd AR L AR U AR Udt dt dt
ρ ρ ρ⎛ ⎞ ⎛ ⎞
+ − − − =⎜ ⎟ ⎜ ⎟⎝ ⎠ ⎝ ⎠
, (3.3)
onde A é a área da seção transversal do tubo, Sj j jL x y= − é o comprimento do
pistão, LSjR é a fração volumétrica de líquido média no pistão, definida pela razão
entre o volume de líquido e o volume total do pistão, LSxjR e LSyjR são as frações de
líquido nas fronteiras jx e jy , respectivamente, definidas como a razão entre a área
de líquido e a área da seção transversal do tubo, Lρ é a massa específica do
líquido, LSxjU e LSyjU são as velocidades do líquido nas fronteiras jx e jy ,
respectivamente, jdx dt e jdy dt são as velocidades das fronteiras jx e jy ,
respectivamente, e t é o tempo.
O primeiro termo da Equação (3.3) representa a variação da massa de líquido
no interior do pistão, enquanto os segundo e terceiro termos representam a
contribuição dos fluxos de massa de líquido através da superfície de controle.
As frações de líquido nas fronteiras do pistão, LSxjR e LSyjR , são ligeiramente
diferentes entre si. No entanto, a consideração dessa diferença tornaria a
modelagem complexa. Por isso, nesse trabalho, é considerado que as frações de
líquido nas fronteiras são iguais entre si e iguais à fração de líquido média ao longo
Capítulo 3 Modelagem matemática 32
do pistão. Essa aproximação significa que as bolhas dispersas estão uniformemente
distribuídas axialmente ao longo do pistão. Dessa forma:
LSxj LSyj LSjR R R= = . (3.4)
Substituindo-se essa definição e dividindo-se por L Aρ , a Equação (3.3) se torna:
SjLSj
dLR
dt( )LSj j j
Sj LSj LSxj LSyj LSj
dR dx dyL R U U R
dt dt dt⎛ ⎞
+ + − − −⎜ ⎟⎝ ⎠
0= . (3.5)
Desenvolvendo-se a Equação (3.5) e reorganizando obtém-se:
Sj LSjLSxj LSyj
LSj
L dRU U
R dt− = − . (3.6)
A Equação (3.6) indica que as velocidades do líquido nas fronteiras do pistão não
são iguais entre si, devido à variação temporal da fração de líquido média no pistão.
Essa equação possui duas incógnitas, LSxjU e LSyjU , e, para ser resolvida, uma
equação adicional deve ser utilizada. Aproximando-se a velocidade média do líquido
no pistão, LSjU , como a média aritmética das velocidades do líquido nas fronteiras:
2
LSxj LSyjLSj
U UU
+= , (3.7)
e resolvendo-se em conjunto com a Equação (3.6), as velocidades nas fronteiras do
pistão são expressas em função da velocidade média do pistão e da taxa de
variação de LSjR :
2
Sj LSjLSxj LSj
LSj
L dRU U
R dt= − , (3.8)
Capítulo 3 Modelagem matemática 33
e:
2
Sj LSjLSyj LSj
LSj
L dRU U
R dt= + . (3.9)
Nas Equações (3.8) e (3.9) está implícita a conservação da massa de líquido no
pistão. A seguir, uma análise semelhante será realizada para a fase de gás.
Gás no pistão
A equação de conservação da massa para a fase de gás no volume de controle da
Figura 3.2 é:
( ) ( )
( )
11 1
1 0
jGSj LSj Sj GBj LSxj GSxj
jGBj LSyj GSyj
dxd A R L A R Udt dt
dyA R U
dt
ρ ρ
ρ
+
⎛ ⎞⎡ ⎤− + − − −⎜ ⎟⎣ ⎦ ⎝ ⎠⎛ ⎞
− − − =⎜ ⎟⎝ ⎠
, (3.10)
onde GSxjU e GSyjU são as velocidades do gás no pistão junto às fronteiras jx e jy ,
respectivamente, GSjρ é a massa específica média do gás ao longo do pistão, GBjρ e
1GBjρ + são as massas específicas do gás no interior das bolhas j e 1j + ,
respectivamente.
O primeiro termo da Equação (3.10) representa a variação da massa de gás
no interior do pistão, enquanto os segundo e terceiro termos representam os fluxos
mássicos de gás que cruzam a superfície de controle. É importante ressaltar que a
massa específica do gás, ao contrário da massa específica do líquido, é variável em
relação ao tempo e ao espaço.
A Equação (3.10) é de difícil resolução, pois a massa específica do gás que
cruza a superfície de controle é diferente da massa específica média ao longo do
pistão. No entanto, essa diferença é muito pequena, devido à ordem de grandeza do
comprimento do pistão em relação à tubulação. Assim, utiliza-se a aproximação:
Capítulo 3 Modelagem matemática 34
1GSj GBj GBjρ ρ ρ+∼ ∼ . (3.11)
Dividindo-se a Equação (3.10) pela área do tubo e pela massa específica do gás na
bolha j , GBjρ , desenvolvendo-se os termos da derivada e considerando-se as
simplificações das Equações (3.4) e (3.11), obtém-se:
( ) ( )
( ) ( )
11 1
1 1 0
Sj LSj GBjLSj Sj Sj LSj
GBj
j jLSj GSxj LSj GSyj
dL dR dR L L R
dt dt dt
dx dyR U R U
dt dt
ρρ
− − + − +
⎛ ⎞ ⎛ ⎞+ − − − − − =⎜ ⎟ ⎜ ⎟
⎝ ⎠ ⎝ ⎠
. (3.12)
Dividindo-se a Equação (3.12) por ( )1 LSjR− e reorganizando os termos, obtém-se:
( ) 01
Sj Sj LSj Sj GBj SjGSxj GSyj
GBjLSj
dL L dR L d dLU U
dt dt dt dtRρ
ρ− + + − − =
−, (3.13)
e evidenciando-se as velocidades do gás no pistão:
( )1 1
1LSj GBj
GSxj GSyj SjGBjLSj
dR dU U L
dt dtRρ
ρ
⎡ ⎤− = −⎢ ⎥
−⎢ ⎥⎣ ⎦. (3.14)
A Equação (3.14) é similar à Equação (3.6) obtida no balanço de massa de líquido
no pistão. Da mesma forma, é necessária uma segunda equação para se calcular as
velocidades do gás nas fronteiras do pistão em função da velocidade de gás média.
Definindo-se a velocidade de gás média no pistão, GSjU , como a média aritmética
das velocidades do gás nas fronteiras:
2
GSxj GSyjGSj
U UU
+= , (3.15)
e combinando-se com a Equação (3.14) obtém-se:
Capítulo 3 Modelagem matemática 35
( )
1 12 1Sj LSj GSj
GSxj GSjGSjLSj
L dR dU U
dt dtRρ
ρ
⎡ ⎤= + −⎢ ⎥
−⎢ ⎥⎣ ⎦, (3.16)
e:
( )
1 12 1Sj LSj GSj
GSyj GSjGSjLSj
L dR dU U
dt dtRρ
ρ
⎡ ⎤= − −⎢ ⎥
−⎢ ⎥⎣ ⎦. (3.17)
As Equações (3.16) e (3.17) garantem a conservação da massa de gás no pistão de
líquido. A seguir serão obtidas as equações para a região da bolha.
Líquido no filme
Na seção anterior as equações de conservação da massa de gás e líquido no
pistão foram obtidas. Os resultados foram expressos através das velocidades que
essas fases têm nas fronteiras que o pistão faz com as bolhas adjacentes. Essas
velocidades serão utilizadas a seguir, nos balanços de massa de líquido no filme e
de gás na bolha alongada.
A Figura 3.3 apresenta um volume de controle (tracejado) definido pelas
fronteiras da bolha alongada e filme de líquido, ou seja, a parede do tubo na direção
radial e as coordenadas 1jx − e jy na direção longitudinal. São apresentados,
também, os fluxos mássicos de líquido e gás que cruzam as fronteiras do volume de
controle.
yj
LBjxj-1
mLxj-1
mGxj-1
mGyj
mLyj
.
.
.
.
Figura 3.3 – Volume de controle definido pelas fronteiras do filme e bolha alongada
Capítulo 3 Modelagem matemática 36
A equação de conservação da massa para a fase líquida do volume de
controle apresentado na Figura 3.3 é:
( ) 11 1 0j j
L LBj Bj L LSyj LSyj L LSxj LSxj
dy dxd AR L AR U AR Udt dt dt
ρ ρ ρ −− −
⎛ ⎞ ⎛ ⎞+ − − − =⎜ ⎟ ⎜ ⎟
⎝ ⎠ ⎝ ⎠, (3.18)
onde LBjR é a fração volumétrica de líquido média na região da bolha e filme,
definida pela razão entre o volume de líquido no filme e o volume total da seção de
comprimento 1Bj j jL y x −= − . Desenvolvendo-se os termos da derivada da Equação
(3.18), dividindo-a por L Aρ e substituindo-se a definição da Equação (3.4) e o
resultado das Equações (3.8) e (3.9) para as velocidades do líquido nas fronteiras do
pistão (dessa forma a equação de conservação da massa de líquido no pistão é
satisfeita), a Equação (3.18) se torna:
1
1 1 11 1
1
2
02
j j LBj Sj LSj jLBj Bj LSj LSj
LSj
Sj LSj jLSj LSj
LSj
dy dx dR L dR dyR L R U
dt dt dt R dt dt
L dR dxR U
R dt dt
−
− − −− −
−
⎛ ⎞⎛ ⎞− + + + − −⎜ ⎟⎜ ⎟ ⎜ ⎟⎝ ⎠ ⎝ ⎠
⎛ ⎞− − − =⎜ ⎟⎜ ⎟
⎝ ⎠
. (3.19)
Reorganizando-se a Equação (3.19) e evidenciando-se os termos comuns obtém-se:
( ) ( ) 1
1
1 11 12 2
j j LBjLSj LBj LSj LBj Bj
Sj LSj Sj LSjLSj LSj LSj LSj
dy dx dRR R R R L
dt dt dtL dR L dR
R U R Udt dt
−−
− −− −
− − − − −
− − = −. (3.20)
A Equação (3.20) é a conservação de massa de líquido no filme e no pistão. Nota-se
que o lado esquerdo dessa equação apresenta somente parâmetros, e o lado direito
apresenta parâmetros cinemáticos.
Capítulo 3 Modelagem matemática 37
Gás na bolha
A equação de conservação da massa de gás para o volume de controle que
envolve a bolha, representado na Figura 3.3 é:
( ) ( )
( ) 11 1
1 1
1 0
jGBj LBj Bj GBj LSyj GSyj
jGBj LSxj GSxj
dyd A R L A R Udt dt
dxA R U
dt
ρ ρ
ρ −− −
⎛ ⎞⎡ ⎤− + − − −⎜ ⎟⎣ ⎦ ⎝ ⎠⎛ ⎞
− − − =⎜ ⎟⎝ ⎠
. (3.21)
Utilizando-se a aproximação da Equação (3.4) para as frações de líquido nas
fronteiras, e dividindo-se a Equação (3.21) por GBj Aρ :
( ) ( )
( ) ( ) 11 1
11 1
1 1 0
Bj LBj GBjLBj Bj Bj LBj
GBj
j jLSj GSyj LSj GSxj
dL dR dR L L R
dt dt dt
dy dxR U R U
dt dt
ρρ
−− −
− − + − +
⎛ ⎞ ⎛ ⎞+ − − − − − =⎜ ⎟ ⎜ ⎟
⎝ ⎠ ⎝ ⎠
. (3.22)
Desenvolvendo-se a Equação (3.22) e evidenciando-se os termos comuns, obtém-
se:
( ) ( )
( ) ( ) ( )
11
1 1
1 1 1 1
11 1 1 0
j j LBjLBj LSj LBj LSj Bj
GBjBj LBj LSj GSjy LSj GSxj
GBj
dy dx dRR R R R L
dt dt dtd
L R R U R Udtρ
ρ
−−
− −
− − + − − − + − +
+ − + − − − =. (3.23)
Substituindo-se, ainda, os resultados das Equações (3.16) e (3.17) para as
velocidades do gás nas fronteiras e reorganizando, obtém-se:
Capítulo 3 Modelagem matemática 38
( ) ( )
( ) ( )
( )
1 1 11
1 11 1
1
2 21 11 1
2 2
11
j j LBj Sj LSj Sj LSjLSj LBj LSj LBj Bj
Sj GBj Sj GBjLSj GSj LSj GSj
GBj GBj
GBjBj LBj
GBj
dy dx dR L dR L dRR R R R L
dt dt dt dt dtL d L d
R U R Udt dt
dL R
dt
ρ ρρ ρ
ρρ
− − −−
− −− −
−
− − − − − − =
⎛ ⎞ ⎛ ⎞− − + + − − −⎜ ⎟ ⎜ ⎟⎜ ⎟ ⎜ ⎟
⎝ ⎠ ⎝ ⎠
− −
. (3.24)
A Equação (3.24) representa a conservação da massa de gás na região da bolha.
Nota-se que o lado esquerdo da Equação (3.24) possui informações geométricas da
célula que são coincidentes com os da Equação (3.20), e o lado direito apresenta
relações cinemáticas para o gás na célula. A seguir as Equações (3.20) e (3.24)
serão combinadas.
Equações combinadas
Os lados esquerdos das Equações (3.20) e (3.24) são iguais. Dessa forma,
essas equações podem ser combinadas, o que resulta em:
( ) ( )
( ) ( ) ( )
1 1
11 11 1
1
1 12
11 1 0
2
LBj LSjGBj Sj GBjLSj LSj LSj LSj Bj
GBj GBj
LSjSj GBjLSj GSj LSj GSj
GBj
R Rd L dR U R U L
dt dt
RL dR U R U
dt
ρ ρρ ρ
ρρ
− −
−− −− −
−
− −− + + +
−+ − − − =
. (3.25)
A Equação (3.25) apresenta as derivadas em relação ao tempo das massas
específicas do gás no interior das bolhas j e 1j − . Porém, como simplificação,
considera-se a diferença entre esses termos pequena. Assim, utiliza-se a
simplificação:
1GBj GBjd ddt dt
ρ ρ− ∼ . (3.26)
Capítulo 3 Modelagem matemática 39
Além disso, a velocidade média do gás no pistão é escrita em função da velocidade
média do líquido no pistão, de acordo com a revisão bibliográfica apresentada no
Capítulo 2. Dessa forma:
DSjGSj LSj
LSj
UU U
R= + , (3.27)
onde DSjU é a velocidade de deslizamento das bolhas dispersas no pistão de líquido.
Substituindo-se essas definições, a Equação (3.25) se torna:
( ) ( ) ( )
( ) ( )
111 1
1
11 1
1
1 1 12 2
1 1 0
LBj LSj LSjGBj Sj SjLSj LSj LSj LSj Bj
GBj GBj GBj
DSj DSjLSj LSj LSj LSj
LSj LSj
R R Rd L LR U R U L
dt
U UR U R U
R R
ρρ ρ ρ
−−− −
−
−− −
−
⎡ ⎤− − −− + + + +⎢ ⎥
⎢ ⎥⎣ ⎦⎛ ⎞ ⎛ ⎞
+ − + − − + =⎜ ⎟ ⎜ ⎟⎜ ⎟ ⎜ ⎟⎝ ⎠ ⎝ ⎠
. (3.28)
O gás presente na bolha é considerado como ideal à temperatura constante. Dessa
forma, sua massa específica e pressão são relacionadas por:
1 1G G
G G
d dPdt P dtρ
ρ= . (3.29)
Substituindo a definição (3.29), a Equação (3.28) é reorganizada como:
( ) ( ) ( )111
1
11
1
1 1 12 2
1 1
LBj LSj LSjGBj Sj SjLSj LSj Bj
GBj GBj GBj
LSj LSjDSj DSj
LSj LSj
R R RdP L LU U L
dt P P P
R RU U
R R
−−−
−
−−
−
⎡ ⎤− − −− = + + +⎢ ⎥
⎢ ⎥⎣ ⎦⎛ ⎞ ⎛ ⎞− −
+ −⎜ ⎟ ⎜ ⎟⎜ ⎟ ⎜ ⎟⎝ ⎠ ⎝ ⎠
. (3.30)
A Equação (3.30) representa o balanço da conservação da massa na célula unitária.
Observa-se que a diferença das velocidades do líquido nos pistões adjacentes é
devida a expansão do gás da bolha entre os pistões (primeiro termo do lado direito)
e a diferença de velocidade das bolhas dispersas nos pistões (segundo e terceiro
Capítulo 3 Modelagem matemática 40
termos do lado direito). Essa equação possui como incógnitas as velocidades do
líquido nos pistões j e 1j − , LSjU e 1LSjU − , e a derivada em relação ao tempo da
pressão do gás no interior da bolha j , GBjP .
3.1.2 Conservação da quantidade de movimento no pistão
A equação de conservação da quantidade de movimento para a fase líquida
do pistão é apresentada a seguir. É importante ressaltar que a quantidade de
movimento relativa à fase de gás foi desprezada, devido à grande diferença entre as
massas específicas do gás e do líquido. A equação de conservação da quantidade
de movimento aplicada à fase líquida do volume de controle definido pelo pistão
(Figura 3.2) é:
( )Sj L LSj Sj LSj LSxj Lxj LSyj LyjdF AR L U U m U mdt
ρ= + −∑ , (3.31)
onde os fluxos de massa de líquido que cruzam a superfície de controle nas
coordenadas jx e jy , respectivamente, são dados por:
jLxj L LSxj LSxj
dxm AR U
dtρ
⎛ ⎞= −⎜ ⎟
⎝ ⎠, (3.32)
e:
jLyj L LSyj LSyj
dym AR U
dtρ
⎛ ⎞= −⎜ ⎟
⎝ ⎠, (3.33)
e o termo do lado esquerdo da Equação (3.31) representa a somatória de forças
atuantes no pistão, que serão discutidas mais à frente. Substituindo-se as Equações
(3.32) e (3.33) na Equação (3.31) encontra-se:
Capítulo 3 Modelagem matemática 41
( ) j
Sj L LSj Sj LSj LSxj L LSxj LSxj
jLSyj L LSyj LSyj
dxdF AR L U U AR Udt dt
dyU AR U
dt
ρ ρ
ρ
⎛ ⎞= + − −⎜ ⎟
⎝ ⎠⎛ ⎞
− −⎜ ⎟⎝ ⎠
∑. (3.34)
O primeiro termo do lado direito da Equação (3.34) é o acúmulo de quantidade de
movimento no interior do volume, enquanto o segundo e terceiro termos são os
fluxos de quantidade de movimento que cruzam as fronteiras. Expandindo-se a
derivada no primeiro termo do lado direito da Equação (3.34), utilizando-se a
simplificação da Equação (3.4) para as frações de líquido, e desenvolvendo-se,
obtém-se:
( )2 2
Sj LSj Sj LSjLSj Sj LSj LSj Sj LSj
L
j jLSj LSxj LSyj LSj LSxj LSyj
F dU dL dRR L R U L U
A dt dt dtdx dy
R U U R U Udt dt
ρ= + + +
⎛ ⎞+ − − −⎜ ⎟
⎝ ⎠
∑. (3.35)
As velocidades do líquido no pistão que aparecem nos termos de fluxo de
quantidade de movimento na Equação (3.35) são as velocidades locais nas
fronteiras jx e jy , enquanto a velocidade do líquido no termo de acúmulo de
quantidade de movimento é a velocidade média. Na seção anterior, a diferenças
entre esses valores foram apresentadas, e as Equações (3.8) e (3.9) foram
desenvolvidas. Substituindo-se essas equações na Equação (3.35), o quarto termo
do lado direito é escrito como:
( )2 2 2 LSjLSj LSxj LSyj Sj LSj
dRR U U L U
dt− = − , (3.36)
e o quinto termo do lado direito é:
2
j j Sj Sj LSj j jLSj LSxj LSyj LSj LSj
LSj
dx dy dL L dR dx dyR U U R U
dt dt dt R dt dt dt⎡ ⎤⎛ ⎞ ⎛ ⎞
− = − +⎢ ⎥⎜ ⎟ ⎜ ⎟⎢ ⎥⎝ ⎠ ⎝ ⎠⎣ ⎦
. (3.37)
Capítulo 3 Modelagem matemática 42
Substituindo-se as Equações (3.36) e (3.37) na Equação (3.35) obtém-se:
Sj LSj SjLSj Sj LSj LSj
L
F dU dLR L R U
A dt dtρ= +∑ 2LSj LSj
Sj LSj Sj LSj
SjLSj LSj
dR dRL U L U
dt dt
dLR U
dt
+ − −
−2
Sj LSj j j
LSj
L dR dx dyR dt dt dt
⎡ ⎤⎛ ⎞− +⎢ ⎥⎜ ⎟
⎢ ⎝ ⎠⎥⎣ ⎦
. (3.38)
Cancelando-se os termos iguais e agrupando-se os termos semelhantes, obtém-se:
2
Sj LSj Sj j j LSjLSj Sj Sj LSj
L
F dU L dx dy dRR L L U
A dt dt dt dtρ⎡ ⎤⎛ ⎞
= + + −⎢ ⎥⎜ ⎟⎝ ⎠⎣ ⎦
∑ . (3.39)
O termo do lado esquerdo da Equação (3.39) é o somatório das forças externas que
atuam no pistão, apresentadas na Figura 3.4 e definidas por:
( ) senSj LS LS LSj Sj L LSj Sjyj xjF P P A DL AR L gτ π ρ θ= − − −∑ . (3.40)
O primeiro termo do lado direito da Equação (3.40) corresponde à diferença das
pressões estáticas exercidas nas faces jy e jx do pistão. O segundo termo do lado
direito é a força de atrito exercida pela parede do tubo no líquido, e o terceiro termo
do lado direito é a força devido ao peso do líquido no volume de controle. Devido
aos resultados da equação de conservação da massa, é interessante que as
pressões que aparecem no primeiro termo do lado direito da Equação (3.40) sejam
relacionadas com as pressões do gás no interior das bolhas adjacentes ao pistão,
GBjP e 1GBjP + . O acoplamento dessas pressões é realizado na seção seguinte.
Capítulo 3 Modelagem matemática 43
PLS j|y
PLS j|x
τ πLSj
SjDL
A
xj
y j
.
.
A
ρθ
LLSj SjAR L gsen
Figura 3.4 – Forças que atuam no volume de controle definido pelo pistão
3.1.3 Acoplamentos entre o pistão e as bolhas adjacentes
Nessa seção, são utilizadas as equações de conservação da massa e
quantidade de movimento aplicadas a volumes de controle infinitesimais, localizados
nas fronteiras entre o pistão j e as bolhas j e 1j + . Dessa forma, é possível
escrever as pressões das fronteiras do pistão de líquido em função das pressões do
gás no interior das bolhas.
Na Figura 3.5 está representado o volume de controle (tracejado) que envolve
o acoplamento entre a frente da bolha j e a traseira do pistão j . Esse volume é
indeformável e suas fronteiras se movem com velocidade jdy dt . Aplicando-se a
equação de conservação da massa a esse volume de controle e anulando-se os
termos de acúmulo obtém-se:
j jL LSyj LSyj L LByj LByj Lyj
dy dyAR U AR U m
dt dtρ ρ
⎛ ⎞ ⎛ ⎞− = − =⎜ ⎟ ⎜ ⎟
⎝ ⎠ ⎝ ⎠, (3.41)
onde LByjR e LByjU são a fração de líquido e a velocidade do filme no início da bolha
alongada.
Capítulo 3 Modelagem matemática 44
myj.
LSj
x j
yj
PGBj
PLS|yj
Figura 3.5 – Representação do volume de controle no acoplamento da traseira do
pistão com a frente da bolha
Aplicando-se a equação de conservação da quantidade de movimento ao
mesmo volume de controle obtém-se:
( ) ( )GBj LS Lyj LSyj Lyj LByj Lyj LSyj LByjyjP P A m U m U m U U− = − = − . (3.42)
O lado esquerdo da Equação (3.42) apresenta a relação entre as pressões no pistão
e na bolha. Porém, de acordo com os resultados experimentais para a variação de
pressão apresentados no Capítulo 2 (Revisão bibliográfica), a pressão na traseira do
pistão é a mesma pressão no interior da bolha, ou seja, não há uma variação brusca
de pressão. Dessa forma, o lado esquerdo da Equação (3.42) deve ser nulo. Isso só
é verdadeiro se o lado direito também for nulo, e LSyj LByjU U= . No Capítulo 2 também
foram mostradas linhas de corrente do líquido ao redor da bolha. Pode-se dizer que,
devido à curvatura do nariz da bolha, o escoamento é como o início de uma
contração suave e sem recirculações. Assim, a variação da velocidade é pequena, e
considera-se que não há salto de velocidade e pressão na traseira do pistão,
portanto:
LS GBjyjP P= . (3.43)
Na Figura 3.6 está representado o volume de controle que envolve a região
do acoplamento entre o pistão j e a bolha 1j + . O volume de controle é
Capítulo 3 Modelagem matemática 45
indeformável e suas fronteiras se movem com velocidade igual a da traseira da
bolha, jdx dt .
mLxj.
LSj
x j
yj PGBj+1
PLS|xj
Figura 3.6 – Representação do volume de controle no acoplamento entre a frente do
pistão e a traseira da bolha
Aplicando-se a equação de conservação da massa ao volume representado
na Figura 3.6 e anulando-se os termos de acúmulo obtém-se:
1 1j j
L LSxj LSxj L LBj LBj Lxj
dx dxAR U AR U m
dt dtρ ρ + +
⎛ ⎞ ⎛ ⎞− = − =⎜ ⎟ ⎜ ⎟
⎝ ⎠ ⎝ ⎠, (3.44)
onde 1LBjR + e 1LBjU + são a fração de líquido e a velocidade do líquido no filme da
bolha 1j + .
A equação de conservação da quantidade de movimento aplicada ao mesmo
volume de controle é:
( ) ( ) ( )1 1 11LS GBj HS HB Lxj LBj Lxj LSxj Lxj LBj LSxjxj xj jP P A F F m U m U m U U+ + ++
− + − = − = − , (3.45)
onde o segundo termo do lado esquerdo representa a diferença entre as pressões
hidrostáticas no pistão e na bolha. Esse termo aparece apenas nessa equação, pois
a variação da seção transversal preenchida com líquido varia bruscamente na
traseira da bolha, o que não acontece no acoplamento na frente da bolha.
Capítulo 3 Modelagem matemática 46
A Equação (3.45) apresenta a relação entre a pressão na frente do pistão j e
na bolha 1j + . Analisando-se os resultados experimentais para a queda de pressão
obtidos na literatura e apresentados no Capítulo 2 nota-se que a diferença de
pressão 1LS GBjxjP P +− não é pequena. Aliás, essa diferença de pressão foi definida
como sendo MixPΔ . Substituindo-se 1LS GBjxjP P +− por MixPΔ na Equação (3.45) obtém-
se:
( ) ( )Mix 1 1Lxj LBj LSxj HS HBxj jP A m U U F F+ +
Δ ⋅ = − − − . (3.46)
A Equação (3.45) pode ser utilizada para relacionar as pressões no pistão e na
bolha à sua frente. Porém, o primeiro termo do lado direito da Equação (3.46) não é
de fácil implementação, pois está relacionado à expansão do filme para o pistão.
Como a modelagem apresentada é unidimensional, alguns efeitos devidos à
geometria tridimensional da expansão não são capturados. Assim, seria necessária
a utilização de constantes de ajuste na Equação (3.46). Na prática, essas constantes
não existem para a geometria da traseira da bolha em movimento, o que torna essa
abordagem inviável.
Para contornar o problema da Equação (3.46), o presente trabalho utilizará a
abordagem apresentada por Taitel e Barnea (1990b). Os autores identificam uma
forma alternativa de calcular a queda de pressão devido à mistura. Apesar do
modelo considerar o escoamento estacionário, pode ser utilizado como uma
aproximação no modelo de seguimento de pistões.
A Figura 3.7 apresenta um volume de controle que envolve o acoplamento
entre a frente do pistão j e a frente da bolha 1j + . Desprezando-se os termos de
acúmulo e os fluxos de quantidade de movimento no filme de líquido, a equação de
conservação da quantidade de movimento é dada por:
( )1 1 1 1 1 1 1 sen 0Bj LS GBj LBj LBj Bj L LBj BjxjF P P A S L AR L gτ ρ θ+ + + + + + += − − − =∑ . (3.47)
Capítulo 3 Modelagem matemática 47
Na Equação (3.47) o segundo e terceiro termos são referentes ao atrito na parede
do filme de líquido e a força gravitacional, 1LBjτ + é a tensão cisalhante na parede e
1LBjS + é o perímetro molhado médio no filme.
LBjx j
yj+1
PGBj+1
PLS|xj
Figura 3.7 – Representação do volume de controle no acoplamento entre a frente do
pistão e a frente da bolha
A Equação (3.47) apresenta uma relação entre as pressões na frente do
pistão j e na bolha 1j + . Essa relação é baseada apenas em forças de atrito e
gravitacionais, que são de compreensão mais fácil que o fenômeno de ressalto
hidráulico. Por isso, a forma da Equação (3.47) foi a escolhida para esse trabalho.
Evidenciando-se o termo da pressão que atua no pistão de líquido, obtém-se:
1 1 1 1 1Mix 1
senLBj LBj Bj L LBj BjLS GBjxj
S L AR L gP P P
Aτ ρ β+ + + + +
+
+Δ = − = . (3.48)
Esse resultado foi demonstrado originalmente por Taitel e Barnea (1990b). Os
resultados das Equações (3.43) e (3.48) devem ser substituídos na Equação (3.40),
o que resulta:
( )
( )1 1 1 1
1 1 sen
Sj GBj GBj LSj Sj LBj LBj Bj
L LSj Sj LBj Bj
F P P A DL S L
A R L R L g
τ π τ
ρ θ
+ + + +
+ +
= − − − −
− +
∑. (3.49)
Capítulo 3 Modelagem matemática 48
Substituindo-se as forças da Equação (3.49), a Equação (3.39) se torna:
( ) ( )1 1 1 11 1 sen
2
GBj GBj LSj Sj LBj LBj BjLSj Sj LBj Bj
L L L
LSj Sj j j LSjLSj Sj Sj LSj
P P DL S LR L R L g
A A
dU L dx dy dRR L L U
dt dt dt dt
τ π τθ
ρ ρ ρ+ + + +
+ +
−− − − + =
⎡ ⎤⎛ ⎞= + + −⎢ ⎥⎜ ⎟
⎝ ⎠⎣ ⎦
. (3.50)
A Equação (3.50) é a equação de conservação da quantidade de movimento do
líquido no pistão. A Equação (3.50) pode ser escrita como:
( )
12 21 1 1 1
1 1
2 2
sen2
LSj Sj BjGBj GBj L LSj Sj LSj L LSj LBj L LBj LBj
Sj j j LSjLSj Sj LBj Bj L Sj LSj L
dU L LP P R L C U C s U
dt D DL dx dy dR
R L R L g L Udt dt dt
ρ ρ ρ
ρ θ ρ
++ + + +
+ +
− = + + +
⎡ ⎤⎛ ⎞+ + + + −⎢ ⎥⎜ ⎟
⎝ ⎠⎣ ⎦
, (3.51)
onde as tensões de cisalhamento na parede do pistão e do filme de líquido, LSjτ e
1LBjτ + , foram substituídas por expressões em função dos fatores de atrito de Fanning
e das velocidades:
21
2
LSjLSj
L LSj
CU
τ
ρ= , (3.52)
e:
11
21
12
LBjLBj
L LBj
CU
τ
ρ
++
+
= . (3.53)
O cálculo dos fatores de atrito é apresentado na seção 3.2.3. O perímetro médio do
filme de líquido, 1LBjS + , foi adimensionalizado através de:
Capítulo 3 Modelagem matemática 49
11
LBjLBj
Ss
Dπ+
+ = . (3.54)
O cálculo dos coeficientes de atrito e do perímetro médio do filme de líquido
adimensional é apresentado na seção 3.2.3.
A Equação (3.51) é a forma final da equação de conservação da quantidade
de movimento aplicada à fase de líquido no pistão. Essa equação indica que a
diferença entre as pressões do gás no interior de duas bolhas consecutivas é devido
a vários efeitos, representados pelo lado direito da equação. O primeiro termo é a
variação da quantidade de movimento do líquido no interior do pistão, ou, também
chamado de inércia do pistão, os segundo e terceiro termos são as forças de atrito
que atuam na parede do pistão e do filme de líquido, respectivamente, o quarto
termo é a força gravitacional no líquido do pistão e filme e o quinto termo é a
variação da quantidade de movimento no pistão devido à variação temporal da
fração de líquido no pistão.
As Equações (3.30) e (3.51) formam um sistema acoplado de duas equações
e duas incógnitas, LSjU e GBjP , para cada célula, que é resolvido numericamente a
cada passo de tempo. O método numérico utilizado é apresentado no Capítulo 4.
3.2 Equações auxiliares
As equações de conservação da massa e da quantidade de movimento,
Equações (3.30) e (3.51) estão acopladas. A sua resolução requer que sejam
conhecidos os deslocamentos das células em cada passo de tempo, requerendo o
cálculo da variação de jx e jy , jdx dt e jdy dt , em função de LSjU e GBjP . Essa
seção apresenta as equações necessárias a esse cálculo. Além disso, são
apresentadas equações auxiliares para o cálculo do fator de atrito e da velocidade
do filme de líquido, que são necessárias para a resolução da equação de
conservação da quantidade de movimento.
Capítulo 3 Modelagem matemática 50
3.2.1 Deslocamento da frente da bolha
A variação temporal da fronteira jy , jdy dt , é a velocidade de translação da
bolha alongada. Como apresentado no capítulo 2 (revisão bibliográfica), a
velocidade de translação da bolha é calculada em função da velocidade do pistão,
LSjU , como:
( )( )0 1 1Tj LSj jU C U C gD h= + + , (3.55)
onde 0C , 1C e jh são constantes. A constante 0C quantifica a influência da
velocidade do pistão, e um valor de 0 ~ 1,2C é aceito para escoamentos turbulentos
na vertical, enquanto 0 ~ 2C é aceito para escoamentos laminares. A constante 1C é
relacionada à velocidade que uma bolha teria se o líquido a sua frente fosse
estagnado. A constante jh é o fator de esteira, que aumenta a velocidade da bolha
quando ela está próxima da bolha precedente, devido à mistura que ocorre na parte
inicial do pistão. Essa constante é modelada por:
exp Sjj w w
Lh a b
D⎛ ⎞
= −⎜ ⎟⎝ ⎠
, (3.56)
onde wa e wb são duas constantes de ajuste experimentais.
Após o cálculo de LSjU , calcula-se a velocidade de translação da bolha TjU
através da Equação (3.55), e, como:
jTj
dyU
dt= , (3.57)
e resolve-se numericamente para o deslocamento da bolha ao longo de um passo
de tempo.
Capítulo 3 Modelagem matemática 51
3.2.2 Deslocamento da traseira da bolha
O deslocamento da traseira da bolha também deve ser conhecido para
completar o movimento da bolha. Essa equação deve respeitar a conservação da
massa. Assim, utiliza-se a Equação (3.20). Evidenciando-se a velocidade de
deslocamento da traseira da bolha, 1jdx dt− , obtém-se:
( ) 1 1
1
1
1 11 1
1
1
2 2
2 2
j LBj Sj LSj Sj LSjGBj GSj Bj
j
GBj GSj
GBj Bj GBj Sj GSj Sj GSjGSj GSj GSj GSj
GBj GBj GBj
GBj GSj
dy dR L dR L dRR R Ldx dt dt dt dt
dt R R
dP L R L R L RR U R U
dt P P PR R
− −
−
−
− −− −
−
−
− − − −= +
−
⎛ ⎞+ + + −⎜ ⎟⎜ ⎟
⎝ ⎠+−
. (3.58)
No presente trabalho a fração de líquido de cada filme, LBjR , é constante ao longo do
tempo (porém cada bolha possui uma fração de líquido diferente), e, assim, pode-se
escrever:
( ) 1 11
1
1 1 1 1
1
1
2 2 2 2
j Bj GBj GBjGBj GSj GSj GSj GSj GSj
j GBj
GBj GSj
Sj GSj Sj GSj GBj Sj LSj Sj LSj
GBj GBj
GBj GSj
dy L R dPR R R U R U
dx dt P dtdt R R
L R L R dP L dR L dRP P dt dt dt
R R
− −−
−
− − − −
−
−
− + + −= +
−
⎛ ⎞+ − −⎜ ⎟⎜ ⎟
⎝ ⎠+−
. (3.59)
A Equação (3.59) também é discretizada e resolvida numericamente a cada passo
de tempo.
3.2.3 Definição do atrito
Na modelagem da equação de conservação da quantidade de movimento no
pistão foram definidos os coeficientes de atrito de Fanning, de acordo com as
Equações (3.52) e (3.53). Esses coeficientes são calculados pela relação de Blasius
como:
Capítulo 3 Modelagem matemática 52
1
0,25
16 Re , se Re 1000
0,079Re , se Re 1000LSj LSj
LSjLSj LSj
C−
−
⎧ <⎪= ⎨>⎪⎩
, (3.60)
e:
0,250,079 ReLBj LBjC −= , (3.61)
onde ReLSj é o número de Reynolds do pistão de líquido, definido como:
Re L LSjLSj
L
U Dρμ
= , (3.62)
e ReLBj é o número de Reynolds do filme de líquido:
Re L LBj HjLBj
L
U Dρμ
= , (3.63)
calculado a partir do diâmetro hidráulico do filme, que é dado por:
4 LBj LBjHj
LBj LBj
R A RD D
s D sπ= = , (3.64)
e da velocidade média do filme de líquido, definida na seção seguinte.
Na equação de conservação da quantidade de movimento também foi
utilizado o perímetro molhado adimensional do filme de líquido, LBjs , que é calculado
de forma aproximada por:
1, se 90º0,5269 0, 2365, se 0ºLBj
LBj
sR
θθ
⎧= ⎨ +⎩
∼∼
, (3.65)
Capítulo 3 Modelagem matemática 53
conforme a equação proposta por Fagundes Netto (1999) para a definição de LBjs no
escoamento horizontal. Essa equação é válida para frações de líquido no filme entre
0,2 e 0,8.
3.2.4 Velocidade do filme de líquido
A velocidade média do filme de líquido é necessária para o cálculo do atrito
no filme. A Figura 3.8 apresenta um volume de controle que engloba a frente da
bolha, desde o pistão de líquido até uma coordenada do filme em que a velocidade é
igual à velocidade média.
UTj
ULSj
ULBj
jésima célula
Figura 3.8 – Representação do volume de controle na interface entre o pistão e a
bolha que o segue
Esse volume de controle é indeformável e se movimenta com velocidade igual
à velocidade de translação da bolha, TjU . Assim, a equação de conservação da
massa para esse volume de controle é:
( ) ( )Lyj L LSj LSj Tj L LBj LBj Tjm AR U U AR U Uρ ρ= − = − (3.66)
Isolando-se a velocidade do filme, LBjU , obtém-se:
( )TLSj
LBj Tj LSj jLBj
RU U U U
R= + − (3.67)
Capítulo 3 Modelagem matemática 54
As Equações (3.57), (3.59), (3.60), (3.61), (3.65) e (3.67) são as equações
não acopladas que devem ser calculadas após o cálculo de LSjU e GBjP em um
tempo novo, e completam o sistema a ser resolvido numericamente para todas as
células do escoamento. No entanto, existem algumas singularidades no modelo, que
resultam da entrada, saída e coalescência das bolhas. A seção a seguir apresenta a
modelagem dessas singularidades.
3.3 Modelagem das singularidades
Diversas singularidades do modelo de seguimento de pistões devem ser
abordadas. Nesta seção é apresentada a modelagem das coalescências entre
bolhas. A entrada e saída de bolhas, bem como as condições de contorno utilizadas
no modelo são apresentas em um capítulo específico, o Capítulo 4.
3.3.1 Coalescência de bolhas
A coalescência de bolhas ocorre no momento em que a frente de uma bolha
toca a traseira da bolha a sua frente, como mostrado na Figura 3.9. Nesse instante,
o comprimento da nova bolha j (representado pelo sobrescrito N ), representada na
Figura 3.9b deve ser calculado em função dos parâmetros antigos (representados
pelo sobrescrito O ). Considera-se que toda a massa de líquido contida no pistão
1j − antigo e nos filmes das bolhas j e 1j + antigos se re-arranja no pistão 1j −
novo e no filme da bolha j nova, ou:
1 1 1O O O N NLSj LBj LBj LSj LBjM M M M M− + −+ + = + , (3.68)
onde LM representa a massa de líquido. Além disso, da geometria das células, tem-
se que:
1 1 1O O O N NSj Bj Bj Sj BjL L L L L− + −+ + = + . (3.69)
Capítulo 3 Modelagem matemática 55
Considerando-se, ainda, que a fração de líquido no pistão 1j − permanece
inalterada:
1 1O NLSj LSjR R− −= , (3.70)
e substituindo-se as massas por L LM ALρ= , onde L é um comprimento, o
comprimento novo da bolha j se torna:
( ) ( )1 1 1 1
1
O O O O O OBj LBj LSj Bj LBj LSjN
Bj N OLBj LSj
L R R L R RL
R R− + + −
−
− + −=
−. (3.71)
A única incógnita na Equação (3.71) é a nova fração de líquido no filme da bolha j .
Considera-se, então, que essa fração de líquido é igual à fração de líquido mínima
entre as frações de líquido antigas nos filmes das bolhas j e 1j + , ou:
( )1Mínimo ,N O OLBj LBj LBjR R R += . (3.72)
Após resolvida a Equação (3.71), o comprimento novo do pistão 1j − é:
1 1 1N O O O NSj Sj Bj Bj BjL L L L L− − += + + − . (3.73)
Capítulo 3 Modelagem matemática 56
LBj-1
LBj-1
O
N
O
N
O
N
O O
N
LSj-1 LBj
LBj
LBj+1 LSj+1
LSj
j+1OjO
jN j+1N
j+2Oj-1O
j-1N
a)
b)
Figura 3.9 – Coalescência de bolhas
Capítulo 4 Método de solução numérica 57
4 MÉTODO DE SOLUÇÃO NUMÉRICA
Este capítulo apresenta o método utilizado para simular o escoamento em
golfadas a partir da modelagem matemática apresentada. Inicialmente é
apresentada a discretização em relação ao tempo das equações acopladas e a sua
implementação em um sistema linear de equações que, depois de resolvido
numericamente, retorna os valores de LSjU e GBjP em um tempo novo para todas as
células do tubo. São discretizadas também as equações auxiliares utilizadas para o
cálculo do deslocamento das frentes dos pistões ( jx ) e das bolhas ( jy ). O algoritmo
geral da simulação também é apresentado.
Na seqüência, são apresentadas as condições iniciais e de contorno, e os
processos realizados na entrada e saída de bolhas e pistões na tubulação. Ênfase é
dada à forma como são calculados os parâmetros de cada nova célula que entra na
tubulação, pois nesse cálculo são incorporados os efeitos intermitentes na
simulação.
Por fim, são apresentados os conceitos das sondas virtuais, e os resultados
de simulações preliminares, que tem como objetivo definir alguns parâmetros do
simulador.
4.1 Discretização das equações acopladas
Esta seção apresenta a discretização, através do método das diferenças
finitas, das equações diferenciais acopladas obtidas na modelagem matemática. As
duas equações principais a serem resolvidas são a conservação da massa de
líquido e gás na célula e da quantidade de movimento no líquido presente no pistão,
repetidas a seguir por conveniência:
( ) ( ) ( )
( ) ( )
111
1
11
1
1 1 12 2
1 1
LBj LSj LSjGBj Sj SjLSj LSj Bj
GBj GBj GBj
LSj LSjDSj DSj
LSj LSj
R R RdP L LU U L
dt P P P
R RU U
R R
−−−
−
−−
−
⎡ ⎤− − −− = + + +⎢ ⎥
⎢ ⎥⎣ ⎦
− −+ −
, (3.30)
Capítulo 4 Método de solução numérica 58
e:
( )
12 21 1 1 1
1 1
2 2
sen2
LSj Sj BjGBj GBj L LSj Sj LSj L LSj LBj L LBj LBj
Sj j j LSjLSj Sj LBj Bj L Sj LSj L
dU L LP P R L C U C s U
dt D DL dx dy dR
R L R L g L Udt dt dt
ρ ρ ρ
ρ θ ρ
++ + + +
+ +
− = + + +
⎡ ⎤⎛ ⎞+ + + + −⎢ ⎥⎜ ⎟
⎝ ⎠⎣ ⎦
. (3.51)
As Equações (3.30) e (3.51) são reescritas de forma simplificada como:
1 11
12 2GBj Bj GBj Sj GSj Sj GSj
LSj LSj DSjGBj GBj GBj
dH L R L R L RU U U
dt H H H− −
−−
⎛ ⎞− = + + + Δ⎜ ⎟⎜ ⎟
⎝ ⎠, (4.1)
e:
( )21 1
12LSj SjGBj GBj LSj Sj LSj LSj Sj Gj j
L
dU LH H R L C U P P I
dt D ρ+ +− = + + Δ + Δ + Δ , (4.2)
onde as pressões no interior das bolhas, GBP , foram substituídas pelo fator de
pressão:
GBjGBj
L
PH
ρ= . (4.3)
A definição da Equação (4.3) é necessária para a estabilidade numérica, pois se
apenas o valor da pressão fosse utilizado, a matriz a ser resolvida seria mal
condicionada, devido à grande diferença na ordem de grandeza de seus elementos.
Além da Equação (4.3), as frações de líquido foram substituídas pelas frações de
gás, que são calculadas como:
1 11 , 1 e 1GBj LSj GSj LSj GSj LSjR R R R R R− −= − = − = − . (4.4)
Os termos novos que surgem nas Equações (4.1) e (4.2) são:
Capítulo 4 Método de solução numérica 59
( ) ( )1
11
1 1LSj LSjDSj DSj DSj
LSj LSj
R RU U U
R R−
−−
− −Δ = − , (4.5)
que agrupa os termos referentes à velocidade do gás no interior do pistão,
1 21 1 1 12 Bj
Sj LBj LBj L LBj
LP C s U
Dρ +
+ + + +Δ = , (4.6)
que representa a queda de pressão devido ao atrito entre a parede do tubo e o
líquido no filme,
( )1 1 senGj LSj Sj LBj Bj LP R L R L gρ θ+ +Δ = + , (4.7)
que é a queda de pressão devido a força gravitacional que atua no líquido presente
no pistão de líquido j e no filme de líquido 1j + , e:
2Sj j j LSj
j Sj LSj L
L dx dy dRI L U
dt dt dtρ
⎡ ⎤⎛ ⎞Δ = + −⎢ ⎥⎜ ⎟
⎝ ⎠⎣ ⎦, (4.8)
que representa a variação da quantidade de movimento no interior do pistão devido
à variação temporal da fração de líquido no pistão.
As Equações (4.1) e (4.2) são discretizadas em relação ao tempo, pois a
discretização em relação ao espaço já foi realizada implicitamente na modelagem
matemática, devido à integração das propriedades ao longo de cada célula. A
discretização em relação ao tempo é realizada através do método de diferenças
finitas, utilizando o esquema semi-implícito de Crank-Nicholson, apenas para as
variáveis principais, LSU e GBP , e as demais variáveis são utilizadas com os seus
valores já conhecidos (em um tempo antigo). Assim, o esquema de Crank-Nicholson
utilizado não é puro. Dessa forma, discretizando-se a Equação (4.1), obtém-se:
1 1 1 1
12 2 2 2
N O N O N OLSj LSj LSj LSj GBj GBj Bj GBj Sj GSj Sj GSj
DjO O OGBj GBj GBj
U U U U H H L R L R L RU
t H H H− − − −
−
⎛ ⎞+ + −− = + + + Δ⎜ ⎟⎜ ⎟Δ ⎝ ⎠
, (4.9)
Capítulo 4 Método de solução numérica 60
onde os sobrescritos N e O indicam que o valor das variáveis é avaliado nos
tempos novo ( t t+ Δ ) e antigo ( t ), respectivamente. Esses índices são utilizados
apenas para identificar as variáveis principais, U e P , e as demais variáveis, que
não apresentam nenhum dos dois índices, devem ser avaliadas no tempo antigo.
Reorganizando-se a Equação (4.9) de forma que as variáveis avaliadas no
tempo novo fiquem isoladas, obtém-se:
1 11 1
1
1 1
1
2
22
Bj GBj Sj GSj Sj GSjN N N OLSj GBj LSj LSjO O O
GBj GBj GBj
OBj GBj Sj GSj Sj GSj GBjO
LSj DjOGBj
L R L R L RU H U U
H t H t H t
L R L R L R HU U
t t t H
− −− −
−
− −
−
⎛ ⎞− + + + + = −⎜ ⎟⎜ ⎟Δ Δ Δ⎝ ⎠
⎛ ⎞− + + + − Δ⎜ ⎟⎜ ⎟Δ Δ Δ⎝ ⎠
. (4.10)
Da mesma forma, a discretização da equação da quantidade de movimento,
Equação (4.2), resulta em:
( )
1 1
1
2 2
12
N O N O N OGBj GBj GBj GBj LSj LSj
LSj Sj
Sj N OLSj LSj LSj Sj Gj j
L
H H H H U UR L
t
LC U U P P I
D ρ
+ +
+
⎛ ⎞+ + −− = +⎜ ⎟⎜ ⎟Δ⎝ ⎠
+ + Δ + Δ + Δ
, (4.11)
que, separando-se os termos avaliados no tempo novo dos termos avaliados no
tempo antigo, resulta em:
( )
1
1 1
24
2 2
LSj Sj SjN N O N OGBj LSj LSj LSj GBj GBj
OLSj Sj LSjO
GBj Sj Gj jL
R L LH U C U H H
t D
R L UH P P I
t ρ
+
+ +
⎛ ⎞− + + + = −⎜ ⎟Δ⎝ ⎠
− + − Δ + Δ + ΔΔ
. (4.12)
As Equações (4.10) e (4.12) referem-se à célula genérica j . Porém, elas
devem ser aplicadas a todas as células no interior do domínio computacional, de 1
até n . Para a célula 1 obtém-se:
Capítulo 4 Método de solução numérica 61
1 1 1 1 0 01 1 0 0
1 1 0
1 1 1 1 0 0 11 1
0
1 1 11 1 1 1 2 1
1 12
2
2 2 ;
2 4
2
N N N OB GB S GS S GSGB LS LS LSO O O
GB GB GB
OO B GB S GS S GS GBLS DO
GB
N N O N OLS S SGB LS LS LS GB GB
O LS SGB
L R L R L RH U U UH t H t H t
L R L R L R HU Ut t t H
R L LH U C U H Ht D
R LH
⎛ ⎞+ + + = + −⎜ ⎟Δ Δ Δ⎝ ⎠
− + + + − ΔΔ Δ Δ
⎛ ⎞− + + + = −⎜ ⎟Δ⎝ ⎠
− + ( )12 1 1
2OLS
S GL
U P P It ρ
⎧⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪
− Δ + Δ + Δ⎪ Δ⎩
. (4.13)
É importante notar que a velocidade do líquido no pistão 0, 0NLSU , aparece do lado
direito da equação, pois é uma condição de contorno do modelo, e deve ser
conhecida a priori. Da mesma forma, para a célula 2 obtém-se:
2 2 2 2 1 11 2 2 1
2 2 1
2 2 2 2 1 1 22 2
1
2 2 22 2 2 2 3 2
23
2
2 2 ;
2 4
2
N N N OB GB S GS S GSLS GB LS LSO O O
GB GB GB
OO B GB S GS S GS GBLS DO
GB
N N O N OLS S SGB LS LS LS GB GB
O LS SGB
L R L R L RU H U UH t H t H t
L R L R L R HU Ut t t H
R L LH U C U H Ht D
R LH
⎛ ⎞− + + + + = −⎜ ⎟Δ Δ Δ⎝ ⎠
− + + + − ΔΔ Δ Δ
⎛ ⎞− + + + = −⎜ ⎟Δ⎝ ⎠
− + ( )2 23 2 2
2OLS
S GL
U P P It ρ
⎧⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪
− Δ + Δ + Δ⎪ Δ⎩
, (4.14)
e assim por diante até a última célula no interior do tubo, célula n , que assume a
forma:
1 11 1
1
1 1
1
2
2 2 ;
2 4
N N N OBn GBn Sn GSn Sn GSnLSn GBn LSn LSnO O O
GBn GBn GBn
OO Bn GBn Sn GSn Sn GSn GBnLSn DnO
GBn
N N O OLSn Sn SnGBn LSn LSn LSn GBn GBn
L R L R L RU H U UH t H t H t
L R L R L R HU Ut t t H
R L LH U C U H Ht D
− −− −
−
− −
−
⎛ ⎞− + + + + = −⎜ ⎟Δ Δ Δ⎝ ⎠
− + + + − ΔΔ Δ Δ
⎛ ⎞− + + = −⎜ ⎟Δ⎝ ⎠
( )
1
1 12 2
N
OO LSn Sn LSnGBn Sn Gn n
L
R L UH P P It ρ
+
+ +
⎧⎪⎪⎪⎪⎪⎨⎪ −⎪⎪⎪
− + − Δ + Δ + Δ⎪ Δ⎩
. (4.15)
Capítulo 4 Método de solução numérica 62
Na Equação (4.15), nota-se a outra condição de contorno do modelo, que é o fator
de pressão na saída do tubo, 1NGBnH + . A definição das condições de contorno é
apresentada nas seções seguintes.
O sistema de equações obtido pela aplicação das Equações (4.13) a (4.15)
pode ser representado através de:
A X B⋅ = (4.16)
onde A é a matriz de coeficientes, X é o vetor de incógnitas e B é o vetor de
termos independentes, que são escritos como:
Capítulo 4 Método de solução numérica 63
1 1 1 1 0 0
1 1 0
1 1 11 1
1 1
1
2 1 0 0
21 4 0 0
,20 0 1
20 0 1 4
B GB S GS S GSO O OGB GB GB
OLS S SLS LS
Bn GBn Sn GSn Sn GSnO O OGBn GBn GBn
OLSn Sn SnLSn LSn
L R L R L RH t H t H t
R L LC Ut D
AL R L R L RH t H t H t
R L LC Ut D
X
− −
−
⎡ ⎤+ +⎢ ⎥Δ Δ Δ⎢ ⎥⎢ ⎥
− +⎢ ⎥Δ⎢ ⎥= ⎢ ⎥
⎢ ⎥⎢ ⎥+ +
Δ Δ Δ⎢ ⎥⎢ ⎥⎢ ⎥− +
Δ⎣ ⎦
=
( )
1 1 1 1 0 0 10 0 1 1
0
1 1 1 11 2 2 1 1
1
1 11
2 2
2 2
,2
ON O O B GB S GS S GS GBLS LS LS DO
GBN O
O OGB LS S LSGB GB S GN
LLS
NGBn O O Bn GBn Sn GSn Sn GSn GB
LSn LSnNLSn
L R L R L R HU U U Ut t t H
H R L UH H P P ItU
BH L R L R L R HU U
t t tU
ρ
− −−
+ − + + + − ΔΔ Δ Δ
⎡ ⎤− + − Δ + Δ + Δ⎢ ⎥ Δ⎢ ⎥
⎢ ⎥ =⎢ ⎥⎢ ⎥ − + + +⎢ ⎥ Δ Δ Δ⎣ ⎦
( )
1
1 1 1
2
2 2
On
DnOGBn
OO N O LSn Sn LSnGBn GBn GBn Sn Gn n
L
UH
R L UH H H P P It ρ
−
+ + +
⎡ ⎤⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥
− Δ⎢ ⎥⎢ ⎥⎢ ⎥
− − + − Δ + Δ + Δ⎢ ⎥Δ⎣ ⎦
(4.17)
Capítulo 4 Método de solução numérica 64
A matriz de coeficientes, A , é tri-diagonal, ou seja, possui elementos não
nulos apenas na diagonal principal e nas duas diagonais secundárias. Por isso, o
sistema é resolvido de forma simples através do algoritmo TDMA (Tridiagonal Matrix
Algorithm). O TDMA é um método direto, e um sistema com 2n equações é
resolvido em apenas 4n operações, onde n é o número de bolhas no interior da
tubulação. No modelo apresentado, o sistema representado pela Equação (4.17) é
resolvido a cada passo de tempo, e resulta nos valores de NLSU e N
GBP para cada
pistão e bolha no domínio computacional.
4.2 Discretização das equações auxiliares
As equações auxiliares são não-acopladas e devem ser resolvidas após o
conhecimento dos valores de NLSjU e N
GBjP , e determinam o deslocamento das frentes
dos pistões e bolhas. Essas equações foram apresentadas no capítulo de
modelagem matemática, e são reescritas a seguir por conveniência:
jTj
dyU
dt= , (4.18)
e:
( ) 1 11
1
1 1 1 1
1
1
2 2 2 2
j Bj GBj GBjGBj GSj GSj GSj GSj GSj
j GBj
GBj GSj
Sj GSj Sj GSj GBj Sj LSj Sj LSj
GBj GBj
GBj GSj
dy L R dPR R R U R U
dx dt P dtdt R R
L R L R dP L dR L dRP P dt dt dt
R R
− −−
−
− − − −
−
−
− + + −= +
−
⎛ ⎞+ − −⎜ ⎟⎜ ⎟
⎝ ⎠+−
. (4.19)
A discretização das equações auxiliares através do método de diferenças finitas
proporciona:
N N Oj Tj jy U t y= Δ + , (4.20)
Capítulo 4 Método de solução numérica 65
para o deslocamento das frentes das bolhas ao longo de um passo de tempo, tΔ , e:
( ) ( ) ( )
( ) ( )
11 1
1 11
1 11 1
1
1
2 2
2 2
O OSj SjO N N O N O
GBj GSj Tj LSj LSj LSj LSjN Oj j O
GBj GSj
O O O O OBj GBj Sj GSj Sj GSjN O O O O O
GBj GBj GSj GSj GSj GSjN N NGBj GBj GBj
OGBj GSj
L LR R tU R R R R
x xR R
L R L R L RP P t R U R U
P P PR R
−− −
− −−
− −− −
−
−
− Δ − − − −= + +
−
⎛ ⎞− + + + Δ −⎜ ⎟⎜ ⎟
⎝ ⎠+−
, (4.21)
para o deslocamento das traseiras das bolhas ao longo de um passo de tempo.
4.3 Condições iniciais e de contorno
O processo numérico de resolução do sistema de equações (4.17) consiste
em uma marcha no tempo a partir de um instante inicial. Portanto, as propriedades
do escoamento e os valores de LSU e GBP devem ser conhecidos para todas as
células no interior da tubulação nesse instante inicial. A princípio, no início da
simulação, poder-se-ia ter diversas células no interior da tubulação. No entanto, isso
traz o inconveniente de como os parâmetros relativos a essas células seriam
conhecidos em 0t = .
Portanto, para tornar a condição inicial mais simples, no presente trabalho
considera-se que, em 0t = , a tubulação está completamente cheia de líquido e a
primeira bolha está posicionada em 0z = pronta para entrar na tubulação. Dessa
forma, só existe uma célula no interior da tubulação, com um pistão do tamanho da
própria tubulação, e a bolha começando a entrar. Assim, os parâmetros dessa única
célula podem ser facilmente calculados através de metodologias para o escoamento
monofásico que está ocorrendo nesse instante.
Essa consideração é compatível com um caso real em que a tubulação está
escoando somente líquido e, em determinado momento, inicia-se o escoamento de
gás.
Quanto às condições de contorno, nota-se que o sistema de equações
apresentado, Equação (4.17), requer o conhecimento da velocidade do líquido no
Capítulo 4 Método de solução numérica 66
pistão 0, 0NLSU , e a pressão na bolha 1j n= + , 1GBnP + . Os valores dessas duas
variáveis devem ser fornecidos para o sistema de equações em todos os passos de
tempo, e podem ser constantes ou variar ao longo da simulação.
A atribuição de valores para 0NLSU e 1GBnP + , a princípio, parece ser uma tarefa
simples. No entanto, existem alguns inconvenientes.
O primeiro deles é que se deve ter sempre um pistão na entrada da tubulação
(o pistão 0) e uma bolha na saída da tubulação (a bolha n+1). No entanto, pela
característica intrínseca de alternância do escoamento, isso não é possível.
Portanto, devem ser tomadas medidas para conviver com essas alternâncias, como
as apresentadas na seção seguinte.
O segundo inconveniente é que o valor de 0NLSU é dependente das
velocidades superficiais, e, consequentemente, das vazões, de líquido e gás na
entrada da tubulação. Dessa forma, a definição de um valor para 0NLSU resulta
indiretamente na definição da vazão de fluidos na entrada do tubo. Portanto, deve-se
ter cuidado nos valores utilizados para essa variável. A próxima seção apresenta a
forma como 0NLSU é calculado.
Em resumo, a Figura 4.1 apresenta as condições iniciais e de contorno da
simulação.
0
t
zLCon
diçã
o de
con
torn
o pa
ra a
vel
ocid
ade:
U(t,
0)=c
onhe
cido
LS
Condição
decontorno
paraa
pressão:P(t,L)=conhecido
GB
Condição inicial (escoamento monofásico):U (0,z)=conhecidoP (0,z)=conhecido
LS
GB Figura 4.1 – Representação esquemática do escoamento ao longo do tubo em
diferentes instantes de tempo e a definição das condições de contorno e inicial
Capítulo 4 Método de solução numérica 67
4.4 Processo de entrada de bolhas e pistões no domínio de cálculo
Os modelos de seguimento de pistões não são projetados para prever a
transição de nenhum outro padrão de escoamento para o padrão de golfadas. Da
mesma forma, no interior da tubulação, não são previstas transições do escoamento
em golfadas para outros padrões. Isso implica que, tanto na entrada da tubulação
como na saída, deve-se garantir a ocorrência do padrão de golfadas. A seguir são
definidas as formas como as bolhas e pistões são inseridas ou retiradas da
tubulação.
A Figura 4.2 apresenta a evolução temporal das células na entrada do tubo.
No instante nt , a bolha 1 está com sua frente, 1y , exatamente em 0z = . Nesse
instante, as seguintes propriedades para a célula 0 devem ser definidas: 0LSU , 0SL ,
0BL , 0LSR e 0GSR (a forma de cálculo desses parâmetros é apresentada na subseção
seguinte). A bolha 1 já faz parte do domínio de cálculo, e o seu deslocamento,
comprimento e demais parâmetros são regidos pela solução do sistema de
equações.
No tempo 1nt + , a bolha 1 se deslocou, e somente o valor de 0LSU deve ser
atualizado (pois é uma condição de contorno para o sistema de equações a cada
passo de tempo). Esse procedimento continua no tempo 2nt + , quando a traseira da
bolha já está dentro do tubo, e, durante esse período, o comprimento do pistão 0
permanece constante.
No tempo 3nt + o pistão 0 continua entrando no tubo, e seu comprimento ainda
é constante. Esse procedimento continua até o tempo 4nt + , quando a frente da bolha
0 chega à coordenada 0z = . Nesse instante, as células 1 até n são renumeradas, e
a célula 0 passa a ser a célula 1. Em seguida, o processo volta para a posição
descrita em nt , quando a nova célula 0 deve ser criada e os valores das
propriedades da célula 0 devem ser definidos.
Capítulo 4 Método de solução numérica 68
0
t
zz=0
tn
tn+1
tn+2
tn+3
tn+4
PGB1PGB0ULS0
y0
y0
y0
y0
y0
x0
x0
x0
x0
x0
LS0
LS0
LS0
LS0
LS0
LB0
LB0
LB0
LB0
LB0
y1
y1
y1
y1
y1
x1
x1
x1
x1
ULS0
ULS0
ULS0
ULS0
ULS1 ULS2
ULS2
ULS2
ULS1
ULS1
ULS1
ULS1
PGB0
PGB0
PGB0
PGB0
PGB2
PGB2
PGB2
PGB2
PGB1
PGB1
PGB1
PGB1
Figura 4.2 – Representação esquemática do processo de entrada de pistão de
líquido e bolhas de gás no domínio computacional
Condições a serem respeitadas na entrada do domínio (pseudo-algoritmo):
• Enquanto a bolha 1 está entrando no tubo ( 1 0y > e 0 0x < )
o Todas as propriedades da célula são obtidas da resolução do sistema de
equações
o 0BL = cte
o 0SL = cte
o 0LSU = informado
• Enquanto o pistão 0 está entrando ( 0 0x > e 0 0y < )
o 1BL = calculado na solução do modelo
Capítulo 4 Método de solução numérica 69
o 0SL = cte
o 0BL = cte (Considera-se que o valor de 0BL é aquele que foi determinado
quando o nariz da bolha 1j = atingiu 0z = . Por esta razão a pressão da
bolha está sendo atualizada, mas o seu comprimento 0BL não)
o 0LSU = informado
• No instante em que a bolha começa a entrar no tubo ( 0 0Ny > e 0 0Oy < )
o Renumeração das células
o Leitura/cálculo das propriedades para a célula 0
4.4.1 Cálculo das variáveis na entrada de uma nova célula
Quando uma nova célula é inserida na entrada da tubulação (como
apresentado na seção 4.4), é necessário o cálculo dos parâmetros geométricos da
célula, isto é: 0SL , 0BL , 0LSR e 0GBR , além da velocidade do líquido no pistão, 0LSU .
Na parte inicial dessa subseção será apresentada a forma como essas variáveis
foram obtidas no presente trabalho, a partir de distribuições baseadas em dados
experimentais. Ao final da subseção são apresentadas formas alternativas para
cálculo desses valores quando não se tem dados experimentais disponíveis e
utilizam-se correlações de fechamento.
Velocidade do líquido no pistão
O valor da velocidade do líquido no pistão, 0LSU , é calculado em função das
velocidades superficiais do gás e do líquido na entrada, 0GJ e 0LJ , como:
( )0 0 0 0 01LS L G R LSU J J U R= + − − , (4.22)
onde 0LSR é o valor da fração de líquido no pistão 0, e 0RU é a velocidade relativa
das bolhas dispersas no pistão em relação à velocidade do líquido, definida como:
Capítulo 4 Método de solução numérica 70
14 3
40 022 senR S
L
gU Rσ ρ θρ
⎛ ⎞Δ= ⎜ ⎟
⎝ ⎠. (4.23)
onde σ é a tensão superficial entre o líquido e o gás e θ é a inclinação da
tubulação. Nota-se das Equações (4.22) e (4.23) que os valores de 0LSR , 0GJ e 0LJ
devem ser conhecidos para o cálculo de 0LSU .
No presente trabalho, os valores de 0GJ e 0LJ são conhecidos
experimentalmente na entrada da tubulação. O valor de 0LSR é obtido através do
processo descrito mais adiante.
Comprimentos do pistão e da bolha
No presente trabalho os comprimentos do pistão e da bolha são conhecidos
experimentalmente na entrada da tubulação.
Fração de líquido no pistão e na bolha
Para se calcular os valores das frações de líquido e gás no pistão e na bolha
na entrada da tubulação, considera-se válido o modelo de estado estacionário
apresentado por Taitel e Barnea (1990). O trabalho apresentado pelos autores
baseia-se em duas equações. A primeira delas é a conservação da massa na célula,
que pode ser arranjada da forma:
( )( )0 0 0 0 0 0 0 00
0 0
1G LS L G R LS TGB
T
J R J J U R UR
Uβ
β− − + + −
= , (4.24)
onde a fração de gás na bolha, GBR , foi explicitada, e a velocidade de translação da
bolha, 0TU , é calculada como:
0 0 0 1T LSU C U C gD= + , (4.25)
Capítulo 4 Método de solução numérica 71
como apresentado no capítulo de revisão bibliográfica e desprezando-se o efeito de
esteira.
Nota-se que a Equação (4.24) é função da fração de líquido no pistão, 0LSR , e
os comprimentos do pistão e da bolha estão agrupados através do fator de
intermitência, definido como:
B
B S
LL L
β =+
. (4.26)
É importante ressaltar que a Equação (4.24), como apresentada em trabalhos
anteriores, é utilizada para se calcular o valor do fator de intermitência, β , em
função das frações de líquido no pistão e no filme já conhecidas. Da forma como foi
implementada no presente trabalho, para se calcular a fração de gás na região da
bolha, pode-se obter resultados incoerentes se as estimativas para a fração de
líquido no pistão e do fator de intermitência não estiverem corretas.
A segunda equação apresentada por Taitel e Barnea (1990) é uma equação
diferencial obtida através da conservação de quantidade de movimento no filme de
líquido, mais conhecida como “modelo de bolha”. Essa equação relaciona a altura do
filme de líquido em função do comprimento da bolha, e é apresentada a seguir com
a notação utilizada originalmente pelos autores:
( )
( ) ( ) ( )( )( )22
1 1 sen
1cos
1
f f G Gi i L G
f G f Gf
f ft L S t b SL G L f G G
f f ff
S S S gA A A Adh
dR dRu u R u u Rdz g v vR dh dhR
τ τ τ ρ ρ β
ρ ρ β ρ ρ
⎛ ⎞− − + + −⎜ ⎟⎜ ⎟
⎝ ⎠=− − −
− − −−
, (4.27)
onde fh é a altura do filme de líquido, que é função do comprimento z ao longo da
bolha, como apresentado na Figura 4.3. Os termos no numerador da Equação (4.27)
referem-se às forças de atrito no líquido, atrito no gás, atrito na interface e força
gravitacional. Os termos no denominador são referentes à variação da pressão
hidrostática do líquido e do gás, e à variação da quantidade de movimento do líquido
Capítulo 4 Método de solução numérica 72
e do gás. No trabalho de Taitel e Barnea (1990) são apresentadas as formas para o
cálculo de cada termo da Equação (4.27) para diversas condições de escoamento.
Todos os termos da Equação (4.27) são avaliados localmente como função
do comprimento z . A integração da Equação (4.27) de 0z = até Bz L= resulta no
perfil da bolha e na altura de filme de líquido média. A fração de gás média na bolha
0GBR é calculada em função da altura do filme de líquido média, conhecendo-se a
geometria da bolha na seção transversal, estratificado no escoamento horizontal ou
anular no escoamento vertical.
Para a integração da Equação (4.27) é necessário conhecer o critério de
parada, ou seja, o comprimento da bolha 0BL . Além disso, essa equação ainda é
função da fração de líquido no pistão, 0LSR .
Dessa forma, têm-se duas incógnitas, 0LSR e 0GBR , e duas Equações, (4.24) e
(4.27). Portanto, pode-se solucionar o sistema através de algum método específico.
Figura 4.3 – Representação da geometria utilizada na integração do perfil da bolha
A partir do exposto acima, o processo de cálculo das variáveis na entrada da
tubulação é representado na Figura 4.4. Quando uma nova célula deve ser inserida
na tubulação, a velocidade superficial do líquido, 0LJ , é obtida dos dados
experimentais. Os valores das variáveis 0GJ , 0BL e 0SL podem ser obtidos de duas
formas. A primeira é através de dados médios experimentais, como o valor de 0LJ . A
segunda forma é através de listas com distribuições independentes para cada
variável, como descrito na subseção seguinte.
Capítulo 4 Método de solução numérica 73
Na seqüência, o valor da fração de líquido, 0LSR , deve ser estimado. A partir
dessa estimativa, e do conhecimento das demais variáveis, calcula-se o valor da
fração de gás na bolha, 0GBR , através do balanço de massa e do modelo de bolha,
Equações (4.24) e (4.27), respectivamente. Se o valor de 0GBR for diferente nos dois
cálculos, um novo valor de 0LSR deve ser estimado. Quando os dois valores de 0GBR
forem iguais (dentro de determinada margem de erro), o processo é finalizado e
obtém-se o valor das variáveis necessárias para a nova célula.
Um problema que pode surgir nessa metodologia é que, como exposto
anteriormente, os valores obtidos para 0GBR podem ser incoerentes ( 0 1GBR > ou
0 0GBR < ). Isso ocorre devido à aleatoriedade na escolha dos comprimentos de
bolhas e pistões, e, consequentemente, do fator de intermitência. Por isso, deve-se
ter cuidado nos resultados obtidos para 0GBR , e, em alguns casos, deve-se descartar
o resultado.
Capítulo 4 Método de solução numérica 74
Figura 4.4 – Fluxograma para cálculo dos parâmetros da nova célula que será
inserida na tubulação
Valores aleatórios
Na subseção anterior foi apresentado o processo de cálculo das variáveis na
entrada de uma nova célula. Foi considerado que os valores médios de 0BL , 0SL e
0GJ eram conhecidos experimentalmente. No entanto, de acordo com resultados de
medidas experimentais, como os apresentados na Figura 4.5 para os comprimentos
de bolha e pistão e velocidade de translação da bolha, nota-se que essas variáveis
apresentam uma grande distribuição em torno do valor médio. Essas distribuições
experimentais são caracterizadas por um valor médio e um desvio padrão.
Estima-se um valor de RLS0
JL0 conhecido dos dados de entrada
LB0, LS0 e JG0 escolhidos aleatoriamente de distribuições
Nova célula deve ser inserida na tubulação
Calcula-se RGB0 através do balanço de massa, Eq. (4.24)
Calcula-se RGB0 através do modelo de bolha, Eq. (4.27)
Calcula-se o valor de UT0, Eq. (4.25)
Valores para a nova célula:
LB0, LS0, RGB0 e RLS0
RGB0 do modelo de bolha é igual ao RGB0 do
balanço de massa? SimNão
Capítulo 4 Método de solução numérica 75
0,00
0,01
0,02
0,03
0,04
0 20 40 60 80 100LB/D
0,00
0,03
0,06
0,09
0,12
0 10 20 30 40LS/D
0
1
2
3
4
1,5 2,0 2,5 3,0UT (m/s)
Figura 4.5 – Representação de distribuições experimentais
O que se propôs para esse trabalho, é que os valores de 0BL , 0SL e 0GJ que
são calculados na entrada da simulação também sejam distribuídos, e, ainda, as
distribuições dessas variáveis devem ser independentes. Além disso, espera-se que
tanto o valor médio como os desvios padrão de cada variável na entrada da
simulação sejam semelhantes ao dos experimentos. Baseado na análise de curvas
como a da Figura 4.5 para diversas condições de escoamento, foi definido que os
comprimentos de bolha e velocidade de translação apresentam distribuição normal,
enquanto os comprimentos de pistão têm distribuição log-normal.
Para que isso ocorra, foram geradas listas de valores para 0BL , 0SL e 0GJ , de
forma que os valores médios e o desvio padrão dessas variáveis sejam iguais aos
dos experimentais. Essas listas são independentes entre si, e seguem uma
distribuição normal para 0BL e 0GJ , e log-normal para 0SL . No momento em que uma
nova célula entra na tubulação, um valor de cada lista é sorteado e atribuído à
célula.
Nos dados experimentais obtém-se facilmente os valores médios e desvios
padrão para os comprimentos de bolha e pistão e para a velocidade de translação
das bolhas, mas não se tem a medida das distribuições para a velocidade superficial
do gás, 0GJ . No entanto, essa variável pode ser relacionada com a velocidade de
translação da bolha, 0TU . A seguir, é apresentada a forma como essas variáveis são
relacionadas, e, mais adiante, como as listas de valores são obtidas.
A variável 0GJ pode ser relacionada com a velocidade de translação da bolha
através de:
( )0 0 0 0 1 .T L GU C J J C g D= + + , (4.28)
Capítulo 4 Método de solução numérica 76
Além disso, a velocidade superficial do líquido é considerada constante, devido à
incompressibilidade do líquido. Portanto, conhecendo-se o coeficiente de variação
(desvio padrão dividido pela média) de 0TU , pode-se calcular o coeficiente de
variação de 0GJ por:
0
JG UT T
G T G
s s UJ U C J
= . (4.29)
Conhecendo-se experimentalmente os valores médios e dos coeficientes de
variação de 0BL , 0SL e 0GJ , é necessária a geração das listas de valores distribuídos
dessas variáveis. A seguir, o procedimento para obtenção dessas listas é
apresentado.
Geração das listas aleatórias distribuídas
Para a obtenção das listas de valores distribuídos das variáveis, foi utilizada a
transformação proposta por Box e Muller (1958), na qual uma lista de dados
aleatórios com distribuição normal pode ser gerada a partir de duas listas
independentes de valores aleatórios com distribuição uniforme (entre 0 e 1). O
programa computacional Microsoft Excel® é capaz de gerar essas listas uniformes,
através da função ALEATÓRIO().
Se 1U e 2U são duas listas geradas pela função citada, uma terceira lista é
calculada por:
2 1cos(2 ) 2 ln( )N U Uπ= − , (4.30)
onde N é uma lista de valores aleatórios com distribuição normal padrão (média 0 e
desvio padrão 1).
Em seguida, deve-se utilizar essa lista com distribuição normal padrão e obter
listas com distribuição normal (para 0BL e 0GJ ) e log-normal (para 0SL ) em função
das médias e desvios padrão experimentais.
Capítulo 4 Método de solução numérica 77
Para a geração de listas normais com média μ e desvio padrão s genéricos,
aplica-se para cada valor da lista normal padrão a transformação:
i iN sϕ μ= + , (4.31)
onde iN são os valores da lista normal padrão e ϕi são os valores da lista com
média e desvio padrão desejados.
Na geração de uma lista log-normal de valores com média μ e desvio padrão
s genéricos, utiliza-se a transformação:
( )0 50%exp lni iNψ σ μ= +⎡ ⎤⎣ ⎦ , (4.32)
onde:
2 2
00 50%ln 1 , exp
2s x σσ μμ
⎡ ⎤ ⎛ ⎞⎛ ⎞= + = −⎢ ⎥ ⎜ ⎟⎜ ⎟
⎝ ⎠⎢ ⎥ ⎝ ⎠⎣ ⎦. (4.33)
e ψi são os valores da lista com distribuição log-normal com média e desvio padrão
desejados.
Após a obtenção das 3 listas independentes para 0BL , 0SL e 0GJ , a cada nova
célula que entra na tubulação, um valor de cada variável é sorteado dessas listas e
utilizado como dado de entrada da simulação.
Formas alternativas para o cálculo dos parâmetros da nova célula
Quando não se tem os dados experimentais na entrada da tubulação, outras
equações devem ser utilizadas para calcular as variáveis 0BL , 0SL , 0GBR e 0LSR . A
seguir algumas dessas equações são apresentadas.
De forma geral, o processo de definição das variáveis ainda segue o
fluxograma descrito na Figura 4.4. No entanto, devem-se obter os valores de 0BL ,
0SL e 0GJ de alguma outra forma. Além disso, uma correlação pode ser utilizada
Capítulo 4 Método de solução numérica 78
para se calcular o valor de 0LSR , ou até, no caso mais simples para escoamento
horizontal, se definir o pistão como não aerado, e RLS0=1.
O valor da velocidade superficial do gás, 0GJ , está relacionado à vazão
volumétrica de gás, que é um dado conhecido. No entanto, é necessário o valor de
0GJ na entrada, e, experimentalmente pode-se conhecer a vazão na entrada ou na
saída da tubulação. Se a vazão for conhecida na saída, deve-se corrigir a velocidade
superficial de gás da saída para a entrada do tubo pela variação da pressão, através
de:
( ) ( ) ( )0 0
( 0)G G G
P z LJ J z J z L
P z=
= = = ==
. (4.34)
As próximas variáveis que devem ser calculadas são os comprimentos da
bolha e do pistão. Uma forma de se obter essas variáveis é através do conhecimento
da freqüência da célula unitária, f, e do fator de intermitência, β. Existem diversas
correlações para o cálculo da freqüência, como apresentado na seção de revisão
bibliográfica. Quanto ao fator de intermitência, trabalhos recentes como o de Rosa
(2006) tem mostrado que existe correlação entre β e a razão entre JG e J.
Sabendo-se os valores de f e β, as seguintes equações cinemáticas podem
ser utilizadas:
( )1, ,T T
U B S BU UL L L L
f fββ
β−
= = = , (4.35)
para o cálculo de LB e LS. Essas equações também são válidas para o modelo de
estado estacionário.
Por fim, devem-se obter os valores das frações de líquido no pistão e gás na
bolha, 0LSR e 0GBR . Pode-se novamente utilizar o processo apresentado na Figura
4.4, com a estimativa de 0LSR e a comparação dos dois valores de 0GBR calculados.
Por outro lado, pode-se utilizar alguma correlação para o cálculo de 0LSR , como
apresentado no capítulo de revisão bibliográfica, ou considerar-se RLS0=1. Dessa
Capítulo 4 Método de solução numérica 79
forma, utiliza-se apenas a equação de conservação da massa para o cálculo de
0GBR .
Existem diversas outras formas de cálculo das variáveis 0BL , 0SL , 0GBR e 0LSR ,
a depender de quais variáveis são conhecidas ou estimadas através de correlações.
É necessária sempre a combinação das variáveis conhecidas com as equações
cinemáticas (4.37) e equações do balanço de massa e modelo de bolha.
Como mais um exemplo, pode-se utilizar uma correlação para 0LSR e se obter
experimentalmente o valor de 0SL . Nesse caso, deve-se iniciar a integração
numérica do modelo de bolha, e a cada incremento da variável z, calcula-se o 0GBR
médio do modelo de bolha. Calcula-se também o valor de 0GBR médio na equação de
conservação da massa, utilizando-se LB=z. Quando os dois valores de 0GBR forem
iguais, o processo finaliza, com o valor de 0GBR obtido e o comprimento da bolha é
igual ao último valor de z incrementado.
No presente trabalho, verificou-se que os valores de 0LSR obtidos das
correlações não são confiáveis. Por isso optou-se pelo cálculo de 0LSR e 0GBR
através das equações de conservação da massa e modelos de bolha. Para a
definição dos comprimentos de bolha e pistão, foi escolhido utilizar os valores
experimentais, após um tratamento para gerar distribuições aleatórias.
4.5 Processo de saída de bolhas e pistões do domínio de cálculo
A Figura 4.6 apresenta o processo de retirada de células do tubo. No tempo
nt , a bolha n está com sua frente, ny , exatamente na posição z L= . Nesse
momento, é considerado que a pressão nessa bolha é igual à pressão de saída,
geralmente a atmosférica, e n atmP P= . A partir desse momento, o sistema de
equações acopladas é resolvido apenas da célula 1 até a 1n − , e não até a célula n .
Isso ocorre, pois não faz sentido calcular a variação da quantidade de movimento do
pistão n , que não pertence mais ao domínio de cálculo. Além disso, a pressão no
interior da bolha deixa de ser incógnita e passa a ser conhecida e constante.
Capítulo 4 Método de solução numérica 80
Como não são calculados os parâmetros da célula n (comprimentos,
velocidades, frações de líquido e gás), também não é calculada a velocidade da
frente da bolha n , ny , e a velocidade da traseira da bolha n , 1nx − . Porém, para que a
bolha n efetivamente deixe o tubo, é necessária a definição da velocidade da
fronteira 1nx − . Para que a bolha n não exerça nenhuma influência no escoamento a
montante, considera-se que a velocidade da sua traseira é igual à velocidade da
frente da bolha 1n − , ou 1 1n nx y− −= , o que implica que o comprimento do pistão 1n −
permanece constante durante a saída da bolha n .
Esse processo é representado na Figura 4.6 nos tempos 1nt + até o tempo 2nt + ,
quando a frente do pistão 1n − chega a z L= . Nesse instante, a última célula é
eliminada do domínio de cálculo, e a célula 1n − passa a ser a célula n , como
representado em '2nt + .
A partir do instante '2nt + até o instante 3nt + , um pistão está deixando o tubo.
Durante esse tempo, o sistema de equações acopladas deve ser resolvido desde a
célula 1 até a célula n . Porém, o modelo numérico necessita que uma bolha sempre
esteja na saída do tubo, para que a condição de contorno ( 1GBnP + ) seja aplicada.
Uma das possibilidades é considerar que a frente do pistão n permanece
estagnada em z L= e há uma ‘bolha’ em z L≥ garantindo a pressão constante.
Nesse caso, o comprimento do pistão n diminui com o tempo até chegar a 0SnL = ,
para que uma nova bolha comece a sair. No entanto, esse procedimento pode levar
à desestabilização do sistema de equações, pois o último termo da diagonal
principal da matriz de coeficientes, que é proporcional ao comprimento do pistão,
apresenta um valor tendendo a zero porque o comprimento do pistão diminui
continuamente.
Em vista deste potencial problema numérico buscou-se uma segunda
alternativa de aplicação da condição de contorno sem que o comprimento do pistão
varie. Considerando-se que o comprimento do pistão n , SnL , permaneça constante
também quando o pistão está saindo do tubo, o comprimento do tubo, L , é
aumentado artificialmente até a frente do pistão ( 'n SnL y L= + ).
Capítulo 4 Método de solução numérica 81
Porém, a condição de contorno para a pressão é rigorosamente definida em
z L= como sendo GB atmP P= e, portanto será necessário conhecer o valor da pressão
'atmP em 'z L= que irá produzir em z L= uma pressão GB atmP P= . O valor de '
atmP deve
ser definido em função de atmP , e das parcelas referentes à queda de pressão devido
ao atrito e ao peso do pistão na região z L> , ou:
( )' 22 sen ,natm atm n L n L Sn n n Sn
CfP P x L U gR x y LD
ρ ρ θ⎛ ⎞= − − + = +⎜ ⎟⎝ ⎠
. (4.36)
De acordo com a Equação (4.36), 'atm atmP P≤ . O valor de '
atmP é artificialmente criado a
fim de que em z L= a pressão seja atmP . Quando a frente da bolha n chega a z L= ,
volta-se ao processo iniciado em nt .
Condições a serem respeitadas na saída do domínio (pseudo-algoritmo):
• Enquanto uma bolha sai do tubo ( ny L> e 1nx L− < )
o O sistema é resolvido até a célula 1n − com a condição de contorno de
pressão na bolha n ( GBn atmP P= )
o 1SnL − = cte
• Enquanto um pistão está saindo do tubo ( nx L> e ny L< )
o O sistema é resolvido até a célula n
o SnL = cte
o Condição de contorno para a pressão na célula 1n + é 'atmP , definida pela
Equação (4.36)
• No instante em que o pistão 1n − chega em z L= ( 1Nnx L− > e 1
Onx L− < )
o Célula n é eliminada do domínio de cálculo
o Célula 1n − passa a ser n
Capítulo 4 Método de solução numérica 82
0
t
zz=L
tn
tn+1
tn+2
tn+2
tn+3
PGBnPGBn-1PGBn-2
PGBn-2
PGBn-2
PGBn-1
PGBn-1
ULSn-1ULSn-2
yn-2
yn-2
yn-2
yn-1
yn-1
xn-3
xn-3
xn-3
xn-2
xn-2
xn-2
xn-2
xn-2
xn-1
xn-1
yn-1 yn
yn
LSn-1
LSn-1
LSn-1
LSn
LSn
yn-1
yn-1
yn
yn
xn-1
xn-1
xn-1
xn
xn
ULSn-2
ULSn-2
ULSn-1
ULSn-1
ULSn-1
ULSn-1
ULSn
ULSn
PGBn-1
PGBn-1
PGBn
PGBn
PGBn
PGBn
Figura 4.6 – Representação esquemática do processo de saída de bolhas e pistões
do domínio computacional
4.6 Processo de início da simulação
Como condição inicial do modelo deve-se definir os valores das propriedades
de todas as células (comprimentos, velocidades, pressões) no interior do tubo. Para
facilitar o processo, considerou-se que no início da simulação o tubo está cheio de
líquido, e a primeira bolha está prestes a entrar no tubo ( 1 0y = ), como definido na
Figura 4.7.
Dessa forma, devem-se conhecer apenas as propriedades da célula 1. A
primeira bolha do tubo comporta-se da mesma forma que uma bolha que está
Capítulo 4 Método de solução numérica 83
iniciando sua entrada no tubo, como descrito na seção 4.4, e as propriedades da
célula 0 devem ser definidas nesse instante de tempo.
Porém, além disso, deve-se definir a velocidade do pistão na célula 1. A
velocidade no pistão 1 é uma função das vazões de gás e líquido na entrada do
tubo. Se a vazão de gás local, na entrada do tubo é conhecida, o problema está
resolvido. Porém, se a vazão de gás na saída é conhecida, a vazão na entrada é
calculada a partir das pressões na saída e na entrada. A pressão na entrada do tubo
pode ser calculada através da queda de pressão por atrito e gravitacional no
escoamento monofásico ao longo do tubo.
No entanto, a queda de pressão por atrito depende da velocidade do
escoamento, que não é conhecida (pois depende da vazão de gás). Para resolver
esse problema, utiliza-se um processo iterativo no qual, através de uma estimativa
inicial para a vazão de gás, calcula-se a velocidade e a queda de pressão ao longo
do tubo. Com a queda de pressão calculada, pode-se corrigir a vazão de gás pela
pressão na entrada. Esse processo é realizado até a convergência da velocidade e
pressão na entrada.
0
t
zz=0 z=L
t=0 PGB1PGB0ULS0
y0 x0LS0LB0 y1 x1
ULS1 ULS1
Figura 4.7 – Representação esquemática da condição inicial
4.7 Algoritmo de simulação
A Figura 4.8 apresenta o algoritmo geral de simulação. No início do processo,
o tubo está cheio de líquido, e a primeira bolha está prestes a entrar no tubo. Nesse
momento, é necessário calcular o comprimento dessa bolha e as frações de líquido
no pistão (que tem comprimento igual ao do tubo) e na bolha. Essa célula é
numerada como 1.
Capítulo 4 Método de solução numérica 84
A cada nova iteração do programa, são verificadas a entrada, saída e
coalescência de bolhas. Se alguma dessas singularidades ocorrer, sub-rotinas
específicas são chamadas. A seguir, as variáveis do tempo antigo são igualadas às
variáveis do tempo novo, e todos os termos necessários à resolução do sistema
acoplado são calculados. Em seguida, o sistema é resolvido, obtendo-se os valores
das velocidades do líquido nos pistões e as pressões do gás nas bolhas no tempo
novo.
As equações auxiliares são resolvidas na seqüência, e os deslocamentos dos
pistões e das bolhas ao longo de um passo de tempo são calculados. Após a
movimentação das células, são armazenados os valores das variáveis nas sondas
virtuais implementadas no modelo, que serão apresentadas na seção seguinte.
Após o movimento das células, é necessário verificar se a primeira bolha
entrou totalmente no tubo. Se isso ocorrer, os parâmetros relativos à nova célula
(célula 0) devem ser calculados.
Por fim, o critério de parada é analisado. Nesse trabalho, utiliza-se o número
de bolhas que saem do tubo como critério, e se ele não é satisfeito, o tempo é
incrementado de um passo de tempo e uma nova iteração inicia. Caso contrário,
todos os valores armazenados nas sondas virtuais são salvos em arquivos de saída
dos resultados e o programa é finalizado.
Capítulo 4 Método de solução numérica 85
Entrada de uma bolha?L >LS0 SIni
Saída de uma bolha?x >Ln-1
Coalescência de bolhas?y>x (j=1, 2,...,n-1)j j
Atualização das variáveisO N
n n-1
t t+ tΔ Renumeraçãodas células
Renumeraçãodas células
Cálculo daspropriedades para
a célula 0
Cálculo daspropriedades para
a célula 1
Salvaresultados FIM
INÍCIO
Cálculo dos termos da matriz
Resolução das equações acopladas(Utiliza o )TDMA
Resolução das equações auxiliares(Deslocamento das frente de bolhas e pistões)
Passagem pelas sondas virtuais(Armazenamento dos dados)
Entrada da célula 0?x >00
Número de bolhas que saíram do tuboé igual ao especificado?
Sim
Sim
Sim
Sim
Sim
Não
Não
Não
Não
Não
Figura 4.8 – Algoritmo geral de simulação
Capítulo 4 Método de solução numérica 86
4.8 Sondas virtuais
Nos modelos de seguimento de pistões é possível implementar diversas
sondas virtuais para monitorar a solução do escoamento em golfadas. No presente
trabalho foram utilizadas as seguintes sondas virtuais:
• Sondas eulerianas: Ao longo do tubo são implementadas 8 sondas eulerianas,
que armazenam o valor de todos os parâmetros de todas as células que
cruzam determinada seção da tubulação. Através desses resultados, obtêm-
se os valores médios de todas as variáveis ao longo do tempo em cada
posição do tubo, bem como as funções densidade de probabilidade (PDF,
Probability Density Functions), que mostram as distribuições das variáveis;
• Sonda lagrangeana (seguindo uma célula): É possível acompanhar uma
determinada célula ao longo de toda a sua passagem pela tubulação,
armazenando-se o valor de todos os parâmetros a cada passo de tempo.
Essa sonda é interessante para avaliar as oscilações típicas do escoamento
em golfadas;
• Sondas eulerianas na entrada e na saída: Duas sondas eulerianas
posicionadas na entrada e na saída do tubo são implementadas. Essas
sondas, além de armazenar as variáveis importantes do escoamento,
calculam e armazenam os valores dos fluxos de massa de líquido que estão
entrando e saindo no tubo. Esses resultados são utilizados em conjunto com
os resultados obtidos pela sonda virtual de massa global, que será
apresentada na seqüência. Além disso, a sonda euleriana na entrada do tubo
mostra as variações de pressão que ocorrem no início da simulação, quando
o tubo ainda está sendo preenchido com bolhas;
• Sonda de massa global: A cada passo de tempo, essa sonda armazena o valor
da massa global de líquido que está presente no tubo. Através de um balanço
de massa de líquido no tubo inteiro, a variação dos valores obtidos através
dessa sonda deve corresponder à diferença de fluxos de massa que entram e
saem do tubo, obtidos nas sondas eulerianas na entrada e saída do tubo;
• Sonda de fotos: É possível armazenar os valores das variáveis para todas as
células no tubo em determinado instante de tempo, como em uma foto. Essa
Capítulo 4 Método de solução numérica 87
sonda é capaz de armazenar várias fotos, em diferentes instantes de tempo, e
é interessante para avaliar a evolução de diversas células ao longo do tubo
em alguns instantes de tempo.
4.9 Simulações preliminares
Nessa seção, os resultados de algumas simulações preliminares são
apresentados. Esses resultados têm como objetivo a definição de parâmetros
importantes para as simulações, como o passo de tempo, as constantes utilizadas
no modelo para o fator de esteira, e a utilização da condição de contorno
intermitente em detrimento da condição de contorno periódica.
As simulações preliminares foram realizadas utilizando-se as configurações
apresentadas na Tabela 4.1. É importante ressaltar que as simulações preliminares
foram realizadas com a condição de entrada de bolhas e pistões aleatória, a não ser
na seção 4.9.3 onde as condições de entrada são analisadas.
Tabela 4.1 – Definição dos parâmetros utilizados nas simulações preliminares Ar e água (A@W) Ar e glicerina (A@G)
Comprimento do tubo (m) 20,098 20,098
Diâmetro do tubo (m) 0,026 0,026
Inclinação do tubo (º) 0 0
Viscosidade do líquido (Pa.s) 0,000855 0,0238
Densidade do líquido (kg/m3) 1000 1190
Velocidade superficial do líquido (m/s) 0,5 0,5
Velocidade superficial do gás na saída (m/s) 0,5 0,5
4.9.1 Análise do passo de tempo
O passo de tempo utilizado na simulação numérica influencia o tempo de
simulação computacional e a qualidade da simulação. Para definição do passo de
tempo a ser utilizado, foram analisados passos de tempo variando entre 0,01s e
0,00001s para o escoamento de ar e água, e armazenados os tempos de simulação
em um PC Intel® Pentium® 4 CPU 3.60GHz com 1.00GB de RAM. Os resultados são
apresentados na Figura 4.9.
Capítulo 4 Método de solução numérica 88
Para analisar a qualidade da simulação, foi definido um resíduo relativo de
massa de líquido na tubulação, que é apresentado a seguir. A equação de
conservação da massa de líquido aplicada a toda a tubulação é:
Sai Entra 0LL L
dM m mdt
+ − = , (4.37)
onde LM é a massa de líquido na tubulação, SaiLm é o fluxo mássico de líquido que
sai da tubulação em z L= e EntraLm é o fluxo mássico de líquido que entra na
tubulação, em 0z = . Discretizando-se a Equação (4.37) através do método das
diferenças finitas e utilizando-se o esquema semi-implícito de Crank-Nicholson,
obtém-se:
Sai Sai Entra Entra
2 2
N O N ON OL L L LL L
Lm m m mM M
tε+ +−
+ − =Δ
, (4.38)
onde o sobrescrito N representa as variáveis em um tempo novo, e O representa as
variáveis em um tempo antigo. Note que a o lado direito da Equação (4.38) não é
nulo, devido a erros numéricos na simulação. Normalmente, quanto menor o passo
de tempo, tΔ , os erros numéricos são menores. Na Equação (4.38) chamou-se esse
erro de resíduo da massa de líquido, Lε . Esse resíduo é absoluto, e tem unidade de
kg s . No entanto, utilizou-se como medida para analisar a qualidade da simulação o
resíduo relativo da massa de líquido, definido por:
RelRef
LL
L
εεε
= , (4.39)
onde RefLε é um resíduo de massa de líquido de referência. O resíduo de referência
deve ser calculado em função de parâmetros conhecidos do escoamento, e decidiu-
se pela utilização de:
RefL L LJ Aε ρ= , (4.40)
Capítulo 4 Método de solução numérica 89
que representa o fluxo médio de massa de líquido na tubulação, e LJ é a velocidade
superficial do líquido ao longo da tubulação.
Portanto, dividindo-se a Equação (4.38) pela Equação (4.40) e reorganizando,
obtém-se:
( ) Sai Sai Sai Sai Entra Entra Entra Entra
Rel 2
N O N N O O N N O OL L L L L L L L L L
LL L
L R R R U R U R U R UJ t J
ε− + − −
= +Δ
, (4.41)
que é o resíduo relativo da massa de líquido na tubulação, onde L é o comprimento
da tubulação, LR é a fração de líquido ao longo de todo o tubo, definida pela razão
do volume de líquido e o volume da tubulação, SaiLR e EntraLR são as frações de
líquido nas áreas de seção transversal na saída e entrada da tubulação, e SaiLU e
EntraLU são as velocidades do líquido na saída e entrada da tubulação.
A cada passo de tempo, o resíduo relativo é calculado, e ao final da
simulação, o seu valor RMS (root mean square) é calculado através de:
2RMS
1
1 N
L LiiN
ε ε=
= ∑ , (4.42)
onde N é o número de passos de tempo utilizado. O valor RMS é utilizado, pois o
valor médio é muito pequeno, devido aos pequenos erros se anularem ao longo da
simulação.
A Figura 4.9 apresenta os resultados obtidos para o resíduo relativo para os
passos de tempo utilizados. Nota-se que quanto menor o passo de tempo, menor é
o resíduo relativo, porém, maior é o tempo necessário para a simulação. Com base
nos resultados apresentados, escolheu-se, para o restante das simulações, o passo
de tempo de 0,0001tΔ = , pois obtém-se um resíduo relativo da ordem de 10-3, com
tempos de simulação na ordem de minutos.
A grosso modo, pode-se dizer que um resíduo relativo da ordem de 10-3
implica que a cada 1 kg s de fluxo de líquido que entra na tubulação, ocorre um erro
numérico de 0,001 kg s . É importante ressaltar que os resultados apresentados
Capítulo 4 Método de solução numérica 90
foram obtidos após 600 bolhas saírem do tubo, o que garante uma população de
bolhas e pistões grande o suficiente para que os valores médios obtidos sejam
representativos.
1,E+01
1,E+02
1,E+03
1,E+04
1,E+05
1,E-05 1,E-04 1,E-03 1,E-02Passo de tempo, s
Tem
po d
e si
mul
ação
(CPU
), s
1,E-04
1,E-03
1,E-02
1,E-01
1,E-05 1,E-04 1,E-03 1,E-02Passo de tempo, s
RM
S do
resí
duo
rela
tivo
Figura 4.9 – Variação do tempo de simulação e do resíduo relativo em função do
passo de tempo utilizado.
4.9.2 Análise das constantes do fator de esteira
Um fator importante para os resultados das simulações são as constantes
utilizadas na equação do fator de esteira, apresentado no capítulo de revisão
bibliográfica e repetido a seguir por conveniência:
exp Sw w
Lh a bD
⎛ ⎞= −⎜ ⎟⎝ ⎠
. (4.43)
A constante wa , que multiplica a função exponencial, está relacionada com a
velocidade máxima que uma bolha tem quando se aproxima da bolha à sua frente.
Por outro lado, a constante wb está relacionada ao tamanho do pistão em que o
efeito de esteira começa a ser importante.
Diversos valores para as constantes wa e wb existem na literatura, porém,
elas foram obtidas em diferentes condições experimentais. Foram testados
diferentes valores para as constantes de esteira no escoamento de ar e água e ar e
glicerina. A Figura 4.10 apresenta os resultados médios ao longo do tubo para os
comprimentos de bolhas e pistões, taxa de coalescência e velocidade da bolha
alongada. Nota-se que a variação das constantes de esteira influenciam
Capítulo 4 Método de solução numérica 91
significativamente na coalescência de bolhas, e, conseqüentemente, no aumento
dos comprimentos de bolhas e pistões.
Na ausência do fator de esteira, os comprimentos praticamente não variam.
Por outro lado, quando o fator de esteira proposto por Grenier (1997) é utilizado,
ocorrem mais coalescências na simulação do que nos resultados experimentais. Por
esse motivo, valores intermediários para as constantes foram testados. Baseado nos
resultados apresentadas na Figura 4.10 adotou-se os valores de 0,4wa = e 1,0wb =
para todas as simulações de ar e água.
Ainda na Figura 4.10 pode-se notar alguns detalhes. Os comprimentos de
bolha na entrada do tubo simulados e experimental são um pouco diferentes. Isso
ocorre porque foi utilizado como condição de entrada a freqüência do escoamento, e
não diretamente o comprimento da bolha. Além disso, os dados para a velocidade
da bolha também apresentam uma grande variação nos dados experimentais, o que
não ocorre na simulação. Nesse caso, acredita-se que o escoamento próximo da
entrada do tubo, nos experimentos, ainda está em desenvolvimento para o padrão
de golfadas. Dessa forma, a velocidade das bolhas é muito alta, pois o processo de
coalescência é intenso. Outro detalhe é que no escoamento em desenvolvimento a
velocidade da bolha não pode ser correlacionada com a velocidade do líquido, como
é feito para o escoamento em golfadas. Como a simulação considera que o
escoamento é em golfadas, os valores de velocidade das bolhas calculados é
diferente dos experimentais.
Para o escoamento de ar e glicerina, espera-se que os valores para o fator de
esteira sejam diferentes, pois as propriedades físicas (densidade e viscosidade) da
glicerina são diferentes da água. A mesma análise apresentada anteriormente foi
realizada para o escoamento de ar e glicerina, e os resultados para os valores
médios dos comprimentos de bolhas e pistões, taxas de coalescências e velocidade
das bolhas são apresentados na Figura 4.11. Nota-se que, nos resultados
experimentais, ocorrem poucas coalescências. Por isso, nas simulações para o
escoamento de ar e glicerina, decidiu-se não utilizar o fator de esteira. Nota-se
também uma diferença na velocidade de bolha entre as simulações e o
experimental. Isso ocorre devido ao valor da constante C0 utilizado na simulação,
que, nesse caso, não se ajustou ao valor experimental.
Capítulo 4 Método de solução numérica 92
10
15
20
25
30
0 200 400 600
L/D
LB/D
0
5
10
15
20
0 200 400 600
L/D
LS/D
0.000
0.001
0.002
0.003
0.004
0.005
0.006
0 200 400 600
L/D
R
aw=0.4, bw= 0.5 (Grenier)aw=0.4, bw= 0.8aw=0.4, bw=1.0aw=0.4, bw=1.2aw=0.3, bw=0.5aw=0.2, bw=0.5aw=0.1, bw=0.5Sem este iraExperim ental
1.00
1.05
1.10
1.15
1.20
1.25
1.30
1.35
1.40
1.45
1.50
0 200 400 600
L/DU
T
Figura 4.10 – Evolução dos comprimentos de bolhas e pistões, taxa de coalescências e velocidade da bolha ao longo do tubo
para diversos valores das constantes de esteira no escoamento de ar e água
Capítulo 4 Método de solução numérica 93
9
12
15
18
21
0 200 400 600
L/D
LB/D
5
10
15
20
0 200 400 600
L/D
LS/D
0.000
0.001
0.001
0.002
0.002
0.003
0.003
0 200 400 600
L/D
R
aw=0.4, bw= 0.5 (Grenier)aw=0.4, bw= 0.8aw=0.4, bw=1.0aw=0.4, bw=1.2aw=0.3, bw=0.5aw=0.2, bw=0.5aw=0.1, bw=0.5Sem este iraExperim ental
1.15
1.20
1.25
1.30
1.35
0 200 400 600
L/D
UT
Figura 4.11 – Evolução dos comprimentos de bolhas e pistões, taxa de coalescências e velocidade da bolha ao longo do tubo
para diversos valores das constantes de esteira no escoamento de ar e glicerina
Capítulo 4 Método de solução numérica 94
4.9.3 Análise da condição de contorno aleatória versus condição periódica
Nesse trabalho duas condições de entrada de células foram consideradas. Na
condição periódica, todas as células que entram no tubo possuem os mesmos
comprimentos de bolhas e pistões, além das mesmas frações de líquido no pistão e
no filme. Na condição de contorno intermitente cada célula que entra no tubo possui
parâmetros diferentes, sorteados através de distribuições independentes para os
comprimentos de bolhas e pistões e velocidades superficiais de gás. A seguir são
apresentadas as comparações de resultados utilizando-se as duas condições de
contorno.
Na Figura 4.12 são apresentadas as variações dos comprimentos de bolha e
pistão, velocidade da bolha e pressão para uma única célula seguida ao longo de
sua passagem ao longo do tubo, nas simulações com as duas condições de entrada.
Nota-se que os comprimentos de bolha e pistão apresentam uma variação bem
maior na simulação com condição de contorno intermitente, principalmente devido à
coalescências da bolha que foi seguida.
A coalescência é percebida no aumento instantâneo do comprimento da bolha
em determinados pontos do tubo. O comprimento do pistão também aumenta, pois
na coalescência, o líquido que estava no pistão entre as bolhas que coalesceram é
redistribuído para os pistões adjacentes, que apresentam um aumento acelerado,
porém não instantâneo.
Através desses resultados mostra-se que a condição de contorno intermitente
introduz uma instabilidade no escoamento e facilita a ocorrência de coalescências.
Além disso, são apresentadas as variação da velocidade da bolha e pressão.
Nota-se que, na média, os resultados são semelhantes. No entanto, na simulação
com a condição de contorno intermitente, essas variáveis apresentam oscilações
grandes. Esse efeito também é resultante da variação na condição de entrada.
De forma geral, pode-se dizer que a condição de contorno intermitente facilita
a ocorrência de coalescências, e, além disso, distribui melhor as coalescências ao
longo do tubo. Isso pode ser mais bem visualizado através da Figura 4.13, que
apresenta a taxa de coalescência ao longo do tubo para as duas condições de
contorno.
Capítulo 4 Método de solução numérica 95
Na condição periódica, como todas as células que entram no tubo são iguais,
as coalescências ocorrem todas em uma mesma posição ao longo do tubo. Por isso
a taxa de coalescência é alta apenas nessa região. Nos resultados para a condição
intermitente, a taxa de coalescência é maior na entrada do tubo, quando existe uma
grande variação nos comprimentos de bolha e pistão, e vai diminuindo ao longo do
tubo, conforme as coalescências vão ocorrendo e as bolhas e pistões começam a ter
comprimentos estáveis. Esse comportamento é coerente com o verificado
experimentalmente.
A Figura 4.14 apresenta as funções densidade de probabilidade para os
comprimentos de bolha utilizando-se as duas condições de contorno. Nota-se que,
apesar dos valores médios serem semelhantes, as distribuições apresentam
grandes diferenças. Isso ocorre também com os comprimentos de pistão
apresentados na Figura 4.14 e com as velocidades de translação da bolha,
apresentadas na Figura 4.14. Esses resultados mostram que a utilização da
condição intermitente é vantajosa em relação à periódica, pois a simulação captura
os efeitos aleatórios do escoamento, o que aumenta a interação entre bolhas.
Pelo exposto acima, no restante desse trabalho, as simulações realizadas
para diversas configurações geométricas e de propriedades dos fluidos, somente a
condição de contorno intermitente foi considerada.
Capítulo 4 Método de solução numérica 96
0
10
20
30
40
50
60
0 200 400 600L/ D
LB/D Periód ico
Intermitente
0
5
10
15
20
25
30
0 200 400 600L/ D
LS/D Periód ico
Intermitente
0,7
0,9
1,1
1,3
1,5
1,7
0 200 400 600L/ D
UT
[m/s
]
Periód icoIntermitente
97
99
101
103
105
107
0 200 400 600L/ D
P [k
Pa]
Periód icoIntermitente
Figura 4.12 – Resultados para uma célula seguida ao longo de sua passagem pelo
tubo com as condições de contorno periódica e intermitente
0,00000,00050,00100,00150,00200,00250,00300,00350,00400,00450,0050
0 200 400 600L/ D
Taxa
de
coal
escê
ncia
ExperimentalIntermitentePeriód ico
Figura 4.13 – Taxa de coalescência ao longo do tubo para as condições de contorno
periódica e intermitente
Capítulo 4 Método de solução numérica 97
0 40 80LB/ D
0.00
0.05
0.10
0.15
0.20PD
F
0 20 40LS/ D
0.00
0.05
0.10
0.15
0.20
0.25
0 0.5 1 1.5UT (m/ s)
05
101520253035
ExperimentalIntermitentePeriodico
Figura 4.14 – Funções densidade de probabilidade para os comprimentos de bolhas
e pistões e para a velocidade da bolha experimentais e simulados utilizando-se as
condições de entrada periódica e aleatória
Capítulo 5 Resultados para o escoamento horizontal 98
5 RESULTADOS PARA O ESCOAMENTO HORIZONTAL
Nesse capítulo são apresentados os resultados obtidos na simulação do
escoamento em golfadas na horizontal. As simulações foram realizadas utilizando-se
as configurações geométricas e os fluidos semelhantes aos utilizados nos
experimentos realizados pelo 2PFG (Two-phase flow group) da FEM/UNICAMP.
Essas configurações são descritas na seção 5.1. Na seção 5.2 são apresentadas
comparações entre os valores médios e as funções densidade de probabilidade ao
longo do tubo dos resultados obtidos na simulação e nos experimentos realizados
pelo 2PFG/FEM/UNICAMP e descritos em Rosa (2006).
5.1 Configurações experimentais
As configurações da simulação foram semelhantes às utilizadas nos
experimentos de Rosa (2006) para o escoamento de ar e água, ar e glicerina, e
nitrogênio e óleo SAE 20-50, que estão representadas na Tabela 5.1. Foram
utilizados 6 casos para cada par de fluidos nos escoamento de ar e água e ar e
glicerina, e 4 casos para o escoamento de nitrogênio e óleo SAE 20-50, com
velocidades superficiais de líquido e gás diferentes. A sigla A@W (Air at water) é
usada para identificar o escoamento de ar e água, enquanto a sigla A@G (Air at
gliceryn) é utilizada para identificar o escoamento de ar e solução de glicerina e
N@O (Nitrogen at oil) é utilizada para o escoamento de nitrogênio e óleo.
Tabela 5.1 – Definição das configurações experimentais Ar e água (A@W) Ar e glicerina (A@G) N2 e Óleo SAE 20-50
Comprimento do tubo (m) 20,098 (773 D) 20,098 (773 D) 26,26 (911,8 D) Diâmetro do tubo (m) 0,026 0,026 0,0288 Inclinação do tubo (º) 0 0 0
Viscosidade do líquido (Pa.s) 0,000855 0,0238 0,230 – 0,360 Densidade do líquido (kg/m3) 1000 1190 880 Estação de medição #1 (m) ~ 0 ~ 0 ~ 0 Estação de medição #2 (m) 3,588 (138 D) 3,588 (138 D) 6,595 (229 D) Estação de medição #3 (m) 9,568 (368 D) 9,568 (368 D) 14,77 (513 D) Estação de medição #4 (m) 16,9 (650 D) 16,9 (650 D) -
As condições de escoamento – velocidades superficiais de líquido e gás na
saída, e freqüência de entrada das bolhas na tubulação – utilizadas nos escoamento
Capítulo 5 Resultados para o escoamento horizontal 99
de ar e água, ar e glicerina e nitrogênio e óleo estão apresentados nas Tabelas 5.2,
5.3 e 5.4, respectivamente.
Tabela 5.2 – Definição das condições de contorno para o escoamento de ar e água JL (m/s) JG (m/s) f (Hz)
A@W#1 0,33 0,64 0,740 A@W#2 0,33 1,31 0,814 A@W#3 0,33 1,62 0,705 A@W#4 0,52 0,52 2,890 A@W#5 0,67 0,67 4,470 A@W#6 0,66 1,30 2,440
Tabela 5.3 – Definição das condições de contorno para o escoamento de ar e
glicerina JL (m/s) JG (m/s) f (Hz)
A@G#1 0,33 0,67 0,870 A@G#2 0,33 1,33 0,780 A@G#3 0,33 1,67 0,914 A@G#4 0,50 0,51 1,760 A@G#5 0,67 0,67 3,850 A@G#6 0,67 1,35 2,720
Tabela 5.4 – Definição das condições de contorno para o escoamento de nitrogênio
e óleo JL (m/s) JG (m/s) f (Hz)
N@O#1 0,17 0,45 2,36 N@O#4 0,34 0,32 5,07 N@O#6 0,34 0,86 2,50 N@O#9 0,68 0,36 7,51
5.2 Comparação entre as simulações e os resultados experimentais
5.2.1 Resultados para o escoamento de ar e água
Nessa seção são apresentados os resultados para o escoamento de ar e
água na horizontal para as condições A@W#2 e A@W#4, e os resultados para as
demais condições são apresentados no Apêndice A. As seções de teste
experimental e numérica para a configuração de ar e água na horizontal estão
mostradas na Figura 5.1.
Capítulo 5 Resultados para o escoamento horizontal 100
Figura 5.1 – Representação das seções de teste experimental e numérica para o
escoamento de ar e água na horizontal
Valores médios ao longo do tubo
As Figuras 5.2 e 5.3 apresentam os valores médios ao longo da tubulação da
simulação e dos experimentos. São apresentadas as variações dos comprimentos
de bolha e pistão, velocidade de translação das bolhas e pressão. A variação dos
comprimentos de bolhas e pistões está ligada a dois efeitos: expansão do gás e
interação (coalescência) de bolhas. O efeito de expansão do gás está relacionado à
queda de pressão ao longo da linha, enquanto a interação está relacionada à
intermitência intrínseca do escoamento. Quando somente o efeito de expansão está
presente, o volume de gás nas bolhas aumenta de forma gradual, e,
conseqüentemente, também aumentam os comprimentos. Como os comprimentos
dos filmes de líquido abaixo das bolhas também aumentam, os comprimentos dos
pistões diminuem, pois os pistões perdem parte do líquido para o filme. Por outro
lado, se os efeitos de interação entre as bolhas são fortes, ocorre um crescimento
grande dos comprimentos de bolha, uma vez que a cada coalescência, a bolha
praticamente dobra o seu comprimento. Nesse caso, os comprimentos de pistões
tendem a aumentar também, pois na coalescência, o líquido do pistão que
desapareceu é redistribuído para os outros pistões. Além disso, a interação entre as
bolhas é afetada pela distância entre duas bolhas adjacentes (ou seja, o
comprimento do pistão que as separa), e, quanto menor o tamanho do pistão, maior
é a interação.
De forma geral, nos resultados para ar e água na horizontal, a queda de
pressão não é grande, diminuindo os efeitos de expansão. Por outro lado, devido à
baixa viscosidade da água, os efeitos de interação entre as bolhas são grandes, o
que gera muitas coalescências. Dessa forma, se os dados na entrada do tubo forem
corretos, espera-se que a simulação capture as coalescências ao longo do tubo.
Capítulo 5 Resultados para o escoamento horizontal 101
Porém, a quantidade de coalescências e a sua distribuição de ocorrência ao longo
do tubo também são importantes.
Nos resultados apresentados a seguir, os comprimentos de bolhas e pistões
simulados na entrada do tubo são semelhantes aos experimentos, com exceção de
uma pequena diferença nos comprimentos de bolhas para o caso A@W#4. No caso
A@W#2 a simulação não capturou adequadamente os efeitos de interação das
bolhas, gerando menos coalescências, o que mostra um crescimento menor dos
comprimentos simulados em relação ao experimental. A pequena ocorrência de
coalescências está relacionada também ao alto comprimento médio de pistões
(~10D), que impede que as bolhas coalescam ao longo do tempo que levam para
percorrer o tubo. Já no caso A@W#4 nota-se um efeito mais pronunciado da
coalescência das bolhas, e o crescimento dos comprimentos de bolhas e pistões ao
longo do tubo, justamente porque o comprimento médio dos pistões é menor (~5D),
e os resultados ficam mais próximos dos experimentais. Quanto à variação das
velocidades das bolhas ao longo do tubo, devido aos efeitos de expansão, a
velocidade da bolha tende a aumentar. Porém, quando ocorrem coalescências, a
velocidade média tende a diminuir, já que na coalescência, a bolha mais rápida
coalesce com a mais lenta, prevalecendo a velocidade da mais lenta. No entanto,
essas variações são muito pequenas tanto nos resultados experimentais como nas
simulações. Por fim, os resultados para a queda de pressão são mais dependentes
das propriedades dos fluidos envolvidos, e não apresentam variação grande para os
diversos casos.
Capítulo 5 Resultados para o escoamento horizontal 102
0 200 400 600L/D
0
40
80
120
L B/D
ExperimentalSlug tracking
0 200 400 600L/D
0
5
10
15
20
L S/D
0 200 400 600L/D
0.0
1.0
2.0
3.0
UT [
m/s
]
0 200 400 600L/D
0
30
60
90
120
P [k
Pa]
Figura 5.2 – Resultados médios para o escoamento horizontal de ar e água na
condição A@W#2
0 200 400 600L/D
0
5
10
15
20
25
L B/D
ExperimentalSlug tracking
0 200 400 600L/D
0
5
10
15
20
L S/D
0 200 400 600L/D
0.0
0.5
1.0
1.5
2.0
UT [
m/s
]
0 200 400 600L/D
0
30
60
90
120
P [k
Pa]
Figura 5.3 – Resultados médios para o escoamento horizontal de ar e água na
condição A@W#4
Capítulo 5 Resultados para o escoamento horizontal 103
Distribuições ao longo do tubo
Nas Figuras 5.4 e 5.5 são apresentadas as funções densidade de
probabilidade dos comprimentos de bolhas e pistões e velocidades de translação
das bolhas. São apresentados os resultados para as estações de medição 2, 3 e 4.
Os resultados são parecidos com os apresentados nos gráficos de médias,
porém pode-se analisar a distribuição das variáveis. No caso A@W2, nota-se que os
comprimentos de bolhas e pistões apresentam uma distribuição bi-modal,
principalmente próximo da entrada do tubo. Essa característica atrapalha o
simulador, que utiliza como condição de entrada uma distribuição normal ou log-
normal baseada no valor médio dos comprimentos experimentais na entrada. Dessa
forma, grande parte dos pistões pequenos que aparecem nos experimentos gera
coalescências. Já na simulação, não são considerados esses pistões pequenos, o
que diminui as coalescências. Já no caso A@W4 o comportamento das distribuições
experimentais é mais próximo de uma distribuição normal, o que torna os dados
utilizados na entrada do simulador mais confiáveis. Dessa forma, as coalescências
geradas na simulação são mais parecidas com as experimentais. Ainda, tanto nos
comprimentos de pistão como de bolha, a ordem de grandeza da distribuição (desvio
padrão) é semelhante na simulação e experimentos, o que reflete diretamente os
dados utilizados como condição de entrada da simulação. Nos resultados para
velocidade de translação das bolhas, ao longo da tubulação, as distribuições
experimentais ficam mais apertadas, como se a interação entre as bolhas
uniformizasse as velocidades. Na simulação esse efeito não ocorre, e a distribuição
na entrada tende a permanecer até a saída do tubo.
Capítulo 5 Resultados para o escoamento horizontal 104
0.000
0.005
0.010
0.015
- Est
ação
3
0 50 100 150 200LB/D
0.000
0.005
0.010
0.015
- Est
ação
2
0.000
0.005
0.010
0.015PD
F - E
staç
ão 4
ExperimentalSlug Tracking
0.000
0.025
0.050
0.075
0.100
- Est
ação
3
0 10 20 30 40 50LS/D
0.000
0.025
0.050
0.075
0.100
- Est
ação
2
0.000
0.025
0.050
0.075
0.100
- Est
ação
40
1
2
3
4
- Est
ação
3
0 1 2 3 4UT (m/s)
0
1
2
3
4
- Est
ação
2
0
1
2
3
4
- Est
ação
4
Figura 5.4 – Resultado das funções densidade de probabilidade para o escoamento horizontal de ar e água na condição
A@W#2
Capítulo 5 Resultados para o escoamento horizontal 105
0.00
0.02
0.04
0.06
0.08
- Est
ação
3
0 10 20 30 40 50LB/D
0.00
0.02
0.04
0.06
0.08
- Est
ação
2
0.00
0.02
0.04
0.06
0.08
- Est
ação
4ExperimentalSlug Tracking
0.00
0.05
0.10
0.15
- Est
ação
3
0 10 20 30 40LS/D
0.00
0.05
0.10
0.15PD
F - E
staç
ão 2
0.00
0.05
0.10
0.15
- Est
ação
40
2
4
6
- Est
ação
3
0 0.5 1 1.5 2UT (m/s)
0
2
4
6
- Est
ação
2
0
2
4
6
- Est
ação
4
Figura 5.5 - Resultado das funções densidade de probabilidade para o escoamento horizontal de ar e água na condição
A@W#4
Capítulo 5 Resultados para o escoamento horizontal 106
5.2.2 Resultados para o escoamento de ar e glicerina
Nessa seção são apresentados os resultados para o escoamento de ar e
glicerina na horizontal para as condições A@G#2 e A@G#4, e os resultados para as
demais condições são apresentados no Apêndice A. As seções de teste
experimental e numérica para a configuração de ar e glicerina na horizontal estão
mostradas na Figura 5.6.
Figura 5.6 – Representação das seções de teste experimental e numérica para o
escoamento de ar e glicerina na horizontal
Valores médios ao longo do tubo
Nas Figuras 5.7 e 5.8 são apresentados os resultados médios para o
escoamento de ar e glicerina. No caso do escoamento de glicerina, a interação entre
as bolhas não é tão forte como na água, devido à maior viscosidade da glicerina, o
que diminui as coalescências. Por outro lado, a queda de pressão ao longo do tubo
é maior, o que aumenta o efeito de expansão do gás. O efeito de expansão do gás é
mais fácil de se simular corretamente, e, dessa forma, quando os comprimentos de
bolhas e pistões simulados na entrada são corretos, eles tendem a permanecer
corretos ao longo do tubo. A velocidade de translação das bolhas também apresenta
pequena variação ao longo do tubo, tanto nos experimentos como nas simulações, à
exceção do caso A@G#2, em que as velocidades simuladas apresentam grande
variação, talvez devido a uma instabilidade numérica. Nota-se, também, que a queda
de pressão simulada tende a ser levemente menor que a experimental, o que gera
alguns erros na variação dos comprimentos devido à expansão do gás.
Capítulo 5 Resultados para o escoamento horizontal 107
0 200 400 600L/D
0
40
80
120
L B/D
ExperimentalSlug tracking
0 200 400 600L/D
0
5
10
15
20
L S/D
0 200 400 600L/D
0.0
1.0
2.0
3.0
UT [
m/s
]
0 200 400 600L/D
0
35
70
105
140
P [k
Pa]
Figura 5.7 – Resultados médios para o escoamento horizontal de ar e solução de
glicerina na condição A@G#2
0 200 400 600L/D
0
4
8
12
16
20
L B/D
ExperimentalSlug tracking
0 200 400 600L/D
0
5
10
15
L S/D
0 200 400 600L/D
0.0
0.5
1.0
1.5
2.0
UT [
m/s
]
0 200 400 600L/D
0
35
70
105
140
P [k
Pa]
Figura 5.8 – Resultados médios para o escoamento horizontal de ar e solução de
glicerina na condição A@G#4
Capítulo 5 Resultados para o escoamento horizontal 108
Distribuições das variáveis ao longo do tubo
Nas Figuras 5.9 e 5.10 são apresentadas as funções densidade de
probabilidade para o escoamento de ar e glicerina na horizontal. Nota-se que nos
resultados experimentais para o caso A@G#2 as distribuições são bi-modais e não
podem ser caracterizadas por curvas normais. Isso gera pouca confiabilidade nos
valores utilizados como condição de entrada, e pode ser responsável pelos maus
resultados médios. No caso A@G#4 as distribuições experimentais são mais bem
comportadas, o que dá mais confiabilidade à condição de entrada. Dessa forma,
para esses casos, as distribuições simuladas são mais coerentes com as
experimentais. A pequena defasagem entre as distribuições experimentais e
simuladas deve-se também à diferença entre as quedas de pressão simulada e
experimental, o que gera um erro no efeito de expansão do gás. Com relação às
distribuições para as velocidades de translação das bolhas, na maioria dos casos os
desvios padrão são parecidos, porém no caso A@G#2 eles são bem diferentes.
Novamente, esse caso foi afetado pelas distribuições experimentais estranhas.
Capítulo 5 Resultados para o escoamento horizontal 109
0.000
0.005
0.010
0.015
0.020
- Est
ação
3
0 50 100 150 200LB/D
0.000
0.005
0.010
0.015
0.020
- Est
ação
2
0.000
0.005
0.010
0.015
0.020PD
F - E
staç
ão 4
ExperimentalSlug Tracking
0.00
0.05
0.10
0.15
0.20
- Est
ação
3
0 10 20 30 40LS/D
0.00
0.05
0.10
0.15
0.20PD
F - E
staç
ão 2
0.00
0.05
0.10
0.15
0.20
- Est
ação
40
1
2
3
4
- Est
ação
3
0 1 2 3 4 5UT (m/s)
0
1
2
3
4
- Est
ação
2
0
1
2
3
4
- Est
ação
4
Figura 5.9 – Resultado das funções densidade de probabilidade para o escoamento horizontal de ar e solução de glicerina na
condição A@G#2
Capítulo 5 Resultados para o escoamento horizontal 110
0.00
0.02
0.04
0.06
0.08
0.10
0.12
- Est
ação
3
0 10 20 30LB/D
0.00
0.02
0.04
0.06
0.08
0.10
0.12
- Est
ação
2
0.00
0.02
0.04
0.06
0.08
0.10
0.12
- Est
ação
4 ExperimentalSlug Tracking
0.00
0.05
0.10
0.15
0.20
- Est
ação
3
0 10 20LS/D
0.00
0.05
0.10
0.15
0.20PD
F - E
staç
ão 2
0.00
0.05
0.10
0.15
0.20
- Est
ação
40
2
4
6
8
10
- Est
ação
3
0 0.5 1 1.5UT (m/s)
0
2
4
6
8
10
- Est
ação
2
0
2
4
6
8
10
- Est
ação
4
Figura 5.10 – Resultado das funções densidade de probabilidade para o escoamento horizontal de ar e solução de glicerina na
condição A@G#4
Capítulo 5 Resultados para o escoamento horizontal 111
5.2.3 Resultados para o escoamento de N2 e óleo SAE 20-50
Nessa seção são apresentados os resultados para o escoamento de
nitrogênio e óleo SAE 20-50 na horizontal para as condições N@O#1 e N@O#4, e
os resultados para as demais condições são apresentados no Apêndice A. As
seções de teste experimental e numérica para a configuração de nitrogênio e óleo
SAE 20-50 na horizontal estão mostradas na Figura 5.11. Nesse caso, a estação de
medição experimental localizada em L=130D não foi utilizada como condição de
entrada porque o escoamento ainda não estava totalmente no padrão de golfadas.
Figura 5.11 – Representação das seções de teste experimental e numérica para o
escoamento de nitrogênio e óleo SAE 20-50 na horizontal
Valores médios ao longo do tubo
As Figuras 5.12 e 5.13 apresentam os resultados médios para o escoamento
horizontal de nitrogênio e óleo SAE 20-50. Como a viscosidade do óleo é muito alta,
espera-se que a interação entre bolhas seja fraca. Além disso, a queda de pressão
deve ser alta, devido ao atrito.
Nos resultados apresentados, a simulação mostra um crescimento de bolhas
e pistões menor que o experimental. Isso pode estar relacionado ao erro que se
percebe na queda de pressão, que é transmitido para a expansão do gás. Por outro
lado, nota-se que tanto na simulação como nos experimentos, os comprimentos de
pistão também aumentam, o que indica coalescência. Essas coalescências são
justificadas porque com a grande expansão do gás as bolhas acabam por coalescer
devido ao seu grande crescimento. Nota-se também um grande erro no cálculo dos
comprimentos na entrada do tubo. Esse erro está relacionado a dificuldade que se
teve de ajustar alguns parâmetros nos dados experimentais na entrada do tubo,
Capítulo 5 Resultados para o escoamento horizontal 112
principalmente a constante C0 para cálculo da velocidade da bolha. Nos
experimentos, foi observado que em cada estação de medição, o valor de C0 é
diferente. Os efeitos de coalescência não foram captados pela simulação devido aos
comprimentos de pistão na entrada serem maiores que os experimentais, e as
coalescências ocorrem com maior freqüência quando os comprimentos de pistão
são menores.
Dessa forma, o crescimento dos comprimentos de bolhas e pistões está mais
relacionado à expansão do gás devido à queda de pressão. Nesse caso, também
devido à alta viscosidade, a queda de pressão ao longo do tubo é alta, o que gera
um grande crescimento do comprimento das bolhas. Em todos os casos a variação
dos comprimentos de bolha e pistão simulados foi parecida com os experimentais,
no entanto, os valores na entrada do tubo não foram muito próximos, o que acabou
se propagando para o resto do tubo. Também nota-se que a velocidade das bolhas
tem uma grande variação ao longo do tubo, tanto nos experimentos como na
simulação. Por fim, a queda de pressão simulada ficou acima da queda de pressão
experimental.
0 200 400L/D
0
10
20
30
L B/D
ExperimentalSlug tracking
0 200 400L/D
0
2
4
6
8
L S/D
0 200 400L/D
0.0
0.5
1.0
1.5
2.0
UT [
m/s
]
0 200 400L/D
0
35
70
105
140
175
P [k
Pa]
Figura 5.12 - Resultados médios para o escoamento horizontal de nitrogênio e óleo
SAE 20-50 N@O#1
Capítulo 5 Resultados para o escoamento horizontal 113
0 200 400L/D
0
4
8
12
L B/D
ExperimentalSlug tracking
0 200 400L/D
0
2
4
6
8
L S/D
0 200 400L/D
0.0
0.5
1.0
1.5
2.0
2.5
UT [
m/s
]
0 200 400L/D
0
60
120
180
240
P [k
Pa]
Figura 5.13 – Resultados médios para o escoamento horizontal de nitrogênio e óleo
SAE 20-50 N@O#4
Distribuições ao longo do tubo
As funções densidade de probabilidade para o escoamento de nitrogênio e
óleo são apresentadas nas Figuras 5.14 e 5.15. De forma geral, nota-se que as
distribuições simuladas e experimentais são parecidas. Porém, nota-se que as
distribuições experimentais para os comprimentos de bolha tendem para uma
distribuição log-normal, enquanto a simulada foi baseada em uma distribuição
normal, o que acaba prejudicando a simulação. Quanto aos comprimentos de pistão
e a velocidade das bolhas, as distribuições simuladas e experimentais são
semelhantes. As distribuições para os comprimentos de pistão apresentam erros
maiores na estação 3, justamente devido às coalescências que não ocoreram na
simulação. Por outro lado, os resultados para as velocidades das bolhas são muito
bons.
Capítulo 5 Resultados para o escoamento horizontal 114
0.00
0.02
0.04
0.06
0.08PD
F - E
staç
ão 3
0 20 40 60LB/D
0.00
0.02
0.04
0.06
0.08
- Est
ação
2ExperimentalSlug Tracking
0.0
0.1
0.2
0.3
- Est
ação
30 5 10 15 20
LS/D0.0
0.1
0.2
0.3
- Est
ação
2
0
1
2
3
4
- Est
ação
3
0 0.5 1 1.5 2 2.5UT (m/s)
0
4
8
12
- Est
ação
2
Figura 5.14 – Resultado das funções densidade de probabilidade para o escoamento horizontal de nitrogênio e óleo SAE 20-50
na condição N@O#1
Capítulo 5 Resultados para o escoamento horizontal 115
0.00
0.05
0.10
0.15
0.20
- Est
ação
3
0 10 20 30 40LB/D
0.00
0.05
0.10
0.15
0.20
- Est
ação
2ExperimentalSlug Tracking
0.0
0.1
0.2
0.3
- Est
ação
30 10 20 30
LS/D0.0
0.1
0.2
0.3
- Est
ação
2
0
2
4
6
8
10
- Est
ação
3
0 1 2 3UT (m/s)
0
2
4
6
8
10
- Est
ação
2
Figura 5.15 – Resultado das funções densidade de probabilidade para o escoamento horizontal de nitrogênio e óleo SAE 20-50
na condição N@O#4
Capítulo 6 Resultados para o escoamento vertical e inclinado 116
6 RESULTADOS PARA O ESCOAMENTO VERTICAL E INCLINADO
Nesse capítulo são apresentados os resultados obtidos nas simulações para
o escoamento vertical e inclinado. Na seção 6.1 são apresentadas as configurações
experimentais, que foram utilizadas também na simulação. A seção 6.2 apresenta a
comparação entre as simulações e os resultados experimentais, na qual são
apresentadas as comparações para os valores médios e as funções densidade de
probabilidade das variáveis ao longo do tubo.
6.1 Configurações experimentais
As configurações experimentais utilizadas nos escoamentos na vertical e
inclinados estão apresentadas na Tabela 6.1, enquanto as Tabelas 6.2, 6.3 e 6.4
apresentam as condições de velocidades superficiais e freqüência de entrada para
cada caso.
Tabela 6.1 – Definição das configurações experimentais
Ar e água (A@WV)
N2 e óleo SAE 20-50 (N@OV)
Ar e água inclinado (A@WI)
Comprimento do tubo (m) 5,811 (223,5 D) 5,587 (194 D) 5,811 (223,5 D) Diâmetro do tubo (m) 0,026 0,0288 0,026 Inclinação do tubo (º) 90 90 60
Viscosidade do líquido (Pa.s) 0,000855 0,2205 – 0,224 0,000855 Densidade do líquido (kg/m3) 1000 880 1000 Estação de medição #1 (m) ~ 0 ~ 0 ~ 0 Estação de medição #2 (m) 4,693 (180,5 D) 4,061 (141 D) 4,693 (180,5 D)
Tabela 6.2 – Condições de contorno para ar e água na vertical JL (m/s) JG (m/s) f (Hz)
A@WV#2 0,33 0,46 1,93 A@WV#5 0,30 1,42 1,91 A@WV#8 0,61 0,83 3,19
A@WV#10 0,88 0,60 4,42
Tabela 6.3 – Condições de contorno para nitrogênio e óleo SAE 20-50 JL (m/s) JG (m/s) f (Hz)
N@OV#2 0,17 0,45 3,09 N@OV#6 0,34 0,56 4,93 N@OV#8 0,18 0,90 3,61
Capítulo 6 Resultados para o escoamento vertical e inclinado 117
Tabela 6.4 – Condições de contorno para ar e água inclinado JL (m/s) JG (m/s) f (Hz)
A@WI#2 0,32 0,68 1,90 A@WI#5 0,59 0,42 2,81 A@WI#8 1,20 0,40 5,10 A@WI#10 1,21 1,40 5,93
6.2 Comparação entre as simulações e resultados experimentais
6.2.1 Resultados para o escoamento de ar e água na vertical
Nessa seção são apresentados os resultados para o escoamento de ar e
água na vertical para as condições A@WV#2 e A@WV#8, e os resultados para as
demais condições são apresentados no Apêndice A. As seções de teste
experimental e numérica para a configuração de ar e água na vertical estão
mostradas na Figura 6.1.
Figura 6.1 – Representação das seções de teste experimental e numérica para o
escoamento de ar e água na vertical
Valores médios ao longo do tubo
Nas Figuras 6.2 e 6.3 são apresentadas as comparações entre as simulações
e os experimentos para os valores médios. No escoamento vertical a queda de
pressão está muito ligada ao termo gravitacional, o que aumenta a expansão do gás.
Dessa forma a queda de pressão simulada fica muito próxima da experimental, pois
a queda de pressão gravitacional não depende dos comprimentos, mas sim da
vazão.
Com relação aos comprimentos de bolha, em todos as condições a
comparação entre simulação e experimento é boa. Nesses casos, os comprimentos
de bolha na entrada do tubo são calculados corretamente e, ao longo do tubo, a
Capítulo 6 Resultados para o escoamento vertical e inclinado 118
simulação descreve bem os fenômenos de expansão do gás no interior das bolhas e
coalescência de bolhas, que resultam no aumento do comprimento médio das
bolhas ao longo do tubo.
Quanto aos comprimentos de pistão, todos os resultados apresentam as
mesmas características. Na entrada do tubo, os comprimentos de pistão médios
simulados têm um valor um pouco abaixo dos experimentais. Além disso, ao longo
do tubo, o grande aumento dos pistões visto nos experimentos não é apresentado
também na simulação, resultando em um erro maior nas proximidades da saída do
tubo. Nota-se uma grande ocorrência de coalescência, por isso o aumento dos
comprimentos de pistões. No entanto, apesar dos comprimentos de bolhas variarem
da mesma forma que os experimentais, os comprimentos de pistão não aumentam
tanto. Isso pode indicar que se deve modelar melhor o que ocorre com as bolhas e
pistões no momento da coalescência, como a correta distribuição das massas de
líquido e gás.
Analisando-se as velocidades de bolhas médias e pressões médias, a
comparação entre a simulação e o experimento é boa em todas as condições. Esse
resultado mostra que essas variáveis não sofrem tanta influência dos comprimentos
de bolhas e pistões.
Capítulo 6 Resultados para o escoamento vertical e inclinado 119
0 100 200L/D
0
5
10
15
20
L B/D
ExperimentalSlug tracking
0 100 200L/D
0
5
10
15
20
L S/D
0 100 200L/D
0.0
0.5
1.0
1.5
2.0
UT [
m/s
]
0 100 200L/D
0
35
70
105
140
P [k
Pa]
Figura 6.2 – Resultados médios para o escoamento vertical de ar e água na
condição A@WV#2
0 100 200L/D
0
10
20
30
L B/D
ExperimentalSlug tracking
0 100 200L/D
0
5
10
15
20
L S/D
0 100 200L/D
0.0
1.0
2.0
3.0
UT [
m/s
]
0 100 200L/D
0
40
80
120
160
P [k
Pa]
Figura 6.3 – Resultados médios para o escoamento vertical de ar e água na
condição A@WV#8
Capítulo 6 Resultados para o escoamento vertical e inclinado 120
Distribuição das variáveis ao longo do tubo
Além dos valores médios das variáveis na entrada e nas proximidades da
saída do tubo, as funções densidade de probabilidade dessas variáveis nesses
pontos também podem ser analisadas. Esses resultados são apresentados nas
Figuras 6.4 e 6.5.
Considerando-se os resultados para os comprimentos de bolha, novamente
percebe-se um bom acordo entre os resultados simulados e experimentais. Porém,
de forma geral, os resultados para todas as condições apresentam bom acordo das
distribuições, inclusive analisando-se a amplitude dos comprimentos de bolha que
ocorrem nos pontos mostrados.
Os resultados para os comprimentos de pistão mostram que, na entrada,
apesar da pequena diferença nos valores médios, as distribuições estão
semelhantes às experimentais. No entanto, próximo da saída do tubo, as
distribuições experimentais ficam mais deslocadas para a direita, e com ocorrência
de comprimentos de pistão maiores que os simulados, fato já observado nos valores
médios.
Quanto aos resultados para as velocidades de bolha, na entrada do tubo as
distribuições são muito parecidas. Porém, nas proximidades da saída do tubo, as
distribuições experimentais tendem a ter desvios maiores que as distribuições
simuladas, apesar dos valores médios simulados estarem em acordo com os
experimentos.
Capítulo 6 Resultados para o escoamento vertical e inclinado 121
0.00
0.04
0.08
0.12
0.16PD
F - E
staç
ão 2
4 8 12 16 20LB/D
0.00
0.04
0.08
0.12
0.16
- Est
ação
1ExperimentalSlug Tracking
0.00
0.05
0.10
0.15
- Est
ação
20 10 20 30 40 50
LS/D0.00
0.05
0.10
0.15
- Est
ação
1
0
1
2
3
4
5
- Est
ação
2
0 0.5 1 1.5 2 2.5UT (m/s)
0
1
2
3
4
5
- Est
ação
1
Figura 6.4 – Resultado das funções densidade de probabilidade para o escoamento vertical de ar e água na condição
A@WV#2
Capítulo 6 Resultados para o escoamento vertical e inclinado 122
0.00
0.05
0.10
0.15
- Est
ação
2
0 10 20 30 40LB/D
0.00
0.05
0.10
0.15
- Est
ação
1ExperimentalSlug Tracking
0.00
0.05
0.10
0.15
- Est
ação
20 10 20 30 40 50
LS/D0.00
0.05
0.10
0.15
- Est
ação
1
0
1
2
3
- Est
ação
2
0 1 2 3 4UT (m/s)
0
1
2
3
- Est
ação
1
Figura 6.5 – Resultado das funções densidade de probabilidade para o escoamento vertical de ar e água na condição
A@WV#8
Capítulo 6 Resultados para o escoamento vertical e inclinado 123
6.2.2 Resultados para o escoamento de N2 e óleo SAE 20-50 na vertical
Nessa seção são apresentados os resultados para o escoamento de
nitrogênio e óleo SAE 20-50 na vertical para as condições N@OV#2 e N@OV#6, e
os resultados para as demais condições são apresentados no Apêndice A. As
seções de teste experimental e numérica para a configuração de nitrogênio e óleo
SAE 20-50 na vertical estão mostradas na Figura 6.6.
Figura 6.6 – Representação das seções de teste experimental e numérica para o
escoamento de nitrogênio e óleo SAE 20-50 na vertical
Valores médios ao longo do tubo
Nas Figuras 6.7 e 6.8 são apresentados os valores médios na entrada e
próximo da saída do tubo. Com relação aos comprimentos de bolha, na entrada do
tubo os valores médios simulados são semelhantes dos experimentais. Próximo da
saída do tubo, existe uma pequena diferença entre os resultados simulados e
experimentais nas condições N@OV#2 e N@OV#6, em que os resultados simulados
apresentam valores menores que os experimentais. Isso mostra que a simulação
não capturou as coalescências que ocorreram nos experimentos.
Analisando-se os resultados para os comprimentos médios dos pistões, nota-
se novamente que, na entrada do tubo, os resultados simulados são bons, enquanto
na saída do tubo há uma diferença em relação aos experimentais. Além disso, nota-
se que nos experimentos, os comprimentos de pistão aumentam ao longo do tubo,
enquanto nas simulações, eles diminuem. Esse resultado pode ser explicado pela
coalescência de bolhas nos experimentos, que tende a aumentar bastante os
comprimentos de bolhas e pistões. Nas simulações, não acontecem tantas
Capítulo 6 Resultados para o escoamento vertical e inclinado 124
coalescências, portanto os comprimentos de bolha não aumentam e os
comprimentos de pistão até diminuem ao longo do tubo.
Os resultados médios para as velocidades de bolhas e pressões mostram,
novamente, que essas variáveis não são tão influenciadas pelos erros nos
comprimentos de pistões e bolhas. Dessa forma, os resultados simulados e
experimentais ficam bem próximos, tanto na entrada do tubo, quanto próximo da
saída.
0 50 100 150L/D
0
5
10
15
20
L B/D
ExperimentalSlug tracking
0 50 100 150L/D
0
2.5
5
7.5
10
L S/D
0 50 100 150L/D
0.0
0.5
1.0
1.5
2.0
UT [
m/s
]
0 50 100 150L/D
0
50
100
150
200
P [k
Pa]
Figura 6.7 – Resultados médios para o escoamento vertical de nitrogênio e óleo SAE
20-50 na condição N@OV#2
Capítulo 6 Resultados para o escoamento vertical e inclinado 125
0 50 100 150L/D
0
5
10
15
L B/D
ExperimentalSlug tracking
0 50 100 150L/D
0
2.5
5
7.5
10
L S/D
0 50 100 150L/D
0.0
1.0
2.0
3.0
UT [
m/s
]
0 50 100 150L/D
0
50
100
150
200
250
P [k
Pa]
Figura 6.8 – Resultados médios para o escoamento vertical de nitrogênio e óleo SAE
20-50 na condição N@OV#6
Distribuição das variáveis ao longo do tubo
As funções densidade de probabilidade são apresentadas nas Figuras 6.9 e
6.10. Para os comprimentos de bolhas e pistões, as distribuições da simulação na
entrada do tubo são muito semelhantes às distribuições experimentais. No entanto,
nas proximidades da saída, as distribuições experimentais apresentam desvio maior
que as simuladas, além de um deslocamento do valor médio para a direita. Esses
resultados são típicos da ocorrência de coalescências ao longo do tubo, que geram
uma distribuição maior das variáveis apresentadas. Nesse caso, a simulação não
capturou essas coalescências de maneira correta. Quanto às velocidades de bolhas,
na condição N@OV#2 as distribuições simuladas e experimental, tanto na entrada
como próximo da saída do tubo, são semelhantes. Isso mostra que a velocidades da
bolha não parece sofrer muita variação para essas condições. Na condição
N@OV#6, logo na entrada do tubo a distribuição simulada já está deslocada em
relação à experimental, o que acaba se repetindo e até aumentando nas
distribuições próximo da saída do tubo.
Capítulo 6 Resultados para o escoamento vertical e inclinado 126
0.00
0.04
0.08
0.12
0.16PD
F - E
staç
ão 2
0 5 10 15 20 25LB/D
0.00
0.04
0.08
0.12
0.16
- Est
ação
1ExperimentalSlug Tracking
0.00
0.05
0.10
0.15
0.20
0.25
- Est
ação
20 10 20
LS/D0.00
0.05
0.10
0.15
0.20
0.25
- Est
ação
1
0
1
2
3
4
- Est
ação
2
0 0.5 1 1.5 2UT (m/s)
0
1
2
3
4
- Est
ação
1
Figura 6.9 – Resultado das funções densidade de probabilidade para o escoamento vertical de nitrogênio e óleo SAE 20-50 na
condição N@OV#2
Capítulo 6 Resultados para o escoamento vertical e inclinado 127
0.00
0.04
0.08
0.12
0.16
- Est
ação
2
0 10 20 30LB/D
0.00
0.04
0.08
0.12
0.16
- Est
ação
1ExperimentalSlug Tracking
0.00
0.05
0.10
0.15
0.20
- Est
ação
20 10 20 30
LS/D0.00
0.05
0.10
0.15
0.20
- Est
ação
1
0.0
0.5
1.0
1.5
2.0
2.5
- Est
ação
2
0 1 2 3 4UT (m/s)
0.0
0.5
1.0
1.5
2.0
2.5
- Est
ação
1
Figura 6.10 – Resultado das funções densidade de probabilidade para o escoamento vertical de nitrogênio e óleo SAE 20-50
na condição N@OV#6
Capítulo 6 Resultados para o escoamento vertical e inclinado 128
6.2.3 Resultados para o escoamento de ar e água inclinado
Nessa seção são apresentados os resultados para o escoamento de ar e
água inclinado para as condições A@WI#2 e A@WI#5, e os resultados para as
demais condições são apresentados no Apêndice A. As seções de teste
experimental e numérica para a configuração de ar e água inclinado estão
mostradas na Figura 6.11.
Figura 6.11 – Representação das seções de teste experimental e numérica para o
escoamento de ar e água inclinado
Valores médios ao longo do tubo
Nas Figuras 6.12 e 6.13 são apresentados os resultados para os valores
médios. Para os comprimentos de bolhas, no caso A@WI#5, e os resultados
simulados são bem próximos dos experimentais. No caso A@WI#2 existe uma
pequena diferença entre a simulação e o experimento na entrada do tubo. Com
relação aos comprimentos de pistão, nas condições A@WI#2 e A@WI#5, os valores
na entrada do tubo apresentam uma pequena diferença com os valores
experimentais. Além disso, ao longo do tubo, o caso A@WI#5 apresenta uma
variação do comprimento médio dos pistões semelhante à experimental, ou seja, o
mesmo erro apresentado na entrada do tubo se repete nas proximidades da saída.
Para o caso A@WI#2 a variação dos comprimentos de pistão ao longo do tubo não
é tão grande como na simulação, o que gera um aumento da diferença entre
simulação e experimento próximo da saída do tubo.
As velocidades de bolhas e pressões apresentam resultados simulados
parecidos com os resultados experimentais.
Capítulo 6 Resultados para o escoamento vertical e inclinado 129
0 100 200L/D
0
10
20
30
L B/D
ExperimentalSlug tracking
0 100 200L/D
0
5
10
15
20
L S/D
0 100 200L/D
0.0
0.5
1.0
1.5
2.0
UT [
m/s
]
0 100 200L/D
0
35
70
105
140
P [k
Pa]
Figura 6.12 Resultados médios para o escoamento inclinado de ar e água na
condição A@WI#2
0 100 200L/D
0
5
10
15
L B/D
ExperimentalSlug tracking
0 100 200L/D
0
5
10
15
20
L S/D
0 100 200L/D
0.0
0.5
1.0
1.5
2.0
UT [
m/s
]
0 100 200L/D
0
50
100
150
P [k
Pa]
Figura 6.13 – Resultados médios para o escoamento inclinado de ar e água na
condição A@WI#5
Capítulo 6 Resultados para o escoamento vertical e inclinado 130
Distribuição das variáveis ao longo do tubo
As funções densidade de probabilidade são apresentadas nas Figuras 6.14 e
6.15. Analisando-se os resultados para os comprimentos de bolhas, nota-se que as
distribuições simuladas e experimentais são bem parecidas. Com relação aos
comprimentos de pistões, em geral, as distribuições para todos os casos na
proximidade da saída do tubo apresentam um deslocamento dos resultados
experimentais em relação aos simulados. Quanto às velocidades de bolhas, nota-se
que já na entrada do tubo as distribuições simuladas e experimentais apresentam
algumas diferenças, tanto nos valores médios como nas distribuições. Essa
característica é notada também na saída do tubo, em que os resultados
experimentais apresentam um alargamento das distribuições, provavelmente devido
à interações entre as bolhas ao longo do tubo.
Capítulo 6 Resultados para o escoamento vertical e inclinado 131
0.00
0.04
0.08
0.12PD
F - E
staç
ão 2
0 10 20 30LB/D
0.00
0.04
0.08
0.12
- Est
ação
1ExperimentalSlug Tracking
0.00
0.04
0.08
0.12
- Est
ação
20 10 20 30 40
LS/D0.00
0.04
0.08
0.12
- Est
ação
1
0
2
4
6
- Est
ação
2
0 1 2 3UT (m/s)
0
2
4
6
- Est
ação
1
Figura 6.14 – Resultado das funções densidade de probabilidade para o escoamento inclinado de ar e água na condição
A@WI#2
Capítulo 6 Resultados para o escoamento vertical e inclinado 132
0.0
0.1
0.2
0.3
0.4
- Est
ação
2
0 5 10 15 20LB/D
0.0
0.1
0.2
0.3
0.4
- Est
ação
1ExperimentalSlug Tracking
0.00
0.05
0.10
0.15
- Est
ação
20 10 20 30 40
LS/D0.00
0.05
0.10
0.15
- Est
ação
1
0
2
4
6
- Est
ação
2
0 0.5 1 1.5 2UT (m/s)
0
2
4
6
- Est
ação
1
Figura 6.15 – Resultado das funções densidade de probabilidade para o escoamento inclinado de ar e água na condição
A@WI#5
Capítulo 7 Conclusões 133
7 CONCLUSÕES
Nesse trabalho, foi apresentado um modelo de seguimento de pistões
considerando-se a influência da inércia do pistão de líquido, fluxos de quantidade de
movimento através das células e influência do filme de líquido que escoa ao lado da
bolha alongada. As condições de contorno e o processo de inserção de células na
tubulação foram trabalhados de forma a estimular a interação entre bolhas
consecutivas, na tentativa de capturar a intermitência intrínseca do escoamento.
Para isso, foi implementado um processo de inserção de células baseado em
distribuições, ou seja, cada célula que entra no domínio computacional apresenta
características diferentes.
A implementação da modelagem com as condições de contorno se mostrou
eficiente na captura dos principais efeitos característicos do escoamento em
golfadas, responsáveis pelas variações dos comprimentos: expansão do gás ao
longo do escoamento devido à queda de pressão e interação entre bolhas. Para
validação da metodologia, foram simulados diversos casos de escoamento
previamente analisados experimentalmente, para escoamentos horizontais,
inclinados e verticais, utilizando-se fluidos com grande variação de viscosidade,
como água, glicerina e óleo.
Os resultados foram analisados em termos de parâmetros médios e
distribuições das variáveis em posições ao longo da tubulação. A evolução desses
parâmetros ao longo da tubulação permite a visualização dos efeitos concorrentes
de expansão do gás e interação entre bolhas.
Os resultados obtidos das simulações para o escoamento em golfadas são
avaliadas em três etapas: cálculo adequado das condições de entrada das células,
cálculo correto da queda de pressão ao longo do tubo e conseqüente expansão do
gás ao longo da tubulação, e correta modelagem da interação entre as bolhas.
Observou-se que na evolução do escoamento ao longo da tubulação, os
efeitos de expansão do gás e interação das bolhas são concorrentes, mas também
interdependentes, uma vez que a interação entre as bolhas também é influenciada
pelo efeito de esteira. O efeito de esteira é função dos comprimentos dos pistões,
que por sua vez dependem das condições de entrada e da expansão do gás.
Capítulo 7 Conclusões 134
A variação dos comprimentos devido à expansão do gás ao longo da
tubulação é um fenômeno de modelagem mais simples, pois depende basicamente
da variação da pressão. A interação entre as bolhas é um fenômeno mais complexo
de se modelar, justamente pela aleatoriedade do fenômeno. Por isso, em alguns
casos, os resultados simulados não se adequaram aos dados experimentais.
De forma geral, se não houver interação entre bolhas, a tendência é de
crescimento do comprimento das bolhas ao longo da tubulação, de acordo com o
aumento no volume do gás. Ocorre, também, pequena diminuição do comprimento
dos pistões, pois parte do líquido que estava nos pistões é transferido para a região
do filme do líquido, que aumenta junto com a bolha.
Quando a interação entre as bolhas provoca várias coalescências, o
comprimento das bolhas apresenta um crescimento muito maior ao longo da
tubulação. Porém, nesse caso, os comprimentos dos pistões também aumentam,
pois a cada coalescência o volume de líquido originário de um pistão é transferido
para os demais pistões.
Observou-se que o escoamento de ar e água (baixa viscosidade) apresenta
maiores taxas de coalescência. Conforme a viscosidade aumenta (no escoamento
de glicerina ou óleo) as coalescências diminuem e a queda de pressão aumenta, o
que intensifica o efeito da expansão do gás.
No escoamento vertical, a queda de pressão é muito maior que no horizontal,
o que gera um efeito maior de expansão do gás e mascara possíveis efeitos de
interação das bolhas. Os resultados obtidos para o escoamento inclinado de ar e
água foram muito semelhantes aos resultados para o escoamento vertical. Isso
mostra que, para a inclinação utilizada (60º) o escoamento já se comporta da
mesma forma que no vertical.
De forma geral, o trabalho apresenta uma contribuição para o
desenvolvimento de simuladores do escoamento no padrão de golfadas, focando
principalmente na interação entre as bolhas. A metodologia utilizada se mostrou
eficiente na caracterização da entrada de células e captura da coalescência.
Alguns pontos que não foram abordados nesse trabalho e podem ser focados
em trabalhos futuros são:
• Estudo aprofundado dos modelos de esteira. Esse efeito é crucial para
modelar corretamente as coalescências. No entanto, um modelo genérico
Capítulo 7 Conclusões 135
para a função de esteira deve ser função de diversas propriedades dos
fluidos e condições de escoamento.
• Estudos mais criteriosos para a geração das listas utilizadas como
condição de entrada da simulação.
• Modelagem do escoamento em tubulações com variações de inclinação,
que é o caso mais comum na indústria.
• Simulação e comparação com dados de tubulações de maiores
comprimentos e diâmetros.
• Modelagem de escoamentos com mudança de fase, ou seja, conforme a
pressão cai ocorre a passagem para a fase gás de parte do líquido. Esse
fenômeno é muito comum, e apresenta um novo efeito para o
crescimento dos pistões e bolhas, concorrente com a expansão do gás e
coalescência.
• Definição de parâmetros para transição do escoamento para outros
padrões, que podem ocorrer ao longo da tubulação.
Referências bibliográficas 136
REFERÊNCIAS BIBLIOGRÁFICAS
Abdul-Majeed G. H., (2000) “Liquid slug holdup in horizontal and slightly inclined two-
phase slug flow”. Journal of Petroleum Science and Engineering, 27, 27-32.
Al-safran, E., Taitel, Y. e Brill, J., (2004) "Prediction of Slug Length Distribution along
a Hilly-Terrain Pipeline Using Slug Tracking Model" ASME J. Energy Resources
Technology, vol. 127, pp. 54-61, 2004.
Andreussi P., Bendiksen K. H., e Nydal O. J., (1993) “Void Distribution in Slug Flow”.
Int. J. Multiphase Flow, 19 (5), 817-828.
Barnea D. e Brauner N., (1985) “Holdup of the liquid slug in two-phase intermittent
flow”. Int. J. Multiphase Flow, 11 (1), 43-49.
Barnea D. e Taitel Y., (1993) "A model for slug lenght distribution in gas-liquid slug
flow" Int. J. Multiphase Flow, 19 (5): 829-838.
Bendiksen K. H., (1984) “An experimental investigation of the motion of long bubbles
in inclined tubes”. Int. J. Multiphase Flow, 10, 467-483.
Bendiksen, K. H., Maines, D., Moe, R. e Nuland, S. (1991) “The dynamic two-fluid
model OLGA; Theory and application” SPE Production Engineering Vol. 6, No. 2,
pp. 171-180.
Benjamin, (1968) “Gravity currents and related phenomena”. J. Fluid Mech., 31 (2),
209-248.
Box, G. e Muller, M., (1958) “A Note on the Generation of Random Normal Deviates”
Annals Math. Statistics, vol. 29, pp. 610-611, 1958.
Bugg, J. D., Mack, K., Rezkallah, K. S. (1998) “A numerical model of Taylor bubbles
rising through stagnant liquids in vertical tubes” International Journal of Multiphase
Flow, Vol. 24, No. 2, pp. 271-281, 1998.
Davis R. M., e Taylor G. I., (1950) “The mechanics of large bubbles rising through
extended liquids and through liquids in tubes”. Proc Roy. Soc., 200A, 375-390.
Dukler A. E., e Hubbard M. G., (1975) “A model for gas-liquid slug flow in horizontal
and near horizontal tubes”. Ind. Eng. Chem. Fundam., 14 (4), 337-347.
Dukler A. E., Moalem Maron D., e Brauner N., (1985) “A physical model for predicting
the minimum stable slug length”. Chem. Eng. Sci., 40, 1379-1385.
Referências bibliográficas 137
Fagundes Netto J. R., (1999), “Dynamique de poches de gaz isolées en écoulement
permanent et non permanent horizontal”. Toulouse: Institut de Mécanique des
Fluides de Toulouse, Institut National Polytechnique de Toulouse, 166p. Tese
(Doutorado).
Fernandes, R. C., Semiat, R., Dukler, A. E. Hydrodynamic model for gas-liquid slug
flow in vertical tubes. AIChE Journal, Vol. 29, No. 6, pp. 981-989, 1983.
Ferré D., (1979) “Ecoulements diphasiques a poches en conduite horizontale”. Rev.
Inst. Fr. Pet., 34, 113-142.
Ferschneider G, (1983) “Ecoulement diphasiques gas-liquide a poches et a
bouchons en conduits”. Revue Inst. Fr. Petrole, 38, 153.
Fortuna, A. O. (2000) Métodos Numéricos para Dinâmica Computacional de Fluidos.
1. ed. São Paulo: EDUSP, 2000. v. 1. 400 p.
Fox, R. W. e McDonald, A. T., (2006) "Introdução à Mecânica dos Fluidos", 6a ed.,
Editora LTC, 2006.
Franklin, E. M. e Rosa, E. S. (2004) “Dynamic slug tracking model for horizontal gás-
liquid flow” Proc. of the 3rd International Symposium on Two-phase Flow
Modelling and Experimentation, Pisa, 22-24 September 2004.
Gomez L. E., Shoham O. e Taitel Y., (2000) “Prediction of slug liquid holdup:
horizontal to upward vertical flow”. Int. J. Multiphase Flow, 26, 517-521.
Gregory G. A., e Scott D. S., (1969) “Correlation of liquid slug velocity and frequency
in horizontal co-current gas-liquid slug flow”. AIChE Journal, 15 (6), 933-935.
Gregory G. A., Nicholson M. K., e Aziz K., (1978) “Correlation of the liquid volume
fraction in the slug for horizontal gas-liquid slug flow”. Int. J. Multiphase Flow, 4,
33-39.
Grenier P., (1997) “Evolution des longueurs de bouchons en écoulement intermittent
horizontal”. Toulouse: Institut de Mécanique des Fluides de Toulouse, Institut
National Polytechnique de Toulouse, 193p. Tese (Doutorado).
Henau, V., Raithby, G. D. (1995) “A transient two-fluid model for the simulation of
slug flow in pipelines – I. Theory” International Journal of Multiphase Flow, Vol. 21,
No. 3, pp. 335-349, 1995.
Heywood N. I., e Richardson J. F., (1979) “Slug flow of air-water mixtures in a
horizontal pipe: determination of liquid holdup by -ray absorption”. Chem. Eng.
Sci., 34, 17-30.
Referências bibliográficas 138
Hill T. J., e Wood D. G., (1990) “A new approach to the prediction of slug frequency”.
65th Annual Technical Conference and Exhibition of the Society of Petroleum
Engineers, September 23-26, New Orleans, USA, SPE 20629, 141-149.
Hill, T. J. e Wood, D. G. (1994) “Slug flow: occurrence, consequences and prediction”
University of Tulsa Centennial Petroleum Engineering Symposium, 29-31 August
1994, Tulsa, Oklahoma.
Ishii, M. (1975) “Thermo-fluid dynamic theory of two-phase flow” Collection de la
Direction des Etudes et Recherches d'Electricite de France, no. 22, Eyrolles,
Paris, France, 1975.
Issa, R. I., Kempf, M. H. W. Simulation of slug flow in horizontal and nearly horizontal
pipes with the two-fluid model. International Journal of Multiphase Flow, Vol. 29,
No. 1, pp. 69-95, 2003.
Malnes D., (1982) “Slug flow in vertical, horizontal and inclined pipes”. Report
IFE/KR/E-83/002, V. Inst. for Energy technology, Kjeller, Norway.
Manolis I. G., (1995) “High Pressure Gas-Liquid Slug Flow”. PhD thesis, Department
of Chemical Engineering and Chemical Technology, Imperial College of Science,
Technology & Medicine, UK.
Mao, Z., Dukler, A. E. An experimental study of gas-liquid slug flow. Experiments in
Fluids, Vol. 8, pp. 169-182, 1989.
Marcano R., Chen T. X., Sarica C., e Brill J. P., (1998) “A Study of Slug
Characteristics for Two-Phase Horizontal Flow”. SPE 39856, 213-219.
Moissis, R. e Griffith, P. (1962) “Entrance effects in a two-phase slug flow” J. Heat
Transfer Vol. 84, pp. 29–39.
Nicklin D. J., Wilkes J. O., e Davidson J. F., (1962) “Two-phase flow in vertical
tubes”. Trans. Inst. Chem. Eng., 40, 61-68.
Nogueira, S., Riethmuler, M. L., Campos, J. B. L. M., Pinto, A. M. F. R. Flow in the
nose region and annular film around a Taylor bubble rising through vertical
columns of stagnant and flowing Newtonian liquids. Chemical Engineering
Science, Vol. 61, No. 2, pp. 845-857, 2006.
Nydal, O. J. e Banerjee, S., (1995) "Object oriented dynamic simulation of slug flow",
Proceedings of the 2nd Int. Conf. Multiphase Flow, Kyoto, vol. 2, pp. IF2_7-14.
Referências bibliográficas 139
Nydal, O. J., Klebert, P., Wangensteen, T. e Kristiansen, O. (2003) “Transient two
phase flow model” Proc. of the 11th BHRG Multiphase Technology '03, San Remo,
Italy: 11-13 June 2003.
Omgba-Essama C., Hanich L., Thompson C. P., and Lezeau P., (2000) “Adaptive
Grid Refinement for Transient Two-Phase Flows”. AMIF-ESF Workshop,
Computing Methods for Two-phase Flows, Centre Paul Langevin, Aussois
(Savoie) France, 12-14 January.
Pauchon, C. L., Dhulesia, H. (1994) “TACITE: A transient tool for multiphase pipeline
and well simulation”, SPE Annual Technical Conference and Exhibition, 25-28
September 1994, New Orleans, Louisiana.
Petalas N., e Aziz K., (1998), “A mechanistic model for multiphase flow in pipes”.
49th Annual Technical Meeting of Petroleum Society of the Canadian Institute of
Mining, Metallurgy & Petroleum. Paper 98-39, Calgary, Alberta June 8-10,
Canada.
Pinto, A. M. F. R., Campos, J. B. L. M. Coalescence of two gas slugs rising in a
vertical column of liquid. Chemical Engineering Science, Vol. 51, No. 1, pp. 45-54,
1996.
Polonsky, S., Shemer, L., Barnea, D. (1999) “The relation between the Taylor bubble
motion and the velocity field ahead of it” International Journal of Multiphase Flow,
Vol. 25, No. 6, pp. 957-975, 1999.
Rodrigues, H. T. (2006) “Simulação do escoamento bifásico líquido-gás intermitente
em golfadas utilizando o modelo de seguimento dinâmico de pistões” Curitiba:
Universidade Tecnológica Federal do Paraná, 115p. monografia (graduação).
Rosa, E. S., (2006) “Análise de Escoamentos em Golfadas de Óleos Pesados e de
Emulsões Óleo-água” Quarto Relatório de Atividades e Resultados alcançados,
Projeto UNICAMP / CENPES – PETROBRAS, Campinas, SP.
Sakaguchi, H., (2001) "Correlations for large bubble lenght, liquid slug lenght, slug
unit lenght and slug period of gas-liquid two-phase slug flow in vertical pipes",
Proc. of 4th International Conference of Multiphase Flow, ICMF 2001.
Shoham, O. (2006) “Mechanistic Modeling of Gas-Liquid Two-phase Flow In Pipes”
Society of Petroleum Engineers, Richardson, TX, 396p.
Taha, T., Cui, Z. F. CFD modelling of slug flow in vertical tubes. Chemical
Engineering Science, Vol. 61, pp. 676-687, 2006.
Referências bibliográficas 140
Taitel, Y. e Barnea, D. (1990a) "Two phase slug flow", Advances in Heat Transfer,
Hartnett J.P. and Irvine Jr. T.F. ed., vol. 20, 83-132, Academic Press (1990).
Taitel, Y., Barnea, (1990b) D. “A consistent approach for calculating pressure drop in
inclined slug flow” Chemical Engineering Science, Vol. 45, No. 5, pp. 1199-1206,
1990b.
Taitel, Y., Barnea, D. Effect of gas compressibility on a slug tracking model. Chemical
Engineering Science, Vol. 53, No. 11, pp. 2089-2097, 1998.
Taitel, Y. e Barnea D. (2000) “Slug Tracking Model for Hilly Terrain Pipelines”, SPE
Journal, Vol. 5, p.102-109 (2000).
Taitel Y., e Dukler A. E., (1977) “A model for slug frequency during gas-liquid flow in
horizontal and near horizontal pipes”. Int. J. Multiphase Flow, 3 (6), 585-596.
Talvy, C. A., Shemer, L., Barnea, D. On the interaction between two consecutive
elongated bubbles in a vertical pipe. International Journal of Multiphase Flow, Vol.
26, No. 12, pp. 1905-1923, 2000.
Théron B., (1989) “Ecoulements diphasiques instationnaires en conduite horizontale”.
PhD thesis, Institut National Polytechnique de Toulouse, France.
Tronconi E., (1990) “Prediction of slug frequency in horizontal two-phase flow”.
AIChE Journal, 36 (5), 701-709.
Ujang, P. M., Lawrence, C. J. e Hewitt, G. F. (2006) “Conservative incompressible
slug tracking model for gas-liquid flow in a pipe” Proc. of the 5th BHRG North
American Conference on Multiphase Technology, Banff, Canada: 31 May-2 June
2006.
van Hout, R., Gulitski, A., Barnea, D. e Shemer, L. (2002) “Experimental investigation
of the velocity field induced by a Taylor bubble rising in stagnant water” Int. J.
Multiphase Flow, Vol. 28 p.579-596, 2002.
Versteeg, H. K. e Malalasekera, W. (1995) “An introduction to computational fluid
dynamics: The finite volume method” Longman Pub Group, 257p., 1995.
Wallis G. B., (1969) “One dimensional two-phase flow”. New York, McGraw-Hill.
Weber M. E., (1981) “Drift in Intermittent Two-Phase Flow in Horizontal Pipes”. Can.
J. Chem. Eng., 59, 398-399.
Woods D. B., e Hanratty T. J., (1996) “Relation of Slug stability to shedding rate”. Int.
J. Multiphase Flow, 22 (5), 809-828.
Referências bibliográficas 141
Zabaras G. J., (2000) "Prediction of Slug Frequency for Gas/Liquid Flows", SPE
Journal 5 (3), 2000. Shell E&P Technology Co.
Zhang H., Wang Q., Sarica C. e Brill J. P. (2003) "A unified mechanistic model for
slug liquid holdup and transition between slug and dispersed bubble flows", Int. J.
of Multiphase Flow, 29, 97-107.
Zheng, G., Brill, J. P. e Taitel, Y. (1994) “Slug flow behavior in a hilly terrain pipeline”
Int. J. Multiphase Flow Vol. 20, No. 1, pp. 63-79, 1994.
Apêndice A – Resultados para diversas condições de escoamento 142
APÊNDICE A – RESULTADOS PARA DIVERSAS CONDIÇÕES DE
ESCOAMENTO
Esse Apêndice apresenta os resultados médios e distribuições obtidas em
diversas configurações experimentais que não foram apresentadas nos capítulos de
resultados. De forma geral, os resultados aqui apresentados são semelhantes aos
do capítulo de resultados.
Ar e água na horizontal
0 200 400 600L/D
0
20
40
60
L B/D
ExperimentalSlug tracking
0 200 400 600L/D
0
10
20
30L S
/D
0 200 400 600L/D
0.0
0.5
1.0
1.5
UT [
m/s
]
0 200 400 600L/D
0
30
60
90
120
P [k
Pa]
Figura A. 1 – Resultados médios para o escoamento horizontal de ar e água na
condição A@W#1
Apêndice A – Resultados para diversas condições de escoamento 143
0 200 400 600L/D
0
40
80
120
160
L B/D
ExperimentalSlug tracking
0 200 400 600L/D
0
5
10
15
20
L S/D
0 200 400 600L/D
0.0
1.0
2.0
3.0
UT [
m/s
]
0 200 400 600L/D
0
35
70
105
140
P [k
Pa]
Figura A. 2 – Resultados médios para o escoamento horizontal de ar e água na
condição A@W#3
0 200 400 600L/D
0
5
10
15
20
25
L B/D
ExperimentalSlug tracking
0 200 400 600L/D
0
5
10
15
20
L S/D
0 200 400 600L/D
0.0
0.5
1.0
1.5
2.0
UT [
m/s
]
0 200 400 600L/D
0
30
60
90
120
P [k
Pa]
Figura A. 3 – Resultados médios para o escoamento horizontal de ar e água na
condição A@W#5
Apêndice A – Resultados para diversas condições de escoamento 144
0 200 400 600L/D
0
20
40
60
L B/D
ExperimentalSlug tracking
0 200 400 600L/D
0
5
10
15
20
L S/D
0 200 400 600L/D
0.0
0.5
1.0
1.5
2.0
2.5
3.0
UT [
m/s
]
0 200 400 600L/D
0
35
70
105
140
P [k
Pa]
Figura A. 4 – Resultados médios para o escoamento horizontal de ar e água na
condição A@W#6
Apêndice A – Resultados para diversas condições de escoamento 145
0.00
0.01
0.02
0.03
0.04
- Est
ação
3
0 20 40 60 80 100LB/D
0.00
0.01
0.02
0.03
0.04
- Est
ação
2
0.00
0.01
0.02
0.03
0.04PD
F - E
staç
ão 4 Experimental
Slug Tracking
0.00
0.03
0.06
0.09
- Est
ação
3
0 10 20 30 40LS/D
0.00
0.03
0.06
0.09PD
F - E
staç
ão 2
0.00
0.03
0.06
0.09
- Est
ação
40
4
8
12
- Est
ação
3
0 0.5 1 1.5 2UT (m/s)
0
4
8
12
- Est
ação
2
0
4
8
12
- Est
ação
4
Figura A. 5 - Resultado das funções densidade de probabilidade para o escoamento horizontal de ar e água na condição AW#1
Apêndice A – Resultados para diversas condições de escoamento 146
0.000
0.005
0.010
0.015
- Est
ação
3
0 100 200 300LB/D
0.000
0.005
0.010
0.015
- Est
ação
2
0.000
0.005
0.010
0.015
- Est
ação
4
ExperimentalSlug Tracking
0.000
0.025
0.050
0.075
0.100
- Est
ação
3
0 20 40 60 80LS/D
0.000
0.025
0.050
0.075
0.100PD
F - E
staç
ão 2
0.000
0.025
0.050
0.075
0.100
- Est
ação
40
1
2
3
- Est
ação
3
0 2 4 6 8 10UT (m/s)
0
1
2
3
- Est
ação
2
0
1
2
3
- Est
ação
4
Figura A. 6 – Resultado das funções densidade de probabilidade para o escoamento horizontal de ar e água na condição
A@W#3
Apêndice A – Resultados para diversas condições de escoamento 147
0.00
0.02
0.04
0.06
0.08
- Est
ação
3
0 10 20 30 40 50LB/D
0.00
0.02
0.04
0.06
0.08
- Est
ação
2
0.00
0.02
0.04
0.06
0.08
- Est
ação
4
ExperimentalSlug Tracking
0.00
0.05
0.10
0.15
- Est
ação
3
0 10 20 30LS/D
0.00
0.05
0.10
0.15PD
F - E
staç
ão 2
0.00
0.05
0.10
0.15
- Est
ação
40
2
4
6
8
- Est
ação
3
0 0.5 1 1.5 2UT (m/s)
0
2
4
6
8
- Est
ação
2
0
2
4
6
8
- Est
ação
4
Figura A. 7 – Resultado das funções densidade de probabilidade para o escoamento horizontal de ar e água na condição
A@W#5
Apêndice A – Resultados para diversas condições de escoamento 148
0.00
0.01
0.02
0.03
0.04
- Est
ação
3
0 20 40 60 80 100LB/D
0.00
0.01
0.02
0.03
0.04
- Est
ação
2
0.00
0.01
0.02
0.03
0.04
- Est
ação
4
ExperimentalSlug Tracking
0.00
0.04
0.08
0.12
- Est
ação
3
0 10 20 30 40LS/D
0.00
0.04
0.08
0.12PD
F - E
staç
ão 2
0.00
0.04
0.08
0.12
- Est
ação
40
1
2
3
4
- Est
ação
3
0 1 2 3UT (m/s)
0
1
2
3
4
- Est
ação
2
0
1
2
3
4
- Est
ação
4
Figura A. 8 – Resultado das funções densidade de probabilidade para o escoamento horizontal de ar e água na condição
A@W#6
Apêndice A – Resultados para diversas condições de escoamento 149
Ar e glicerina na horizontal
0 200 400 600L/D
0
10
20
30
40
L B/D
ExperimentalSlug tracking
0 200 400 600L/D
0
5
10
15
20
L S/D
0 200 400 600L/D
0.0
0.5
1.0
1.5
2.0
UT [
m/s
]
0 200 400 600L/D
0
35
70
105
140P
[kPa
]
Figura A. 9 – Resultados médios para o escoamento horizontal de ar e solução de
glicerina na condição A@G#1
Apêndice A – Resultados para diversas condições de escoamento 150
0 200 400 600L/D
0
50
100
150
200
L B/D
ExperimentalSlug tracking
0 200 400 600L/D
0
5
10
15
20
L S/D
0 200 400 600L/D
0.0
1.0
2.0
3.0
4.0
UT [
m/s
]
0 200 400 600L/D
0
35
70
105
140
P [k
Pa]
Figura A. 10 – Resultados médios para o escoamento horizontal de ar e solução de
glicerina na condição A@G#3
0 200 400 600L/D
0
4
8
12
16
20
L B/D
ExperimentalSlug tracking
0 200 400 600L/D
0
5
10
15
L S/D
0 200 400 600L/D
0.0
0.5
1.0
1.5
2.0
2.5
UT [
m/s
]
0 200 400 600L/D
0
35
70
105
140
P [k
Pa]
Figura A. 11 – Resultados médios para o escoamento horizontal de ar e solução de
glicerina na condição A@G#5
Apêndice A – Resultados para diversas condições de escoamento 151
0 200 400 600L/D
0
10
20
30
40
L B/D
ExperimentalSlug tracking
0 200 400 600L/D
0
5
10
15
L S/D
0 200 400 600L/D
0.0
0.5
1.0
1.5
2.0
2.5
3.0
UT [
m/s
]
0 200 400 600L/D
0
35
70
105
140
P [k
Pa]
Figura A. 12 – Resultados médios para o escoamento horizontal de ar e solução de
glicerina na condição A@G#6
Apêndice A – Resultados para diversas condições de escoamento 152
0.00
0.01
0.02
0.03
0.04
- Est
ação
3
0 20 40 60 80LB/D
0.00
0.01
0.02
0.03
0.04
- Est
ação
2
0.00
0.01
0.02
0.03
0.04PD
F - E
staç
ão 4 Experimental
Slug Tracking
0.00
0.05
0.10
0.15
0.20
- Est
ação
3
0 5 10 15 20 25LS/D
0.00
0.05
0.10
0.15
0.20
- Est
ação
2
0.00
0.05
0.10
0.15
0.20
- Est
ação
40
2
4
6
8
10
- Est
ação
3
0 0.5 1 1.5UT (m/s)
0
4
8
12
- Est
ação
2
0
2
4
6
8
- Est
ação
4
Figura A. 13 – Resultado das funções densidade de probabilidade para o escoamento horizontal de ar e solução de glicerina na
condição A@G#1
Apêndice A – Resultados para diversas condições de escoamento 153
0.000
0.005
0.010
0.015
0.020
- Est
ação
3
0 100 200 300LB/D
0.000
0.005
0.010
0.015
0.020
- Est
ação
2
0.000
0.005
0.010
0.015
0.020
- Est
ação
4 ExperimentalSlug Tracking
0.00
0.05
0.10
0.15
- Est
ação
3
0 10 20 30 40 50LS/D
0.00
0.05
0.10
0.15PD
F - E
staç
ão 2
0.00
0.05
0.10
0.15
- Est
ação
40
1
2
3
4
- Est
ação
3
0 2 4 6 8UT (m/s)
0
1
2
3
4
- Est
ação
2
0
1
2
3
4
- Est
ação
4
Figura A. 14 – Resultado das funções densidade de probabilidade para o escoamento horizontal de ar e solução de glicerina na
condição A@G#3
Apêndice A – Resultados para diversas condições de escoamento 154
0.00
0.05
0.10
0.15
- Est
ação
3
0 10 20 30LB/D
0.00
0.05
0.10
0.15
- Est
ação
2
0.00
0.05
0.10
0.15
- Est
ação
4 ExperimentalSlug Tracking
0.00
0.05
0.10
0.15
0.20
0.25
- Est
ação
3
0 5 10 15LS/D
0.00
0.05
0.10
0.15
0.20
0.25PD
F - E
staç
ão 2
0.00
0.05
0.10
0.15
0.20
0.25
- Est
ação
40
2
4
6
- Est
ação
3
0 0.5 1 1.5 2UT (m/s)
0
2
4
6
- Est
ação
2
0
2
4
6
- Est
ação
4
Figura A. 15 – Resultado das funções densidade de probabilidade para o escoamento horizontal de ar e solução de glicerina na
condição A@G#5
Apêndice A – Resultados para diversas condições de escoamento 155
0.00
0.02
0.04
0.06
- Est
ação
3
0 20 40 60LB/D
0.00
0.02
0.04
0.06
- Est
ação
2
0.00
0.02
0.04
0.06
- Est
ação
4 ExperimentalSlug Tracking
0.0
0.1
0.2
0.3
- Est
ação
3
0 5 10 15 20LS/D
0.0
0.1
0.2
0.3PD
F - E
staç
ão 2
0.0
0.1
0.2
0.3
- Est
ação
40
1
2
3
- Est
ação
3
0 1 2 3UT (m/s)
0
1
2
3
- Est
ação
2
0
1
2
3
- Est
ação
4
Figura A. 16 – Resultado das funções densidade de probabilidade para o escoamento horizontal de ar e solução de glicerina na
condição A@G#6
Apêndice A – Resultados para diversas condições de escoamento 156
Nitrogênio e óleo SAE 20-50 na horizontal
0 200 400L/D
0
20
40
60
L B/D
ExperimentalSlug tracking
0 200 400L/D
0
5
10
15
20
L S/D
0 200 400L/D
0.0
1.0
2.0
3.0
4.0
UT [
m/s
]
0 200 400L/D
0
60
120
180P
[kPa
]
Figura A. 17 – Resultados médios para o escoamento horizontal de nitrogênio e óleo
SAE 20-50 N@O#6
Apêndice A – Resultados para diversas condições de escoamento 157
0 200 400L/D
0
4
8
12
16
20
L B/D
ExperimentalSlug tracking
0 200 400L/D
0
5
10
15
20
L S/D
0 200 400L/D
0.0
1.0
2.0
3.0
4.0
UT [
m/s
]
0 200 400L/D
0
60
120
180
240
P [k
Pa]
Figura A. 18 – Resultados médios para o escoamento horizontal de nitrogênio e óleo
SAE 20-50 N@O#9
Apêndice A – Resultados para diversas condições de escoamento 158
0.00
0.01
0.02
0.03
0.04PD
F - E
staç
ão 3
0 50 100 150 200LB/D
0.00
0.01
0.02
0.03
0.04
- Est
ação
2ExperimentalSlug Tracking
0.00
0.02
0.04
0.06
0.08
0.10
- Est
ação
30 10 20 30 40 50
LS/D0.00
0.02
0.04
0.06
0.08
0.10
- Est
ação
2
0
1
2
3
4
- Est
ação
3
0 2 4 6UT (m/s)
0
1
2
3
4
- Est
ação
2
Figura A. 19 – Resultado das funções densidade de probabilidade para o escoamento horizontal de nitrogênio e óleo SAE 20-
50 na condição N@O#6
Apêndice A – Resultados para diversas condições de escoamento 159
0.00
0.04
0.08
0.12
0.16
- Est
ação
3
0 10 20 30 40LB/D
0.00
0.04
0.08
0.12
0.16
- Est
ação
2ExperimentalSlug Tracking
0.00
0.02
0.04
0.06
0.08
0.10
- Est
ação
30 10 20 30 40
LS/D0.00
0.02
0.04
0.06
0.08
0.10
- Est
ação
2
0
1
2
3
4
5
- Est
ação
3
0 0.5 1 1.5 2 2.5 3 3.5UT (m/s)
0
1
2
3
4
5
- Est
ação
2
Figura A. 20 – Resultado das funções densidade de probabilidade para o escoamento horizontal de nitrogênio e óleo SAE 20-
50 na condição N@O#9
Apêndice A – Resultados para diversas condições de escoamento 160
Ar e água na vertical
0 100 200L/D
0
20
40
60
L B/D
ExperimentalSlug tracking
0 100 200L/D
0
5
10
15
20
25
L S/D
0 100 200L/D
0.00.51.01.52.02.53.03.54.0
UT [
m/s
]
0 100 200L/D
0
50
100
150P
[kPa
]
Figura A. 21 – Resultados médios para o escoamento vertical de ar e água na
condição A@WV#5
Apêndice A – Resultados para diversas condições de escoamento 161
0 100 200L/D
0
5
10
15
L B/D
ExperimentalSlug tracking
0 100 200L/D
0
5
10
15
20
L S/D
0 100 200L/D
0.0
1.0
2.0
3.0
UT [
m/s
]
0 100 200L/D
0
40
80
120
160
P [k
Pa]
Figura A. 22 – Resultados médios para o escoamento vertical de ar e água na
condição A@WV#10
Apêndice A – Resultados para diversas condições de escoamento 162
0.00
0.01
0.02
0.03
0.04
0.05PD
F - E
staç
ão 2
0 40 80 120LB/D
0.00
0.01
0.02
0.03
0.04
0.05
- Est
ação
1ExperimentalSlug Tracking
0.00
0.04
0.08
0.12
- Est
ação
20 10 20 30 40
LS/D0.00
0.04
0.08
0.12
- Est
ação
1
0.0
0.5
1.0
1.5
- Est
ação
2
0 1 2 3 4 5UT (m/s)
0.0
0.5
1.0
1.5
- Est
ação
1
Figura A. 23 – Resultado das funções densidade de probabilidade para o escoamento vertical de ar e água na condição
A@WV#5
Apêndice A – Resultados para diversas condições de escoamento 163
0.0
0.1
0.2
0.3
0.4
- Est
ação
2
0 10 20 30LB/D
0.0
0.1
0.2
0.3
0.4
- Est
ação
1ExperimentalSlug Tracking
0.00
0.04
0.08
0.12
0.16
- Est
ação
20 10 20 30 40
LS/D0.00
0.04
0.08
0.12
0.16
- Est
ação
1
0
1
2
3
4
5
- Est
ação
2
0 1 2 3UT (m/s)
0
1
2
3
4
5
- Est
ação
1
Figura A. 24 – Resultado das funções densidade de probabilidade para o escoamento vertical de ar e água na condição
A@WV#10
Apêndice A – Resultados para diversas condições de escoamento 164
Nitrogênio e óleo SAE 20-50 na vertical
0 50 100 150L/D
0
5
10
15
20
25
30
L B/D
ExperimentalSlug tracking
0 50 100 150L/D
0
2.5
5
7.5
10
L S/D
0 50 100 150L/D
0.0
1.0
2.0
3.0
UT [
m/s
]
0 50 100 150L/D
0
50
100
150
200
250P
[kPa
]
Figura A. 25 – Resultados médios para o escoamento vertical de nitrogênio e óleo
SAE 20-50 na condição N@OV#8
Apêndice A – Resultados para diversas condições de escoamento 165
0.00
0.02
0.04
0.06
0.08PD
F - E
staç
ão 2
0 20 40 60LB/D
0.00
0.02
0.04
0.06
0.08
- Est
ação
1ExperimentalSlug Tracking
0.0
0.1
0.2
0.3
- Est
ação
20 5 10 15
LS/D0.0
0.1
0.2
0.3
- Est
ação
1
0.0
0.4
0.8
1.2
- Est
ação
2
0 1 2 3 4UT (m/s)
0.0
0.4
0.8
1.2
- Est
ação
1
Figura A. 26 – Resultado das funções densidade de probabilidade para o escoamento vertical de nitrogênio e óleo SAE 20-50
na condição N@OV#8
Apêndice A – Resultados para diversas condições de escoamento 166
Ar e água inclinado
0 100 200L/D
0
2.5
5
7.5
10L B
/D
ExperimentalSlug tracking
0 100 200L/D
0
5
10
15
20
L S/D
0 100 200L/D
0.0
1.0
2.0
3.0
UT [
m/s
]
0 100 200L/D
0
40
80
120
160
P [k
Pa]
Figura A. 27 – Resultados médios para o escoamento inclinado de ar e água na
condição A@WI#8
Apêndice A – Resultados para diversas condições de escoamento 167
0 100 200L/D
0
10
20
30
L B/D
ExperimentalSlug tracking
0 100 200L/D
0
5
10
15
20
L S/D
0 100 200L/D
0.0
1.0
2.0
3.0
4.0
5.0
UT [
m/s
]
0 100 200L/D
0
40
80
120
160
P [k
Pa]
Figura A. 28 – Resultados médios para o escoamento inclinado de ar e água na
condição A@WI#10
Apêndice A – Resultados para diversas condições de escoamento 168
0.0
0.4
0.8
1.2PD
F - E
staç
ão 2
0 5 10 15LB/D
0.0
0.4
0.8
1.2
- Est
ação
1ExperimentalSlug Tracking
0.00
0.04
0.08
0.12
0.16
- Est
ação
20 10 20 30 40 50
LS/D0.00
0.04
0.08
0.12
0.16
- Est
ação
1
0
1
2
3
4
5
- Est
ação
2
0 2 4 6UT (m/s)
0
1
2
3
4
5
- Est
ação
1
Figura A. 29 – Resultado das funções densidade de probabilidade para o escoamento inclinado de ar e água na condição
A@WI#8
Apêndice A – Resultados para diversas condições de escoamento 169
0.00
0.05
0.10
0.15
- Est
ação
2
0 10 20 30 40 50LB/D
0.00
0.05
0.10
0.15
- Est
ação
1ExperimentalSlug Tracking
0.00
0.04
0.08
0.12
- Est
ação
20 10 20 30 40
LS/D0.00
0.04
0.08
0.12
- Est
ação
1
0.00
0.25
0.50
0.75
1.00
- Est
ação
2
0 2 4 6UT (m/s)
0.00
0.25
0.50
0.75
1.00
- Est
ação
1
Figura A. 30 – Resultado das funções densidade de probabilidade para o escoamento inclinado de ar e água na condição
A@WI#10
Apêndice B - Resumo de correlações para cálculo das variáveis de fechamento
170
APÊNDICE B – RESUMO DE CORRELAÇÕES PARA CÁLCULO DAS
VARIÁVEIS DE FECHAMENTO
Esse apêndice está conectado ao exposto no capítulo 2, de revisão
bibliográfica. Aqui é apresentado um resumo das correlações existentes na literatura
para o cálculo das equações de fechamento para a velocidade de translação das
bolhas, fração de líquido no pistão e freqüência da célula unitária. As Tabelas B.1,
B.2 e B.3 apresentam o resumo das correlações.
Apêndice B - Resumo de correlações para cálculo das variáveis de fechamento
171
Tabela B. 1 – Valores das constantes para o cálculo da velocidade da bolha alongada
Autores 0C /DU gD h Observações
Davies e Taylor (1950) - 0,328 - Experimental, vertical.
Nicklin et al. (1962) 1,2 0,351 - Experimental, vertical.
Moissis e Griffith (1962) - - ( )8exp 1,06 SL D− Experimental, vertical.
Benjamin (1968) - 0,542 - Analítico, horizontal.
Dukler e Hubbard (1975) ( )1,022 0,021ln ReM+ 0 - Analítico, horizontal,
ReM JDρ μ= .
Ferré (1979) 1,10; 2, 261,30;2, 26 8, 281,02; 8, 28
M
M
M
FrFr
Fr
≤⎧⎪ < <⎨⎪ ≥⎩
0, 44; 2, 260;2,26 8, 283; 8, 28
M
M
M
FrFr
Fr
≤⎧⎪ ≤ <⎨⎪ ≥⎩
- Experimental, horizontal,
MFr J gD= .
Weber (1981) - 0,56
1,760,54Eo
− - Experimental, vertical,
2Eo gDρ σ= Δ .
Bendiksen (1984) ( )21,05 0,15 sen ; 3,51,2; 3,5
L
L
FrFr
θ⎧ + <⎪⎨
≥⎪⎩
0,54cos 0,35sen ; 3,50,35sen ; 3,5
L
L
FrFr
θ θθ
+ <⎧⎨ ≥⎩
- Experimental, diversas
inclinações,
L LFr J gD= .
Dukler et al. (1985) 1,225 0 - Analítico, horizontal e
vertical.
Théron (1989) ( )20, 231,3 0,13 senβ− +
Γ 0,80,5 cos 0,35senβ β⎛ ⎞− + +⎜ ⎟
⎝ Γ ⎠
- Experimental, diversas
inclinações, 3,5crítFr = ,
( )101 cos M crítFr FrθΓ = +.
Apêndice B - Resumo de correlações para cálculo das variáveis de fechamento
172
Taitel e Barnea (1993) - -
Stab
5,5exp 6 SLL
⎛ ⎞−⎜ ⎟
⎝ ⎠
Experimental, vertical,
10 15stabL D = ∼
Manolis (1995) 1,033; 2,861, 216; 2,86
M
M
FrFr
<⎧⎨ ≥⎩
0, 477; 2,860; 2,86
M
M
FrFr
<⎧⎨ ≥⎩
- Experimental, horizontal.
Woods e Hanratty (1996) 1,1; 3,11, 2; 3,1
M
M
FrFr
<⎧⎨ ≥⎩
0,52; 3,10; 3,1
M
M
FrFr
<⎧⎨ ≥⎩
- Experimental, horizontal.
Grenier (1997) - - ( )0, 4exp 0,5 SL D− Experimental, horizontal.
Petalas e Aziz (1998) ( ) 0,0311,64 0,12sen ReMβ −+ - - Experimental, horizontal e
vertical.
Fagundes Netto (1999) - - 0
crít
1 exp 0,16S S
S
L LL D
υ⎛ ⎞ ⎛ ⎞− −⎜ ⎟ ⎜ ⎟
⎝ ⎠⎝ ⎠
Experimental, horizontal,
0 0, 22υ = , crít 6,3SL D= .
Apêndice B - Resumo de correlações para cálculo das variáveis de fechamento
173
Tabela B. 2 – Modelos para a fração de líquido no pistão
Autores LSR Observações
Gregory et al. (1978) ( )
11,391 8,66J−
⎡ ⎤+⎣ ⎦ Experimental, horizontal.
Malnes (1982) ( ) 11 0,251 1 83 M LFr Bo−− −− + Experimental, horizontal, MFr J gD= , 2
L LBo gD ρ σ= .
Fershneider (1983) ( )
2221 Lj ABo gDβ ρ ρ
−−−⎡ ⎤+ Δ⎢ ⎥⎣ ⎦
Experimental, horizontal, 25A = , 0,1β = .
Barnea e Brauner
(1985) ( )
20,1 1,2ˆ1 0,058 0,605 Re 0,725S MBo Fr⎡ ⎤− −⎣ ⎦ Analítico, vertical.
Andreussi (1993) ( ) ( )0 1 1MF F Fr F+ + Experimental, horizontal, ( )( )20 máx 0;2,6 1 2 0,025F D⎡ ⎤= −⎣ ⎦ ,
( )3
41 2400 1 sen 3F Boθ −= − .
Marcano et al. (1998) ( ) 121,001 0,0179 0,0011J J−
+ + Experimental, horizontal.
Gomez et al. (2000) ( )6exp 0,45 2,48 10 ReSθ −⎡ ⎤− + ⋅⎣ ⎦ Experimental, diversas inclinações, 0 1,57 radianosθ< < .
Abdul Majeed (2000) ( )1,009 C J A− ⋅ Experimental, horizontal e levemente inclinado, 0,006 1,3377 G LC μ μ= + ,
1, se 01 sen , se 0
AA
θθ θ
= ≤⎧⎨ = − >⎩
.
Zhang et al. (2003)
( )
1
13,16
sm
L G
Tgρ ρ σ
−⎡ ⎤⎢ ⎥+⎢ ⎥−⎣ ⎦
Analítico, diversas inclinações,
( )( )212 4
L LF T F Fssm S
e S
u u j uf dT jC L
ρ αρ
− −⎡ ⎤= +⎢ ⎥
⎣ ⎦.
Apêndice B - Resumo de correlações para cálculo das variáveis de fechamento
174
Tabela B. 3 – Modelos para a freqüência da célula
Autores f Observações
Gregory e Scott (1969) ( ) 1,20,0266 2,02 D Frλ +⎡ ⎤⎣ ⎦ Experimental, horizontal, LJ Jλ = , Fr J gD= .
Heywood e Richardson
(1979) ( ) 1,02
0,0434 2,02 D Frλ +⎡ ⎤⎣ ⎦ Experimental, horizontal.
Correlação Shell,
apud Zabaras (2000) ( ) ( )min
20,1 0,0641,17
sl sg slFr Fr Fr Frg N A N N ND
⎧ ⎫⎡ ⎤+ + −⎨ ⎬⎢ ⎥⎣ ⎦⎩ ⎭
Experimental, horizontal,
e sl sgFr L Fr GN J gD N J gD= = ,
( ) ( )min
0,81 2,340,048 e 0,73
sl slFr Fr FrN N A N= = .
Zabaras (2000) ( ) ( ) 1,20,250,0226 0,836 2,75sen 19,75LJ gD J Jθ+ +⎡ ⎤⎣ ⎦
Experimental, diversas inclinações.
Sakaguchi et al. (2001) 10,317 1,61 0,333 3,040,5641,38 0,166 2
16100 0,087G G GL L L
L L L
J J JD J DD JJ J J gD
μ ρρ ρμ σ μ ρ
−− −−⎡ ⎤⎛ ⎞ ⎛ ⎞ ⎛ ⎞ ⎛ ⎞⎛ ⎞⎛ ⎞ ⎛ ⎞⎛ ⎞⎢ ⎥+⎜ ⎟ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟⎜ ⎟⎜ ⎟ ⎜ ⎟⎜ ⎟ ⎜ ⎟⎢ ⎥⎝ ⎠ ⎝ ⎠⎝ ⎠ ⎝ ⎠⎝ ⎠ ⎝ ⎠ ⎝ ⎠⎝ ⎠⎣ ⎦, experimental,
vertical.
Taitel e Dukler (1977) Método numérico transiente baseado na discretização espacial e temporal das equações de conservação da quantidade de
movimento no líquido e gás do escoamento estratificado. No processo transiente é calculado o tempo necessário para que uma
onda do escoamento estratificado encoste-se à parede superior do tubo, formando um pistão. Segundo os autores, o tempo
decorrido na solução numérica é proporcional à freqüência do escoamento.
Tronconi (1990) 0,61 G G
L G
uh
ρρ
Analítico, horizontal, Gu e Gh são a velocidade e altura do gás no
escoamento estratificado em equilíbrio.
Apêndice B - Resumo de correlações para cálculo das variáveis de fechamento
175
Hill e Wood (1990) ( )2.680.275 103600
LRL GJ JD
⋅+⎛ ⎞ ⋅⎜ ⎟⎝ ⎠
Experimental, horizontal, LR é a fração de líquido no escoamento
estratificado em equilíbrio.
Hill e Wood (1994) ( )( )
9.91209 ' 0.20524 '
0.3
24.729 0.00766 24.721 13600 1 0.05
H HL G
G
e e J JD J D
− + + +⎛ ⎞⎜ ⎟ −⎝ ⎠
Experimental, horizontal, ( )( )' 1 1 0.068e LH Jα= − − e LR é a fração de líquido no escoamento estratificado em equilíbrio.
Apêndice C - Modelos de Dois Fluidos e Drift Flux
176
APÊNDICE C – MODELOS DE DOIS FLUIDOS E DRIFT FLUX
Modelo de dois fluidos (two-fluid model)
O modelo de Dois Fluidos foi apresentado inicialmente por Ishii (1975) e será
apresentado aqui de forma sucinta, como mostrado em Shoham (2006). O modelo
trata o líquido e o gás como se escoassem de forma separada na tubulação, de
acordo com a Figura C.1. A maior hipótese utilizada é a de escoamento
unidimensional, fazendo-se as velocidades no líquido e no gás possuir um perfil
uniforme.
UL
UG
dz
z
AL
AG
Figura C. 1 – Esquema do escoamento para o modelo de Dois Fluidos
No modelo, aplicam-se as equações de conservação da massa e quantidade
de movimento a um volume de controle diferencial na direção axial, de acordo com a
Figura C.1. As equações de conservação da massa para as fases de líquido e gás
são:
( ) ( )L L L L L LA A U At z
ρ ρ∂ ∂+ = Γ
∂ ∂ (C.1)
e:
( ) ( )G G G G G GA A U At z
ρ ρ∂ ∂+ = Γ
∂ ∂, (C.2)
Onde ρL e ρG são as massas específicas de líquido e gás, respectivamente, AL e AG
são as seções transversais preenchidas com líquido e gás, respectivamente, UL e UG
Apêndice C - Modelos de Dois Fluidos e Drift Flux
177
são as velocidades médias de líquido e gás, respectivamente, e LΓ e GΓ
representam a transferência de massa entre as fases de líquido e gás por unidade
de comprimento da tubulação. Por definição 0LΓ > quando líquido é gerado, e,
devido à conservação da massa, L GΓ = −Γ . As equações de conservação da
quantidade de movimento para o líquido e gás são:
( ) ( )
( )
2
ˆ
sen
L L LL L LL
L L LL L I I L L IL
A UA UAU
t zP A AS S A g P
z z
ρρ
τ τ ρ θ
∂∂+ − Γ =
∂ ∂∂ ∂
= − + − − +∂ ∂
, (C.3)
e:
( ) ( )
( )
2
ˆ
sen
G G GG G GG
G G GG G I I G G IG
A UA UAU
t zP A AS S A g P
z z
ρρ
τ τ ρ θ
∂∂+ − Γ =
∂ ∂∂ ∂
= − − − − +∂ ∂
, (C.4)
onde Lτ e Gτ são as tensões de cisalhamento na parede nas fases de líquido e gás,
e Iτ é a tensão na interface, LS , GS e IS são os perímetros do líquido, gás e
interface, ˆGU U= se 0LΓ > e ˆ
LU U= se 0LΓ < , e ILP e IGP são as pressões do
líquido e do gás na interface.
Dessa forma, o modelo de Dois Fluidos possui 4 equações e 4 incógnitas: LU ,
GU , LA (ou GA ) e IP . O domínio computacional é discretizado e as Equações (C.1),
(C.2), (C.3) e (C.4) são resolvidas numericamente utilizando, geralmente, o método
de volumes finitos.
Bendiksen et al. (1991) apresentam a implementação do modelo de Dois
Fluidos no programa computacional OLGA, que foi um dos primeiros programas para
cálculo de parâmetros transientes no escoamento bifásico, e até hoje é largamente
utilizado na indústria petrolífera. O modelo apresentado por Bendiksen et al. (1991)
apresenta uma pequena particularidade para os escoamentos estratificados e
anulares que podem conter pequenas partículas de líquido dispersas na região do
Apêndice C - Modelos de Dois Fluidos e Drift Flux
178
gás. Nesse caso, uma nova equação de conservação da massa é aplicada apenas
às partículas, porém, a equação de conservação da quantidade de movimento
continua sendo aplicada à mistura gás/partículas de líquido.
Em seguida, Henau e Raithby (1995) apresentaram um modelo de Dois
Fluidos aplicado especificamente para o escoamento em golfadas. A partir dos
resultados de LU , GU e LA para cada posição axial do modelo de Dois Fluidos, os
autores resolvem um sub-modelo de escoamento em golfadas (baseado no modelo
estacionário de Dukler e Hubbard, 1975). Para isso, é necessária a utilização de
relações constitutivas para a tensão cisalhante na interface, coeficiente de arrasto
entre o líquido e o gás e para o coeficiente da força de massa virtual.
Mais recentemente, Issa e Kempf (2003) apresentaram resultados para um
modelo de Dois Fluidos aplicado à formação do escoamento em golfadas a partir de
um escoamento estratificado. Depois de formado, o escoamento em golfadas
continua a ser simulado, porém, sem a utilização de um sub-modelo como em
Henau e Raithby (1995). O modelo só é válido para escoamento horizontais ou
levemente inclinados e, segundo os autores, os resultados das simulações são bons
em comparação com resultados experimentais.
De forma geral, os modelos de Dois Fluidos apresentam alguns problemas
quanto à estabilidade e convergência dos resultados, principalmente se em algum
local da tubulação a fração do líquido ou do gás torna-se nula. Por isso, deve-se
realizar um tratamento cuidadoso das equações resultantes. Além disso, testes de
malha revelam que para uma solução independente da malha, esta deve ser muito
refinada, o que torna a simulação inviável para tubulações de grandes extensões.
Por fim, o modelo é muito sensível às equações constitutivas utilizadas.
Modelo de drift flux
O modelo de Drift Flux surgiu como alternativa ao modelo de Dois Fluidos, na
modelagem de escoamento em que a diferença entre as velocidades das fases é
pequena, característica de escoamentos dispersos, como representado na Figura
C.1. O modelo de Drift Flux gera uma equação a menos que o modelo de Dois
Fluidos, ou seja, 3 equações (se considerássemos o escoamento com transferência
Apêndice C - Modelos de Dois Fluidos e Drift Flux
179
de calor, o modelo geraria 4 equações, contra 6 do modelo de Dois Fluidos). A
seguir é apresentado o modelo de acordo com o apresentado em Shoham (2006).
G
dz
z
c
Figura C. 2 – Esquema do escoamento para o modelo de Drift Flux
A Figura C.1 apresenta o escoamento utilizado para a apresentação do
modelo, e G é o fluxo mássico total por unidade de área ( L L L G G GG R U R Uρ ρ= + ), e c
é a concentração de gás no escoamento ( G G Mc R ρ ρ= ). Aplicando-se a equação de
conservação da massa para a fase de gás no volume de controle infinitesimal da
Figura C.1, obtém-se:
( ) ( )G G G G G GA A U At z
ρ ρ∂ ∂+ = Γ
∂ ∂. (C.5)
É possível também aplicar a equação de conservação da massa para a mistura:
( ) ( ) 0M A AGt z
ρ∂ ∂+ =
∂ ∂, (C.6)
onde Mρ é a massa específica da mistura ( M G G L LR Rρ ρ ρ= + ). A Equação (C.5) pode
ser rearranjada em função de G e c como:
( )1 M MG G
M M M
AcUc G ct z A z
ρρ ρ ρ
∂ Γ∂ ∂+ + =
∂ ∂ ∂, (C.7)
onde c é a concentração mássica de gás na mistura e MGU é a velocidade difusiva
do gás, definida como a velocidade do gás relativa a uma superfície que se move no
centro de massa ( MG G MU U G ρ= − ).
Apêndice C - Modelos de Dois Fluidos e Drift Flux
180
A equação de conservação da quantidade de movimento também é escrita
para a mistura, da forma:
( ) ( )2 2
senL L L G G GL L G G M
A U A UAG PS S Ag At z z
ρ ρτ τ ρ θ
∂ +∂ ∂+ = − − − −
∂ ∂ ∂, (C.8)
que pode ser rearranjada em função de G e c como:
( ) 22 sen
1MG M L L G G MM
AG AG c PAU S S Ag At z z c z
ρ τ τ ρ θρ
∂ ⎛ ⎞∂ ∂ ∂⎛ ⎞+ + = − − − −⎜ ⎟ ⎜ ⎟∂ ∂ ∂ − ∂⎝ ⎠⎝ ⎠. (C.9)
O modelo consiste em se resolver as Equações (C.6), (C.7) e (C.9) ao longo do tubo,
em uma malha unidimensional, para as incógnitas G , c e P .
Segundo apresentado por Pauchon et al. (1994), o programa computacional
TACITE utiliza o modelo de Drift Flux para simular escoamentos bifásicos em
diversos padrões. Para isso, o programa possui um identificador do padrão de
escoamento e aplica equações de fechamento pertinentes ao padrão detectado.
Omgba et al. (2000) apresenta um modelo de Drift Flux que utiliza uma malha
adaptativa, que é refinada automaticamente de acordo com a necessidade do
modelo. Os autores não citam quais os modelos utilizados para cada padrão de
escoamento ou equações constitutivas.
Os modelos de Drift Flux possuem resolução semelhante aos de Dois Fluidos.
Por isso, apresentam as mesmas deficiências com relação à estabilidade e
refinamento da malha computacional.