ads:
UNIVERSIDADE DO ESTADO DO RIO DE JANEIRO
INSTITUTO POLITÉCNICO
Nelson Machado Barbosa
Resolução numérica de equações diferenciais parciais hiperbólicas não lineares: um
estudo visando a recuperação de petróleo
Nova Friburgo
2010
ads:
Livros Grátis
http://www.livrosgratis.com.br
Milhares de livros grátis para download.
Nelson Machado Barbosa
Resolução numérica de equações diferenciais parciais
hiperbólicas não lineares: um estudo visando a
recuperação de petróleo
Dissertação apresentada como requisito parcial
para a obtenção do título de Mestre, ao
Programa de Pós-Graduação em Modelagem
Computacional, do Instituto Politécnico, da
Universidade do Estado do Rio de Janeiro.
Orientador: Prof. Luiz Nélio Henderson Guedes de Oliveira
Nova Friburgo
2010
ads:
Nelson Machado Barbosa
Resolução numérica de equações diferenciais parciais
hiperbólicas não lineares: um estudo visando a
recuperação de petróleo
Dissertação apresentada como requisito parcial
para a obtenção do título de Mestre, ao
Programa de Pós-Graduação em Modelagem
Computacional, do Instituto Politécnico, da
Universidade do Estado do Rio de Janeiro.
Orientador: Prof. Luiz Nélio Henderson Guedes de Oliveira
Aprovado em 26 de fevereiro de 2010.
Banca Examinadora:
_______________________________________________
Luiz Nélio Henderson G. de Oliveira, D.Sc. - orientador
IPRJ/UERJ
_______________________________________________
Maurício Kischinhevsky, D.Sc.
UFF
_______________________________________________
Mikhail Vishnevsk Petrovich, PhD.
UENF
DEDICATÓRIA
A DEUS, aos meus
pais e
irmã e a minha namorada Rozeane.
AGRADECIMENTOS
Primeiramente agradeço ao meu Senhor Deus, por tudo de bom que me deste, desde o
ínicio da vida, pois sem ele nada consigo e nada sou.
Aos meus amados pais Ernane e Fátima, que sempre me motivaram e incentivavam ao
estudo, com certeza sempre serão lembrados com amor, pois foram à base fundamental para
formar o homem que sou.
A minha namorada e futura esposa Rozeane, que sempre esteve inseparadamente ao meu
lado nesses últimos seis anos acadêmicos, demonstrando dedicação e compreensão que foram
essenciais para minha formação.
Ao meu orientador Luiz Nélio Henderson o qual sempre esteve disposto a me receber,
agradeço também pela valiosa e indispensável orientação aonde sempre agiu com
competência demonstrando-se sempre claro, objetivo e incentivador, qualidades que foram
fundamentais para a realização desse trabalho.
A minha família que sempre esteve presente na minha vida, em especial a minha avó
Nilza pela sua serenidade sempre me acalmando em momentos difícieis, minha irmã Salete
que trouxe a vida minha afiliada Ana Beatriz alegrando mais ainda a minha vida e o meu
primo Bruce (em memória), que infelizmente nos deixou no ano de 2006, mas sempre será
lembrado com carinho e saudades pela sua forma de viver alegre e brincalhão.
À Faperj, que me confiou uma bolsa através do programa (Bolsa Faperj Nota 10), que foi
fundamental para as viagens em congressos, colóquios, etc.
Aos meus amigos e professores do IPRJ, em especial o pessoal da casa15 (Newton,
Felipe, Sandro, Wagner e Agnaldo) e os das casas vizinha (Juan, Julliany, Pojukan, Alex,
Thiago, etc) o pessoal da administração, principalmente da Biblioteca (Tereza, Rogério e
Juliana), guarnição e todos outros do IPRJ não mencionados.
Aos professores Maurício Kischinhevsky e Mikhail Vishnevsk Petrovich por aceitarem a
fazer parte da banca examinadora.
Ao IPRJ, pela infra-estrutura e organização que oferece aos estudantes, o qual não só
agradeço e também parabenizo a todos que zelam pela forma de gerenciar esse patrimônio.
E novamente agradeço a Deus por dispor-me de saúde, sabedoria e força de vontade,
pois sem nenhum desses requisitos não estaria obtendo esse título.
RESUMO
BARBOSA, Nelson Machado. Resolução numérica de equações diferenciais parciais
hiperbólicas não lineares: um estudo visando a recuperação de petróleo. 2010. 140 f.
Dissertação (Mestrado em Modelagem Computacional) – Instituto Politécnico,
Universidade do Estado do Rio de Janeiro, Nova Friburgo, 2010.
O processo de recuperação secundária de petróleo é comumente realizado com a injeção de
água no reservatório a fim de manter a pressão necessária para sua extração. Para que o
investimento seja viável, os gastos com a extração têm de ser menores do que o retorno
financeiro obtido com o petróleo. Para tanto, tornam-se extremamente importantes as
simulações dos processos de extração. Neste trabalho são estudados os problemas de Burgers
e de Buckley-Leverett visando o escoamento imiscível água-óleo em meios porosos, onde o
escoamento é incompressível e os efeitos difusivos (devido à pressão capilar) são
desprezados. Com o objetivo de incorporar conhecimento matemático mais avançado, para
em seguida utilizá-lo no entendimento do problema estudado, abordou-se com razoável
profundidade a teoria das leis de conservação. Foram consideradas soluções fracas que,
fisicamente, podem ser interpretadas como ondas de choque ou rarefações, então, para que
fossem distinguidas as fisicamente admissíveis, foi utilizado o princípio de entropia, nas suas
diversas formas. Inicialmente consideramos alguns exemplos clássicos de métodos numéricos
para uma lei de conservação escalar, os quais podem ser vistos como esquemas conservativos
de três pontos. Entre eles, o método de Lax-Friedrichs (LF) e o método de Lax-Wendroff
(LW). Em seguida, um esquema composto foi testado, o qual inclui na sua formulação os
métodos LF e LW (chamado de LWLF-4). Respeitando a condição CFL, foram obtidas
soluções numéricas de todos os problemas tratados aqui. Com o objetivo de validar tais
soluções, foram utilizadas soluções analíticas oriundas dos problemas de Burgers e Buckley-
Leverett. Também foi feita uma comparação com os métodos do tipo TVD’s com limitadores
de fluxo, obtendo resultado satisfatório. Vale à pena ressaltar que o esquema LWLF-4, pelo
que nos consta, nunca foi antes utilizado nas resoluções das equações de Burgers e Buckley-
Leverett.
Palavras-chave: Recuperação secundária de petróleo; Equações hiperbólicas não lineares;
Problemas de Burgers e Buckley-Leverett; Método composto LWLF-k.
ABSTRACT
The secondary recovery of petroleum is usually performed with injection of water through an
oil reservoir to keep the oil pressure for the exploration. In order to make the exploration
profitable, the extraction cost must be less than the financial return, which means that the
simulation of the exploration process is extremely relevant. In this work, the Burgers- and-
Buckley-Leverett problems are studied seeking a two-phase displacement in porous media.
The flow is considered incompressible and capillary effects are ignored. In order to analyze
the problem, it was necessary to use the theory of conservation law in a spatial variable. Weak
solutions, which can be understood as shock or rarefaction waves, are studied with the
entropy condition, so that only the physically correct solutions are considered. Some classical
numerical methods, which can be seen as conservative schemes of three points, are studied,
among them the Lax-Friedrichs (LF) and Lax-Wendroff (LW) methods. A composite scheme,
called LWLF-k, is tested using LF and LW methods, being respected the CFL condition, with
satisfactory results. In order to validate the numerical schemes, we consider analytical
solutions of the Burgers-and-Buckley-Leverett equations. Was also made a comparison with
TVD’s methods with flux limiters, obtaining satisfactory results. We emphasize that to the
best of our knowledge, the LWLF-4 scheme has never been used to solve the Buckley-
Leverett equation.
Keywords: Secondary recuperation process of petroleum; Nonlinear hyperbolic equations;
Burgers and Buckley-Leverett problems; LWLF-k Composite Scheme.
LISTA DE TABELAS
Tabela 1- Limitadores de fluxo...............................................................................................105
Tabela 2 -Valores da saturação no ponto de choque da comparação dos métodos para a razão
de viscosidade c=0.1...............................................................................................................130
Tabela 3 -Valores da saturação no ponto de choque da comparação dos métodos para a razão
de viscosidade c=0.1...............................................................................................................130
LISTA DE FIGURAS
Figura 1.1-Esboço de um Reservatório de Petróleo sob injeção de água.................................18
Figura 2.1- Experimento de Darcy............................................................................................22
Figura 2.2- Volume elementar
V
com fronteira
V
................................................................24
Figura 2.3- Gráfico típico de permeabilidade relativa água/óleo..............................................34
Figura 2.4- Gráfico da função de fluxo fracionária para o modelo de Stone em distintos
valores para a razão da viscosidade
ow
c
/
.......................................................................37
Figura 2.5-Gráfico da função de fluxo fracionária para o modelo de Brooks – Corey.............37
Figura 2.6- Distribuição da Saturação da água a partit da equação de Buckley-Leverett........41
Figura 2.7- Solução descontínua para a equação de Buckley-Leverett....................................41
Figura 2.8- Determinação geométrica da saturação no ponto de choque.................................43
Figura 2.9- Gráfico da função de fluxo para o problema de Burgers.......................................45
Figura 3.1- Inclinação da reta característica.............................................................................50
Figura 3.2- Retas características que não se cruzam.................................................................51
Figura 3.3- Um típico leque de rarefação..................................................................................52
Figura 3.4- Duas retas características que se cruzam................................................................53
Figura 3.5- Instante do primeiro choque...................................................................................56
Figura 3.6- Os gráficos das funções
xS
0
e
xS
0
..................................................................59
Figura 3.7- Características definidas pela equação (3.35)........................................................59
Figura 3.8: Divisão do domínio de interesse em três conjuntos disjuntos................................60
Figura 3.9- Solução em cada subconjunto considerado............................................................61
Figura 3.10- Curva de descontinuidade....................................................................................63
Figura 3.11- Esboço da solução fraca do problema do exemplo 3.1........................................66
Figura 3.12- Suporte da função teste
....................................................................................67
Figura 3.13- Curva fechada em
2
........................................................................................69
Figura 3.14- Retas características do problema do exemplo 3.2...............................................71
Figura 3.15- Uma solução contínua para o problema do exemplo 3.2.....................................71
Figura 3.16- Uma solução descontínua para o problema do exemplo 3.2 ...............................72
Figura 3.17- Uma típica descontinuidade admissível...............................................................74
Figura 3.18- Descontinuidades admissíveis..............................................................................77
Figura 3.19- Descontinuidades não admissíveis.......................................................................77
Figura 3.20- Função de fluxo de Burgers para o dado inicial 1................................................79
Figura 3.21- Retas características para o problema de Burgers para o referido dado inicial 1.80
Figura 3.22- Gráfico da solução analítica do problema de Burgers para o dado inicial 1........81
Figura 3.23- Retas características para o problema de Burgers para o referido dado inicial 2.83
Figura 3.24- Gráfico da solução analítica do problema de Burgers para o dado inicial 2........84
Figura 3.25- Gráfico da solução analítica do problema de Burgers para o dado inicial 3........85
Figura 3.26- Solução do problema de Riemann, f convexa e
RL
SS
..................................86
Figura 3.27- Solução do problema de Riemann, f convexa e
RL
SS ..................................86
Figura 3.28- Choque como solução do problema de Riemann, f sigmoidal..........................87
Figura 3.29- Rarefação como a solução do problema de Riemann, f sigmoidal....................88
Figura 3.30- Rarefação seguida de um choque, f sigmoidal..................................................88
Figura3.31- Um choque seguido de rarefação, f sigmoidal....................................................89
Figura 4.1- Grade de células centradas.....................................................................................93
Figura 4.2- Domínio de dependência para um esquema de três pontos....................................96
Figura 4.3- Representação geométrica da condição CFL.........................................................97
Figura 4.4- Ilustração do domínio de dependência para o esquema Lax-Friedrichs..............101
Figura 5.1- Gráfico da Função
Sf
para o problema de Burgers........................................112
Figura 5.2- Gráfico da Função
Sf
para o problema de Buckley-Leverett, c=0.1..............113
Figura 5.3- Gráfico da Função
Sf
para o problema de Buckley-Leverett, c=0.5..............113
Figura 5.4- Gráfico da Função
Sf
para o problema de Buckley-Leverett, c=1.................114
Figura 5.5- Comparação dos esquemas
LW
e
LF
para o problema de Buckley-Leverett para
o passo de tempo
st 5
e razão de viscosidade
5.0
c
.....................................................115
Figura 5.6- Comparação dos esquemas
LW
e
LF
para o problema de Burgers para o passo
de tempo
st 20
e razão de viscosidade
1
c
....................................................................115
Figura 5.7- Comparação dos esquemas
LW
,
LF
e
4
LWLF
para o problema de Buckley-
Levertt para o passo de tempo
st 5
no tempo fixo
st 360
e razão de viscosidade
1
c
........................................................................................................................................116
Figura 5.8- Comparação dos esquemas
LW
,
LF
e
4
LWLF
para o problema de Buckley-
Levertt para o passo de tempo
st 5
no tempo fixo
st 900
e razão de viscosidade
5.0
c
....................................................................................................................................117
Figura 5.9- Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 150 células com
0067.0
x
................................118
Figura 5.10- Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 350 células com 0029.0
x ................................118
Figura 5.11- Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 450 células com 002.0
x ..................................119
Figura 5.12- Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 400 células com
01.0
x
....................................120
Figura 5.13- Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 800 células com
005.0
x
..................................120
Figura 5.14- Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 1000 células com
004.0
x
................................121
Figura 5.15- Comparação do esquema composto 4
LWLF com a solução analítica do
Problema de Burgers para uma grade de 1200 células com
0033.0
x
..............................121
Figura 5.16- Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 300 células com
0067.0
x
................................122
Figura 5.17- Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 400 células com
005.0
x
..................................122
Figura 5.18- Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 600 células com
0033.0
x
................................123
Figura 5.19- Aplicação do método
4
LWLF
para o problema de Buckley-Leverett com 150
células, cmx 67,66
,
st 10
e razão de viscosidade igual a
1
c
..................................124
Figura 5.20- Aplicação do método
4
LWLF
para o problema de Buckley-Leverett com 150
células, cmx 67,66
,
st 10
e razão de viscosidade igual a
5.0
c
..............................124
Figura 5.21- Aplicação do método
4
LWLF
para o problema de Buckley-Leverett com 250
células,
cmx 40
,
st 7
e razão de viscosidade igual a
1
c
........................................125
Figura 5.22- Aplicação do método
4
LWLF
para o problema de Buckley-Leverett com 250
células,
cmx 40
,
st 7
e razão de viscosidade igual a
5.0
c
.....................................125
Figura 5.23- Aplicação do método
4
LWLF
para o problema de Buckley-Leverett com 250
células,
cmx 40
,
st 7
e razão de viscosidade igual a
1.0
c
.....................................126
Figura 5.24- Comparação do método
4
LWLF
com o limitador de fluxo Superbee, para
1
c
........................................................................................................................................127
Figura 5.25- Comparação do método
4
LWLF
com o limitador de fluxo Superbee, para
5.0
c
.....................................................................................................................................127
Figura 5.26- Comparação do método
4
LWLF
com o limitador de fluxo Superbee, para
1.0
c
.....................................................................................................................................128
Figura 5.27- Comparação do método 4
LWLF com o limitador de fluxo Van Leer, para
1
c
........................................................................................................................................128
Figura 5.28- Comparação do método 4
LWLF com o limitador de fluxo Van Leer, para
1.0
c
.....................................................................................................................................129
Figura 5.29- Comparação do método
4
LWLF
com o limitador de fluxo Minmod, para
1
c
........................................................................................................................................129
Figura 5.30- Comparação do método
4
LWLF
com o limitador de fluxo Minmod, para
1.0
c
.....................................................................................................................................130
Figura 6.1- Mudança de variável............................................................................................134
Figura 6.2- Função de fluxo fracionária da água (Modelo de Juanes e Patzek).....................136
Figura 6.3- Função de fluxo fracionária do gás (Modelo de Juanes e Patzek).......................137
LISTA DE SÍMBOLOS
Índice inferior
relativo ao fluido
, (
w
o
,
com
óleoo
e águaw
)
i
relativo a i-ésima célula na discretização espacial
D
indica que a variável está sendo adimensionalizada
Alfabeto latino
A área euclidiana
c fator de compressibilidade do fluido
ow,
c razão de viscosidade
dA
elemento infinitesimal de área
dV
elemento infinitesimal de volume
DS
derivada lagrangeana da saturação
SE
função de entropia
f função de fluxo fracionária do fluido
ow,
g
magnitude da aceleração devido à gravidade
k
permeabilidade efetiva
K condutividade hidráulica
FIS
K tensor de difusão física
r
k permeabilidade relativa da fase
ow,
x
permeabilidade absoluta da rocha
w
SJ função
J
de Leverett
l distância percorrida pelo o fluido no experimento de Darcy
F
L
abreviação do método de Lax-Friedrichs
W
L abreviação do método de Lax-Wendroff
kLWLF
abreviação do método composto
n
vetor normal
RL
N
número de Rapoport e Leas
p
pressão global
p pressão da fase
ow,
c
p pressão capilar
q
taxa de massa por unidade de volume injetado em V
r
parâmetro do limitador de fluxo
e
S
saturação efetiva
o
S saturação do fluido óleo
w
S ou
S
, saturação do fluido água
*
w
S valor da saturação da água no ponto de choque
t
variável tempo
TVDs
abreviação do método do tipo TVDs com limitadores de fluxo
u velocidade de Darcy para o fluido
ow,
T
u
velocidade total (
w
u
+
o
u
)
T
u
~
velocidade total dividido pela a porosidade do meio
U
variável utilizada na discretização numérica
*
s
v
velocidade do fluido água no ponto de choque
V
elemento de volume
V
fronteira do elemento de volume
x
variável espacial
Z
altura piezométrica do experimento de Darcy
Letras gregas
representa um incremento da função que o acompanha
h
diferença de altura manométrica
t
passo de tempo para os problemas hiperbólicos convectivos
x
largura da célula na malha, na direção x
curva arbitrária
mobilidade total
mobilidade para a fase
, (
w
o
,
com óleoo
e águaw
)
viscosidade
densidade por unidade de volume da fase
, (
w
o
,
com
óleoo
e
águaw
)
porosidade do meio
ângulo de contato da interface rocha e água
tensão superficial dos fluidos
valor constante no eixo da abscissa no domínio unidimensional
conjunto de funções teste de
2
r
limitador de fluxo
coeficiente de difusão maior que zero
k
operador definido por (k-1) etapa de
w
L
seguido de uma aplicação do método
F
L
função de fluxo de entropia
fluxo numérico
fluxo numérico entrópico
somatório
domínio do reservatório
SUMÁRIO
1. INTRODUÇÃO...................................................................................................................17
2. MODELAGEM MATEMÁTICA PARA O ESCOAMENTO EM MEIOS
POROSOS ...............................................................................................................................20
2.1 Modelagem matemática para o escoamento bifásico.....................................................20
2.1.1 Escoamento bifásico........................................................................................................20
2.1.2 Lei de Darcy.....................................................................................................................21
2.1.3 Equação de balanço de massa..........................................................................................24
2.1.4 Desenvolvimento matemático para uma lei de conservação............................................26
2.2 Equação de Buckley-Leverett..........................................................................................30
2.2.1 Dedução da equação.........................................................................................................31
2.2.2 Solução semi-analítica.....................................................................................................37
2.2.3 Localização da frente de água..........................................................................................38
2.2.4 Solução analítica de Welge para a equação de Buckley-Leverett....................................40
2.3 Equação de Burgers..........................................................................................................43
3. LEIS DE CONSERVAÇÃO...............................................................................................46
3.1 Formulação........................................................................................................................46
3.2 A formação de ondas de choque e ondas de rarefação..................................................48
3.3 Soluções descontínuas e a condição de Rankine-Hugoniot...........................................58
3.4 Perda de unicidade e a condição de entropia.................................................................69
3.5 Solução analítica para o problema de Burgers..............................................................78
3.5.1 Dado inicial 1...................................................................................................................78
3.5.2 Dado inicial 2...................................................................................................................81
3.5.3 Dado inicial 3...................................................................................................................84
3.6 Análise do problema de Riemann....................................................................................85
3.7 O teorema de existência e unicidade de Kruzkhov........................................................89
4. ABORDAGENS NUMÉRICAS PARA LEIS DE CONSERVAÇÃO ...........................93
4.1 Esquemas conservativos...................................................................................................93
4.2 A condição CFL.................................................................................................................95
4.3 Método de Lax-Friedrichs................................................................................................98
4.4 Método de Lax-Wendroff...............................................................................................101
4.5 Métodos TVDs com limitadores de fluxo......................................................................103
4.6 Teorema de Lax-Wendroff.............................................................................................105
4.7 A desigualdade de entropia de discreta.........................................................................106
4.8 Esquema composto..........................................................................................................108
5. RESULTADOS..................................................................................................................111
5.1 Aspectos gerais................................................................................................................111
5.2 Condições de injeção e extração....................................................................................114
5.3 Comparações dos métodos LF, LW e LFLW-4............................................................114
5.4 Aplicação do método LWLF-4 ao problema de Burgers.............................................117
5.4.1 Resultados para o dado inicial 1....................................................................................117
5.4.2 Resultados para o dado inicial 2....................................................................................119
5.4.3 Resultados para o dado inicial 3....................................................................................122
5.5 Aplicação do método LWLF-4 no problema de Buckley-Leverett.............................123
5.6 Comparações do método LWLF-4 com os métodos tipo TVDs com limitadores de
fluxo na captura do choque..................................................................................................126
6. CONCLUSÕES E PERSPECTIVAS..............................................................................131
6.1 Conclusões........................................................................................................................131
6.2 Perspectivas.....................................................................................................................132
REFERÊNCIAS....................................................................................................................138
17
1. INTRODUÇÃO
O petróleo é um recurso natural abundante nos dias atuais, considerado uma fonte de
energia não renovável. O petróleo bruto possui em sua composição uma cadeia de
hidrocarbonetos, cujas frações leves formam os gases e as frações pesadas o óleo cru. A
distribuição percentual desses hidrocarbonetos é o que define os diversos tipos de petróleo
existentes no mundo. A sua pesquisa envolve elevados custos e complexidades de estudos. No
entanto, ele continua a ser a principal fonte de energia, servindo como base para a fabricação
dos mais variados produtos, dentre os quais se destacam: benzinas, óleo diesel, gasolina,
alcatrão, polímeros e até medicamentos, dentre outros. Hoje é a principal renda de muitos
paises, incluindo os do Oriente Médio e o Brasil.
Em virtude dessa importância as companhias petrolíferas continuam buscando a sua
exploração, com a descoberta de novas jazidas e com a explotação adequada dos reservatórios
já descobertos, preocupadas sempre em manter ou aumentar suas reservas de hidrocarbonetos.
Para tais explorações tornam-se indispensável simulações, as quais são ferramentas cruciais
no gerenciamento de reservatórios, sendo importante em todas as fases, desde a descoberta do
campo até o abandono. A simulação de petróleo tem como meta a predição da produção de
modo que novas técnicas possam ser implementadas, a fim de otimizar a recuperação do
petróleo. Tal recuperação pode ser classificada como primária, secundária e terciária.
No estágio de recuperação primária a pressão no reservatório é muito alta, de modo
que o óleo ou o gás é produzido naturalmente devido à diferença de pressão, e o seu fim
ocorre quando há equilíbrio entre a pressão do óleo e a atmosférica. Ao final desse estágio,
geralmente, ainda existe de 70% a 85% de óleo no reservatório.
Na recuperação secundária explora-se a imiscibilidade dos fluidos água e óleo, ou seja,
água (por ser um fluido mais econômico) é injetada através de poços de injeção e o óleo é
recuperado nos poços de produção. Este processo serve para manter elevada a pressão do
reservatório, bem como para deslocar o óleo na direção do poço produtor. Sabe-se que esse
processo não é totalmente eficiente, pois ao final de tal estágio ainda existirá
aproximadamente 50% ou mais de óleo no reservatório. Essa quantidade se mantém devido à
tensão superficial, de fato uma grande quantidade de óleo ficará pressa nos poros menores e a
água não será capaz de removê-los. Além disso, como o óleo é (geralmente) muito viscoso e
18
pesado, a água terá uma grande mobilidade no interior do reservatório. Depois de certo tempo
os poços de produção recuperarão essencialmente água denotando o fim desse estágio.
Nesse trabalho será dada ênfase a essa etapa de recuperação, procurando solucionar os
problemas de Burgers e Buckley-Leverett. Os referidos problemas são clássicos em
engenharia de petróleo, representados por uma equação hiperbólica de transporte não linear,
onde objetivo é capturar o choque e a verificação do comportamento da saturação da água no
reservatório durante o processo de injeção.
Figura 1.1 - Esboço de um Reservatório de Petróleo sob injeção de água.
A recuperação terciária se baseia em novas técnicas envolvendo processos químicos e
térmicos. Nessa etapa são injetados diferentes materiais, normalmente não presentes nos
reservatórios com o objetivo de alcançar uma maior miscibilidade do óleo.
A simulação computacional é uma poderosa ferramenta, a qual é usada na simulação
de reservatórios de petróleo. Ao contrário de métodos puramente analíticos, a simulação
computacional possibilita o uso de um grande número de grau de liberdade no tratamento dos
fenômenos, bem como a incorporação de efeitos não lineares, não simétricos e
multidisciplinares inerentes ao problema, envolvendo, em geral, a solução numérica de
equações diferenciais parciais (EDP) através da aplicação de diversos métodos numéricos.
Entre eles, método de diferenças finitas (Peacemann, 1977, Aziz e Settari, 1979), método dos
elementos finitos (Chen, 2006), método dos volumes finitos e o método dos elementos finitos
mistos (Chen, 2006). Neste trabalho, entre outros, será considerado o método de diferenças
finitas em um esquema composto para a resolução numérica para os problemas cujos modelos
são descritos por equações diferenciais parciais hiperbólicas não lineares.
19
O esquema composto abordado aqui é um esquema conservativo, baseado em dois
métodos de diferenças finitas: O método de Lax-Friedrichs (LF), que é um método de
primeira ordem no tempo e no espaço e difusivo em regiões de descontinuidades, e o método
Lax-Wendroff (LW) que é um método de segunda ordem oscilatório em regiões descontínuas.
A composição entre eles fornece o método Composto chamado LWLF-4 sugerido por Liska e
Wendroff (1998), que foi testado para os problemas de Buckley-Leverett e Burgers,
respeitando a condição CFL, com resultados satisfatórios.
Com o objetivo de avaliar o desempenho do método composto, no presente trabalho
foram feitas comparações numéricas com as soluções analíticas dos problemas de Burgers e
Buckley-Leverett. Também foi feita uma comparação com o método do tipo TVD’s com
limitadores de fluxo. Os métodos TVDs, assim como o esquema composto, foram
desenvolvidos com o objetivo de evitar as oscilações espúrias tipicamente encontradas nas
soluções numéricas fornecidas pelos esquemas clássicos, apesar de ambos os métodos serem
eficientes para a captura do choque, dependendo da escolha do limitador de fluxo, o método
composto se mostrou superior e mais capaz na captura das descontinuidades.
Com o objetivo de incorporar conhecimentos matemáticos mais avançados, para em
seguida utilizá-lo no entendimento dos problemas estudados, o presente trabalho trata de
elaborar um estudo da teoria das Leis de Conservação, sendo estudadas as soluções fracas
que, fisicamente, podem ser interpretadas como ondas de choques ou rarefação, então, para
que possam ser distinguidas as fisicamente admissíveis, foi considerado o princípio de
entropia, em suas diversas formas. Logo, a teoria das equações diferenciais parciais (EDP)
hiperbólicas referente às Leis de Conservação não lineares e aspectos numéricos dessas leis
são tratadas em detalhes neste trabalho.
20
2. MODELAGEM MATEMÁTICA PARA O
ESCOAMENTO EM MEIOS POROSOS
Neste capítulo são vistas as equações matemáticas que modelam o escoamento de
fluidos em meios porosos dando ênfase ao processo de recuperação secundária de petróleo.
Tal processo é feito injetando água em alguns poços (denominados de poços de injeção) com
o objetivo de obter o óleo nos poços de produção. Assim, se o meio poroso é suficientemente
permeável, a mistura água-óleo ou água-óleo-gás deve escoar na direção desses poços de
produção, geralmente localizados estrategicamente ao longo do reservatório.
2.1 Modelagem matemática para o escoamento bifásico
O modelo matemático clássico para o processo de injeção de água em reservatório de
petróleo considera que o meio poroso encontra-se totalmente preenchido por uma mistura
formada por duas fases imiscíveis (água e óleo), as quais são supostamente incompressíveis e
não possuem espécies químicas em comum. Assim, tratam-se as fases água e óleo como se
fossem componentes químicos individuais, mesmo sabendo-se que nas condições típicas de
um reservatório ambas as fases apresentam mais de uma espécie química, principalmente o
óleo que é constituído por vários hidrocarbonetos.
2.1.1 Escoamento bifásico
Dois tipos de escoamento podem ocorrer quando dois fluidos escoam juntos em um
meio poroso: o deslocamento miscível e o imiscível. No primeiro caso, os fluidos são solúveis
um no outro não havendo uma interface bem definida entre os mesmos. No segundo tipo, que
será tratado nesse trabalho, há uma interface bem definida entre os fluidos. Nesse tipo de
escoamento imiscível, utilizam-se outras propriedades de interesse, as saturações das fases.
Considerando o escoamento bifásico água e óleo, define-se a saturação da água (denotada por
w
S ), como sendo a razão do volume da água contida em um dado volume do material poroso
pelo volume de poros, ou seja:
21
.
porosovolume
águadevolume
S
w
(2.1)
Análogo para a saturação do óleo, denotada por
o
S , temos
.
porosovolume
óleodovolume
S
o
(2.2)
Assim nota-se que as saturações das fases são grandezas (percentuais) adimensionais, tendo
seus valores pertencentes ao intervalo
1,0 .
Desde que o óleo e a água preencham totalmente o espaço vazio do meio poroso, então
temos das definições (2.1) e (2.2) a seguinte relação:
,1
S
onde
.,ow
(2.3)
O restante desse capítulo será dedicado à obtenção das equações diferenciais utilizadas
na descrição desse tipo de processo. As equações do escoamento bifásico são obtidas a partir
das equações fundamentais da Continuidade e de uma equação constitutiva conhecida como
Lei de Darcy.
2.1.2 Lei de Darcy
A Lei de Darcy, descoberta em 1856 pelo o engenheiro Henry Darcy, no contexto de
uma única fase, estabelece a relação básica entre a velocidade macroscópica do fluido e o
gradiente da pressão. Na verdade Darcy estabeleceu através de trabalhos empíricos
(experimentais) que a velocidade de um escoamento em um meio poroso é diretamente
proporcional à diferença da altura manométrica (Darcy, 1856) podendo então relacionar estas
duas grandezas através da introdução de uma constante:
,
l
h
Ku
(2.4)
22
onde
u
é a velocidade do escoamento,
h
é a diferença das alturas manométricas,
l
é a
distancia percorrida pelo fluido e
K
a constante de proporcionalidade denominada
condutividade hidráulica, sendo uma propriedade do material poroso e do fluido que escoa em
seu interior. A Figura 2.1 adiante ilustra o experimento de Darcy.
Figura 2.1 – Experimento de Darcy
Escrevendo a equação (2.4) na forma diferencial obtemos:
.
dl
dh
Ku
(2.5)
A pressão em um dado ponto por onde o fluido escoa pode ser relacionada com a
altura manométrica deste ponto através da expressão:
,Zhgp
(2.6)
onde
é a densidade ou massa específica do fluido,
g
é a magnitude da aceleração da
gravidade e
Z
a elevação em relação a um plano de referência.
23
Isolando
h
na equação (2.6) tem-se
.
1
gZ
p
g
h
(2.7)
Substituindo (2.7) em (2.5) obtemos:
.
gZ
p
dl
d
g
K
u
(2.8)
A constante
gK /
da equação (2.8) é válida apenas para o fluido água, entretanto, vários
experimentos com diversificados fluidos levaram a uma generalização da Lei de Darcy que
passa a ser expressa por:
,
gZ
p
dl
dk
u
(2.9)
onde k é a permeabilidade e
é a viscosidade do fluido.
Considerando o fluido incompressível, a massa específica
é constante, então a
equação (2.9) pode ser escrita como
.gZp
dl
dk
u
(2.10)
A permeabilidade do meio poroso
k
é um parâmetro que descreve a capacidade de
fluidos escoarem através da rocha que constitui o reservatório. Entretanto no caso de um
sistema bifásico é necessário introduzir o conceito de permeabilidade relativa
r
k , devido ao
fato de uma fase interferir no escoamento da outra. As permeabilidades relativas são definidas
de maneira que a seguinte relação seja satisfeita:
,1
Sk
r
(2.11)
24
,
k
k
k
r
(2.12)
onde
k
é a permeabilidade absoluta e
k é a permeabilidade efetiva da fase
.
A Lei de Darcy aplicada a escoamento bifásicos incompressíveis pode ser escrita para
cada fase da seguinte forma:
,gZp
dl
dk
ku
r
(2.13)
onde
denota a massa específica do fluido da fase
.
Introduzindo a mobilidade
/k na equação (2.13) e generalizando para problemas
multidimensionais obtemos
.Zgpku
(2.14)
2.1.3 Equação de balanço de massa
Apresentaremos inicialmente a derivação da equação que governa a conservação da
massa para um sistema monofásico. Para isso, considere um elemento de volume
V
com
fronteira V
como mostra a Figura 2.2.
Figura 2.2- Volume elementar
V
com fronteira
V
25
A conservação da massa implica que a taxa de massa acumulada em
V
é igual à taxa do fluxo
de massa através de suas fronteiras, somada à quantidade de massa injetada em V .
Então seja
n
o vetor normal ao volume elementar,
x
a porosidade do meio e
q
a taxa de massa por unidade de volume injetada em
V
. Logo, a conservação da massa pode
ser descrita pela equação
..
VVV
qdxdAnudx
dt
d
(2.15)
Aplicando o Teorema da Divergência,
VV
dxudAnu
..
e trocando
dtd /
com a integração espacial, pois supostamente
V
não varia com o tempo,
VV
dx
t
dxtx
dt
d
,
a equação (2.15) fica da seguinte forma:
..
VV
dxqudx
t
(2.16)
A equação acima é válida para qualquer volume de controle
V
. Então, fazendo
0
V
obtemos uma equação diferencial também conhecida como a Equação da Continuidade.
.. qu
t
(2.17)
Para obtermos a respectiva equação para um fluxo bifásico, primeiramente teremos
que assumir que as duas fases estão escoando simultaneamente, não havendo transferência de
26
massa entre as fases. Como os fluidos estão escoando juntos, eles estão dividindo o mesmo
volume poroso. Assim sendo, as equações de conservação de massa para cada uma das fases
deve considerar apenas a fração do volume poroso que contém massa dessa fase em questão.
Essa fração é a saturação da fase. Consequentemente, podemos escrever a equação de
conservação da massa para cada fase
de um fluxo bifásico como sendo:
..
qu
t
S
(2.18)
2.1.4 Desenvolvimento matemático para uma lei de conservação
Desenvolvendo a derivada temporal da equação da continuidade (2.18) para cada uma
das fases obtemos:
..
qu
t
S
t
S
t
S
(2.19)
Dividindo a equação (2.19) por
e reescrevendo-a para cada uma das fases temos para o
fluido água a seguinte equação:
w
w
w
www
w
w
w
w
qu
t
S
t
S
t
S
.
(2.20)
e para o óleo
.
.
o
o
o
ooo
o
o
o
o
qu
t
S
t
S
t
S
(2.21)
Somando as equações (2.20) e (2.21), temos:
27
.
..
w
w
o
owo
wo
w
w
wo
o
oo
o
o
w
ww
o
oo
qq
t
SS
t
SS
t
S
t
S
t
Suu
(2.22)
Substituindo a relação (2.3),
1
ow
SS
, e desenvolvendo os gradientes do lado esquerdo da
equação (2.22) tem-se:
.
....
w
w
o
ow
w
wo
o
o
w
ww
w
ww
o
oo
o
oo
qq
tt
S
t
S
uuuu
(2.23)
Como tanto a variação temporal da massa específica do fluido quanto da porosidade
dependem apenas da pressão no reservatório, a equação (2.23) toma a forma
.
.
w
w
o
ow
w
wo
o
o
o
o
ow
w
w
ow
qq
t
p
pt
p
p
S
t
p
p
S
p
p
u
p
p
u
uu
(2.24)
Com o objetivo de simplificar a equação (2.24), introduziremos algumas definições, entre elas
a velocidade total
T
u , fluxo total
Q e a compressibilidade da fase
c , como mostram
respectivamente às equações (2.25) a (2.27):
,
owT
uuu
(2.25)
,Q
qq
o
o
w
w
(2.26)
28
.
1
c
p
(2.27)
Assim, a equação (2.24) torna-se:
.
1
. Q
t
p
p
cScSpcucuu
oowwoowwT
(2.28)
Uma vez que o meio poroso é considerado nesse trabalho como incompressível, ou seja,
0
p
(2.29)
e como as fases também são consideradas incompressíveis temos que,
0
c
(2.30)
então a equação (2.28) pode ser reescrita na seguinte forma:
.. Qu
T
(2.31)
Como a mobilidade total (
T
) é a soma das mobilidades das fases e devido à Lei de Darcy
(2.14), sendo desprezados os efeitos gravitacionais, a equação (2.31) pode ser modificada,
chegando-se à equação que descreve o comportamento da pressão em todo o reservatório,
como segue:
.. Qpk
T
(2.32)
Para obter a equação da lei de conservação que descreve o comportamento da
saturação da água no reservatório, devemos primeiramente somar e subtrair a velocidade do
29
fluxo de cada fase, lembrando que os efeitos gravitacionais estão sendo desprezados. Logo,
desta soma teremos:
pkuu
owow
(2.33)
e da subtração obtemos:
.pkuu
owow
(2.34)
Explicitando o operador gradiente, a equação (2.33) também pode ser escrita como
.
T
T
k
u
p
(2.35)
Análogo para a equação (2.34), temos:
.
T
ow
k
uu
p
(2.36)
Igualando as equações (2.35) e (2.36) chega-se a relação:
.
T
T
o
T
T
w
ow
uuuu
(2.37)
Assim, obtemos uma relação para velocidade do fluido de cada fase, como mostram as
relações seguintes:
T
T
w
w
uu
e
.
T
T
o
o
uu
(2.38)
Em seguida, definimos o fluxo fracionário como:
30
,
T
f
(2.39)
onde a descrição matemática do fluxo fracionária com seus respectivos modelos será vista
com mais detalhes nos próximos capítulos. Continuando a dedução, notamos que as
velocidades das fases podem ser reescritas da forma
Tww
ufu
e
.
Too
ufu
(2.40)
Utilizando a equação da continuidade (2.18) e introduzindo nela a relação (2.40) para o fluido
água, obtemos a equação para leis de conservação que descreve o comportamento do campo
de saturação desse fluido.
..
w
w
Tww
w
q
uSf
t
S
(2.41)
Neste trabalho, estaremos interessados no problema de transporte convectivo,
unidimensional e sem difusão, com a porosidade constante, onde não haverá termo de fonte
ou sumidouro e nem a ação da gravidade. Portanto, nessas condições, a equação (2.41)
resume-se a
.0
~
x
Sf
u
t
S
ww
T
w
(2.42)
onde
T
u
~
é a velocidade total dos fluidos dividida pela porosidade do reservatório.
2.2 Equação de Buckley-Leverett
O problema de Buckley-Leverett é clássico na engenharia de exploração de petróleo,
sendo utilizado na recuperação secundária. A equação que modela esse problema tem o
objetivo de fornecer o perfil da saturação da água ao longo do reservatório com o passar do
tempo, constituindo uma equação diferencial parcial não linear, parabólica quando os efeitos
31
difusivos são considerados e hiperbólica quando esses efeitos são desconsiderados. A seguir,
será feita uma análise para o problema de Buckley-Leverett.
2.2.1 Dedução da equação
Para deduzir a equação de Buckley-Leverett, partiremos da equação (2.31), porém,
para um problema unidimensional, ou seja:
.Q
x
u
T
(2.43)
A equação de Darcy (2.14) para um problema unidimensional, onde são desprezados os
efeitos gravitacionais, fica (para cada fase) com as formas:
,
x
p
ku
w
ww
(2.44)
.
x
p
ku
o
oo
(2.45)
Introduziremos o conceito da pressão capilar como sendo a diferença entre as pressões de
cada fase, ou seja,
.
wowc
ppSp
(2.46)
Empiricamente, a pressão capilar é descrita em função da saturação da água. Assim a
velocidade do fluido água pode ser modelado como:
.
x
p
k
x
p
ku
c
w
o
ww
(2.47)
Isolando a diferença parcial da pressão do óleo na equação (2.47) e em seguida substituindo
na equação (2.45), obtemos uma expressão para a velocidade do fluido óleo.
32
.
x
p
k
u
x
p
k
u
ku
c
o
w
woc
w
w
oo
(2.48)
Como
owT
uuu
, segue:
.1
x
p
kuu
c
o
w
o
wT
(2.49)
Introduzindo a equação do fluxo fracionário da água (2.39) na equação (2.49), temos:
x
p
k
f
u
u
c
o
w
w
T
(2.50)
ou
.
x
p
fkufu
c
owTww
(2.51)
Substituindo a equação (2.51) em (2.18) para a fase água, com 0
w
q , vem:
.
t
S
x
p
kfu
x
wc
wo
ow
wT
(2.52)
Aplicando a regra da cadeia para a pressão capilar
x
S
Sp
x
S
S
p
x
p
w
wc
w
w
cc
(2.53)
e substituindo na equação (2.52) obtemos:
33
.0
t
S
x
S
pk
xx
fu
ww
c
wo
owwT
(2.54)
É notório que, como 0/
wcc
dSdpp , o termo
c
wo
ow
FIS
pkK
(2.55)
tem na equação (2.54) o caráter de uma dispersão física, esse termo é conhecido como tensor
de dispersão física. Logo, na equação de Buckley-Leverett a capilaridade é responsável pela
introdução de um termo difusivo.
A equação (2.54) pode ser adimensionalizada definindo-se:
,
L
x
x
D
(2.56)
.
0
t
T
D
dt
L
u
t
(2.57)
Substituindo as equações (2.56) e (2.57) em (2.55) obtemos:
,0
1
D
w
D
w
FIS
TD
w
t
S
x
S
K
xLux
f
(2.58)
onde o quociente
LuK
TFIS
/
pode ser escrito como:
.
/1
1
w
c
ow
w
TT
FIS
dS
dp
k
LuLu
K
(2.59)
Com a finalidade de obter a equação de Buckley-Leverett, será necessário introduzir as
definições para as permeabilidades relativas de cada fase. A permeabilidade relativa é uma
34
função das saturações dos fluidos, e é expressa por uma curva como a da Figura 2.3 que pode
ser aproximada por um modelo exponencial, como mostra a seguir:
,
1
0
a
orwi
wiw
rwrw
SS
SS
kk
(2.60)
,
1
1
0
b
orwi
wor
roro
SS
SS
kk
(2.61)
onde
or
S
é a saturação do óleo residual,
wi
S
é a saturação da água irredutível,
0
rw
k
é a
permeabilidade relativa da água na saturação residual do óleo,
0
ro
k
é a permeabilidade relativa
do óleo na saturação residual da água e os expoentes
a
e
b
são constantes reais não
negativos. Os modelos de potência representam com bastante precisão os valores de
permeabilidade relativa medidos experimentalmente.
Figura 2.3 – Gráfico típico de permeabilidade relativa água/óleo
Podemos também utilizar o conceito da função J de Leverett (Poston et al., 1970;
Bear-Bachmat, 1991), adimensional:
0 1
0
1
Saturação da Água
Permeabilidade Relativa
Permeabilidade Relativa x Saturação da Água
Permeabilidade Relativa do Óleo
Permeabilidade Relativa da Água
35
.
cos
kSp
SJ
wc
w
(2.62)
Onde
é a tensão superficial dos fluidos e
é o ângulo de contato da interface rocha e
água. Então, isolando a pressão capilar da expressão acima e devirando ambos os lados com
relação à saturação da água, obtemos:
.cos
ww
c
dS
dJ
kdS
dp
(2.63)
Substituindo a expressão (2.63) em (2.59), vem:
.cos
/1
1
1
1
0
wow
a
orwi
wiw
w
rw
TT
FIS
dS
dJ
kSS
SSkk
LuLu
K
(2.64)
Logo, definimos:
w
a
orwi
wiw
w
w
a
orwi
wiw
ow
w
dS
dJ
SS
SS
f
dS
dJ
SS
SS
S
1
1
1/1
1
(2.65)
como sendo o número de Rapoport e Leas, o qual fornece uma razão entre as forças viscosas
(convectivas) e capilares (difusivas). Quanto maior o seu valor, mais convectivo é o fluxo,
mais detalhes em (Lake, 1989). Considerando
cos
0
rw
Tw
RL
k
Lu
k
N
(2.66)
obtemos a seguinte relação:
.
RL
w
T
FIS
N
S
Lu
K
(2.67)
36
Substituindo-se em (2.58), resulta finalmente a equação de Buckley-Leverett em sua forma
adimensionalizada,
.
1
D
w
w
DRLD
w
D
w
x
S
S
xNx
f
t
S
(2.68)
A equação de Buckley-Leverett é composta por uma parte convectiva (parte esquerda da
equação) e uma parte difusiva (parte direita). Ela é estritamente não linear e parabólica, porém
em nosso estudo estaremos interessados na parte convectiva desprezando os efeitos difusivos.
Logo, a equação (2.68) se torna estritamente hiperbólica não linear, como é mostrado na
equação abaixo,
.0
D
w
D
w
x
f
t
S
(2.69)
Aplicando a regra da cadeia para a função de fluxo fracionário, obtemos uma nova versão
para a equação de Buckley-Leverett, sem os índices de adimensionalização.
.0
x
S
f
t
S
w
w
w
(2.70)
Existem vários modelos para a função de fluxo fracionário para o problema de
Buckley-Leverett. A função de fluxo é não linear e não convexa, ela representa a transferência
da quantidade da saturação que ocorre em determinadas regiões do escoamento em certos
intervalos de tempo. Para o nosso estudo utilizaremos o modelo de Stone (Sethian, 1983),
esse modelo além da saturação do fluido molhante (no caso a água) também depende da razão
viscosidade entre os fluidos sendo dada por:
.
1
2
2
2
w
o
w
w
w
w
SS
S
Sf
(2.71)
37
Pode ser vista em (Flensberg, 2005) a aplicação de um outro modelo para a função de fluxo
fracionário para o problema de Buckley-Leverett, o modelo de Brooks e Corey, mais detalhes
em (Brooks, 1964), o qual também se mostrou satisfatório na resolução do problema
proposto, sendo dependente da saturação residual, como mostra a equação a seguir:
,
11
2
2
4
4
eee
e
w
SSS
S
Sf
(2.72)
onde
6.0/
rwwe
SSS
é chamada de saturação efetiva.
Adiante, nas Figuras 2.4 e 2.5, podemos ver os gráficos para a função de fluxo
fracionário em função da saturação para os respectivos modelos.
Figura 2.4 – Gráfico da função de fluxo Figura 2.5 – Gráfico da função de fluxo
fracionária para o modelo de Stone em fracionária para o modelo de Brooks-Corey
distintos valores para a razão de viscosidade (Brooks, 1964).
ow
c
/
.
2.2.2 Solução semi-analítica
Para uma frente constante de saturação
w
S
, da equação (2.70) temos:
.0
dx
x
S
dt
t
S
dS
t
w
x
w
w
(2.73)
Assim,
0 0.2 0.4 0.6 0.8 1
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
Saturação (S)
Função de Fluxo para o Modelo de Stone's
Função de Fluxo f(S)
c=1
c=0.5
c=0.2
0 0.2 0.4 0.6 0.8 1
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
Saturação (S)
Função de Fluxo para o Modelo de Brooks e Corey [7]
Função de Fluxo f(S)
38
.
/
/
*
ww
t
w
x
w
sw
Sf
xS
tS
dt
dx
v
(2.74)
Portanto,
..| tSfx
wwS
w
(2.75)
A relação (2.75) possui uma solução fisicamente incoerente na frente do avanço da
água, a qual será vista em mais detalhes em seguida. Para que a solução seja selecionada de
forma correta é necessário utilizar a teoria do choque, que será tratada no Capítulo 3. A
velocidade correta no ponto de saturação de choque
*
w
S
é determinada através da aplicação do
balanço de massa, ou, equivalentemente, pela a aplicação do princípio de Rankine-Hugoniot
(Capítulo 3). A saturação da água no ponto de choque é obtida igualando-se essa velocidade
com aquela calculada usado (2.74). A equação (2.75) é utilizada para obter as posições das
saturações da água maiores que
*
w
S
. Não é difícil verificar que a aplicação do princípio da
entropia também conduz à saturação correta no choque, ou seja:
• Se o choque fosse fixado para
*
ww
SS
, teríamos soluções com três valores distintos para a
saturação na mesma posição, na frente de avanço, ocorrendo uma incoerência física.
• Se o choque fosse fixado para
*
ww
SS , a condição de entropia na descontinuidade seria
violada.
Em seguida, será visto com mais detalhes o avanço de água e a solução analítica de
Welge (1952) para a captura do choque mencionado.
2.2.3 Localização da frente de água
O processo de injeção de água aumenta a pressão no interior do reservatório e
impulsiona o óleo para os poços de produção. Entretanto, devido a tensão superficial e a
diferença entre os valores da viscosidade existente entre os fluidos imiscíveis água e óleo,
após uma certa quantidade de água injetada é observado que o óleo residente resiste ao
deslocamento e tende para um estado de imobilidade, o qual permite que a fase água avance
em grande quantidade na direção dos poços de produção, prejudicando a recuperação do
petróleo. Logo, com o objetivo de prever esse avanço da água, baseando na equação de
39
Buckley-Leverett não adimensionalizada, encontraremos uma equação que modela o avanço
de água no reservatório.
Introduziremos a relação (2.74) na equação de Buckley-leverett não
adimensionalizada, obtendo uma relação análoga a relação (2.74), porém com a presença da
velocidade total dos fluidos, como mostra a equação (2.76).
.
~
wwT
Sftu
dt
dx
(2.76)
A equação (2.76) fornece, a cada instante de tempo, a velocidade da frente de água com
valores de saturação constante que avança no meio poroso na direção dos poços de produção.
Logo, integrando a equação (2.76) em relação ao tempo, obtemos
dttuSfdt
dt
dx
t
Tww
t
00
~
(2.77)
ou de forma equivalente,
,
11
0
00
AdttuSfdt
tu
Sfxtx
t
Tww
t
T
ww
(2.78)
onde
é a área da seção do meio poroso transversal ao escoamento, supostamente constante
ao longo da direção
x
. Seja
dtutdV
T
(2.79)
a diferencial do volume total da mistura que escoa através da área
. Então a equação (2.78)
pode ser reescrita como,
tV
V
ww
tdVSfxtx
0
1
0
(2.80)
40
ou seja,
,0
ww
Sf
tV
xtx
(2.81)
onde
0VtVtV representa o volume total de fluido acumulado no intervalo de
tempo que vai de zero a
t
. Suponha que
00 x e
00 V , de modo que a equação (2.81)
pode ser reescrita como:
.
ww
Sf
tV
tx
(2.82)
A equação (2.82) informa a localização da frente de avanço de água em função do valor da
saturação
w
S
, para 0
t .
2.2.4 Solução analítica de Welge para a equação de Buckley-Leverett
Considerando o modelo de Stone para a função de fluxo fracionário da água, descrito
na equação (2.71) e tomando a sua derivada temos
.
1
12
2
2
2
w
o
w
w
ww
o
w
ww
SS
SS
Sf
(2.83)
A partir das equações (2.82) e (2.83), chegamos à seguinte equação algébrica:
.0121
2
2
2
ww
o
w
w
o
w
w
SS
tV
xSS
(2.84)
A equação (2.84) define uma curva no plano
w
Sx , para 0
t . A Figura 2.6 mostra o
esboço dessa relação para o caso particular onde
2/1/
ow
e
.1/
tV
41
Figura 2.6 – Distribuição da Saturação da água a partit da equação de Buckley-Leverett.
Como já tínhamos mencionado antes, o gráfico do avanço da água mostra que a
solução obtida determina três valores distintos para a saturação da água relacionados com um
único ponto
x
. Esta incoerência física tem origem no fato de não ter sido incluído no modelo
apresentado no presente trabalho o efeito das forças difusivas decorrentes da pressão capilar,
ver equação (2.68). Para eliminar esses valores múltiplos da saturação, consideramos uma
solução descontínua, como mostra a Figura 2.7.
Figura 2.7 – Solução descontínua para a equação de Buckley-Leverett
Na Figura (2.7),
c
x é chamada abscissa de corte e
*
w
S é a saturação superior de corte ou valor
de saturação no ponto de choque.
Na década de cinquenta, Welge (1952) mostrou como determinar o valor da saturação
no ponto de choque e consequentemente o valor da abscissa de corte a partir da geometria da
42
função de fluxo fracionária da água
w
f . Considerando a descontinuidade na abscissa de corte,
é notável que o volume de fluido no intervalo
c
x,0 em um determinado instante de tempo
t
é dado por:
.
0
dxStV
c
x
w
(2.85)
Integrando por partes, a equação (2.85) pode ser escrita como segue,
.
*
1
*
w
S
wc
dStxSxtV
w
(2.86)
Substituindo a equação (2.82) no último termo do lado direito da equação (2.86), temos,
www
S
wc
dSSf
tV
SxtV
w
*
0
*
(2.87)
ou seja,
.1
**
1
*
*
wwwwcw
S
wwwc
fSftVSxdSSftVSxtV
w
(2.88)
Fazendo 1
w
S no modelo de Stone para a função de fluxo fracionário, obtemos
11
w
f .
Assim, podemos reescrever a equação (2.88) da seguinte forma equivalente
1
*
wwwc
SftVSxtV
(2.89)
ou seja,
43
.
*
*
tV
S
Sf
x
w
ww
c
(2.90)
Por outro lado, considerando a equação do avanço da saturação da água (2.82) para
c
xx
temos
.
*
tVSfx
wwc
(2.91)
Igualando as equações (2.90) e (2.91) obtemos uma expressão para o valor de saturação no
ponto de choque, como mostra a relação adiante.
.
*
*
*
ww
ww
w
Sf
Sf
S
(2.92)
Logo podemos ver que o valor da saturação no ponto de choque depende rigorosamente do
modelo da função de fluxo e a sua derivada. Welge (1952) também mostra essa relação
geometricamente; de fato, observando a Figura 2.8 notamos que o valor da saturação no ponto
de choque no qual a reta que passa pela origem 0
w
S e 0
w
f é tangente ao gráfico da
função de fluxo, exatamente no ponto de coordenadas
*
ww
SS e
*
www
Sff .
Figura 2.8 – Determinação geométrica da saturação no ponto de choque
44
Tomando a equação da derivada da função de fluxo (2.83) e substituindo na equação
(2.92), verificamos que a saturação no ponto de choque depende apenas da razão das
viscosidades dos fluidos, como mostrado na seguinte equação,
.
/1
/
*
ow
ow
w
S
(2.93)
Essa relação é de extrema importância para o nosso estudo, pois será comparado com a
solução numérica no problema de Buckley-Leverett.
Será visto no próximo capítulo, em mais detalhes que, a equação de Buckley-Leverett
é uma lei de conservação não linear e a solução descontínua mostrada na Figura 2.7 é uma
solução fraca dessa lei, a qual satisfaz a condição de entropia.
2.3 Equação de Burgers
A equação de Burgers também é muito utilizada na matemática e em engenharia. Essa
equação, assim como a de Buckley-Leverett, é uma lei de conservação hiperbólica não linear.
A equação de Burgers difere da equação de Buckley-Leverett apenas na função de fluxo, a
qual é mostrada a seguir,
.
2
1
2
ww
SSf
(2.94)
A equação (2.94) representa o modelo da função de fluxo fracionário para o problema de
Burgers, ela é não linear e convexa, como é visto na Figura 2.9.
45
Figura 2.9 – Gráfico da função de fluxo para o problema de Burgers
Introduzindo a equação do fluxo do modelo de Burgers na equação de Lei de conservação do
tipo mostrado na equação (2.70) obtemos:
.0
2
1
2
x
S
t
S
w
w
(2.95)
Aplicando a regra da cadeia para a equação (2.95) temos:
.0
x
S
S
t
S
w
w
w
(2.96)
A equação acima é a forma mais geral para o problema de Burgers.
-1 -0.5 0 0.5 1
0
0.05
0.1
0.15
0.2
0.25
0.3
0.35
0.4
0.45
0.5
Saturação (S)
Função de Fluxo para o Modelo de Burgers
Função de Fluxo f(S)
46
3. LEIS DE CONSERVAÇÃO
3.1 Formulação
Nesse capítulo serão estudadas as equações diferencias do tipo
,0.
SfS
t
(3.1)
onde
x , e
:f é uma função não linear. Tais equações são EDP’s hiperbólicas e
ocupam um lugar de destaque na física, pois descrevem fenômenos onde há conservação de
alguma grandeza tal como conservação de massa, energia, momento, sendo por isto objeto de
estudo importante no mundo científico. Pode-se dizer que a equação (3.1) é uma lei física de
caráter local que surge na ciência aplicada oriunda de leis de conservação globais ou integrais
originalmente escritas como
..
dsnfSdx
dt
d
(3.2)
A equação (3.2) afirma que a taxa de variação da propriedade S ao longo de toda uma
região
N
(conjunto limitado e conexo) é igual ao fluxo total que atravessa a fronteira de
, denotada por
, onde a aplicação f representa localmente o fluxo da propriedade e
n
é
o vetor normal exterior a
. Se a região
não varia com o tempo, pode-se reescrever a
equação (3.2) como segue
..
dsnfdxS
t
(3.3)
Note que uma lei de conservação na sua forma integral (equação (3.2)) não exige que
o fluxo f seja diferenciável. Mas supondo-se que isso ocorra, o teorema da divergência diz
que:
47
...
dxfdsnf
(3.4)
Logo, a equação (3.3) pode ser escrita na seguinte forma
.0.
dxfS
t
(3.5)
Se o integrando
fS
t
. é supostamente contínuo em
, como a região
é tomada de
forma arbitrária, então da igualdade na equação (3.5) obtém-se imediatamente a lei de
conservação local descrita pela equação (3.1). É enfatizado mais uma vez que a passagem da
formulação global para a formulação local exige que a função incógnita S e o fluxo f
possuam um certo grau de regularidade.
Neste capítulo são tratados problemas de valor inicial, também chamados de
problemas de Cauchy, associados às leis de conservação:
xxStxS
txSfS
xt
0
0,
0,0
(3.6)
A equação hiperbólica no problema descrito na equação (3.6) se apresenta numa forma
denominada conservativa. Entretanto, se
S
for suave, aplicando a regra da cadeia esse
problema pode ser escrito em sua forma não conservativa, como mostra a equação adiante,
xxStxS
txSfSS
xt
0
0,
0,0
(3.7)
O problema de Cauchy resume uma série de problemas releventes da engenharia e da
física aplicada, por exemplo, o problema relacionado com a equação de Burgers pode ser
visto como um protótipo unidimensional da conhecida equação de Euler, usada para descrever
a conservação do momento linear de um fluido incompressível e sem viscosidade. Neste caso,
a incógnita
txSS , é interpretada como sendo a velocidade unidimensional de uma
48
partícula de fluido que passa pela posição
x
no instante
t
. Na mecânica dos meios contínuos,
a equação de Burgers (2.96) é muitas vezes escrita na seguinte forma,
,0
Dt
DS
(3.8)
onde xStDtD
/// é denotada derivada material ou Lagrangeana. Assim, do ponto
de vista da mecânica Lagrangeana, a equação de Burgers (2.96) descreve o fato de que num
escoamento incompressível e sem viscosidade a velocidade de uma dada partícula fixada se
mantém constante ao longo do caminho percorrido por ela no campo de escoamento.
O problema de Buckley-Leverett com condição de Cauchy é dado por:
0,
0,
0,
0,0
xseS
xseS
xS
txSfSS
R
L
xt
(3.9)
onde 10
LR
SS .
3.2 Formação de ondas de choque e ondas de rarefação
Seja
txS , uma função que depende apenas de uma variável espacial e do tempo.
Além disso, considere que a função de fluxo
txSff , é suficientemente diferenciável.
Definição 3.1 Diz-se que
txSS , é uma solução clássica do problema de Cauchy para uma
lei de conservação, se
txS , é continuamente diferenciável no semiplano superior aberto
0e;,
22
txtx e satisfaz a EDP e a condição inicial descritas na
equação (3.6).
Proposição 3.1 Se
txSS , é uma solução clássica do problema de valor inicial indicado na
equação (3.6), entãoS satisfaz à fórmula implícita
49
.0,;,,
0
txttxSfxStxS
(3.10)
Demonstração: Considerando que a condição inicial
0
S é uma função diferenciável, sabe-se,
do teorema da função implícita, que a expressão em (3.10) pode ser resolvida explicitamente
para S como função de
x
e
t
(para
t
pequeno). Uma vez que isso é possível, aplicando a
regra da cadeia, obtêm-se suas derivadas parciais, como mostram as equações abaixo:
,
1
0
0
tSfS
SSf
S
t
(3.11)
.
1
0
0
tSfS
S
S
x
(3.12)
As equações (3.11) e (3.12) satisfazem a EDP
0
SfSS
xt
. Evidentemente, a solução
mostrada na equação (3.10) cumpre também a condição inicial em (3.6), o que conclui a
demonstração.
Definição 3.2 Uma curva característica da EDP
0
SfSS
xt
é uma solução da EDO
.,ttxSftx
dt
d
(3.13)
Proposição 3.2 Uma solução clássica do problema de Cauchy para uma lei de conservação é
constante ao longo de qualquer curva característica definida pela equação (3.13).
Demonstração: Supondo que
ttxSS , é a solução da equação
0
SfSS
xt
, pela a
regra da cadeia, tem-se imediatamente que
0,,,
,,,
ttxSttxSfttxS
ttxStx
dt
d
ttxSttxS
dt
d
tx
tx
(3.14)
com isso é demonstrado o seguinte resultado.
50
Considera-se uma curva característica arbitrária parametrizada pelo tempo, dado por
txx , onde (de acordo com a condição inicial)
00
, SttxS
t
é a solução do
problema (3.6), para 0
t e
,0
x . Visto que S é constante ao longo de
qualquer curva característica, então, por (3.13), a curva
txx satisfaz a EDO
0
/ Sfdtdx
e a condição inicial
0x
. Consequentemente, as curvas características
do problema (3.6), para o caso estudado, são retas da forma:
tSftx
0
(3.15)
para
cte
.
Como denotado acima, o valor constante de S
ao longo da reta característica que
emana de um ponto arbitrário
0,
do plano
t
x
é dado por
0
S . Se essa reta forma um
ângulo
com o eixo das abscissas (veja Figura 3.1), então o seu coeficiente angular é
.
1
tan
0
Sf
(3.16)
Figura 3.1 – Inclinação da reta característica
Assim, se
51
121020
SfSf
(3.17)
então, tem-se
21
tantan
. Isto implica que duas retas características quaisquer nunca se
cruzam, como mostra a Figura 3.2.
Figura 3.2 – Retas características que não se cruzam
Considerando a equação de Burgers, onde
SSf
, a condição descrita na equação
(3.17) torna-se
.
121020
SS
(3.18)
Assim, utilizando considerações relativas ao problema físico modelado pela a equação de
Burgers, em vista da relação (3.18) acima, pode-se dizer que a partícula que estava localizada
na posição
1
, no instante de tempo inicial 0
t , possui velocidade
10
S menor ou igual à
velocidade
20
S da partícula que se encontrava na posição
2
, no instante de tempo inicial
0
t
. Uma vez que
1020
SS e pelo o fato de tais partículas viajarem com velocidades
constantes ao longo das retas características que emanam, respectivamente, de
1
e
2
,
fisicamente é de se esperar que elas nunca colidam, ou seja, que os caminhos no plano
t
x
nunca se cruzem ou se choquem, o que está inteiramente de acordo com a Figura 3.2. Observe
que, se
0
S
for uma função não decrescente, ou seja, se
0
0
S
, para todo S , então a
relação (3.18) deve ocorrer sempre.
52
Curvas características decorrentes da relação descrita na equação (3.18) são
denominadas de curvas de rarefação. Neste caso, as soluções associadas a este tipo de
problema são chamadas de onda de rarefação. Na Figura 3.3 é esboçado um leque de curvas
de rarefação.
Figura 3.3 – Um típico leque de rarefação
Por outro lado, se existirem duas posições
1
e
2
, com
21
, nas quais ocorrem
,
1020
SfSf
(3.19)
então existirão retas características emanando de
1
e
2
, respectivamente, e possuindo
coeficientes angulares
1
tan
e
2
tan
tais que
21
tantan
. Portanto, essas características
necessariamente se cruzam em um ponto, como mostrado na Figura 3.4.
53
Figura 3.4 – Duas retas características que se cruzam
Voltando ao problema particular da equação de Burgers, na situação esboçada na
Figura 3.4, pode-se dizer que existe um instante
0
t
onde a partícula que estava na posição
1
, no tempo inicial
0
t
, colide com a partícula que se encontrava na posição
2
neste
mesmo instante inicial. Entretanto, isto é intuitivamente esperado, pois em vista da equação
(3.19), a velocidade da partícula que estava inicialmente atrás é maior do que a velocidade da
partícula que se encontrava a frente, ou seja,
.para
211020
SS
(3.20)
Quando a relação (3.19) é satisfeita, diz-se que as curvas características planas
apresentam choques. Observa-se que, se
0
S é uma função decrescente, ou seja, se
,0
0
S
(3.21)
então as curvas características da equação de Burgers exibirão choques.
Para começar a entender as possíveis complexidades associadas com a resolução de
um problema hiperbólico não linear, é considerada a condição descrita na equação (3.20) e
analisada a solução do problema de Burgers sobre um ponto de choque de duas retas
54
características que emanam de valores distintos
1
e
2
, como mostra na Figura 3.4. Por um
lado, desde que a solução desse problema (velocidade) se mantenha constante ao longo das
retas características, então
S
sobre a característica que emana do ponto
1
é dada por
10
S
e a solução ao longo da reta que emana
2
é
20
S . Por outro lado, como essas retas se
cruzam em um ponto, então, seguindo essa argumentação, conclui-se que neste ponto
2010
SS , o que contradiz a desigualdade descrita na equação (3.20). Este é o caso
típico onde a solução do problema não linear é uma onda de choque.
Com a argumentação desenvolvida no último parágrafo, pode-se observar que a
solução clássica considerada na Proposição 3.2 nem sempre representa, pelo menos com a
regularidade exigida até aqui, a realidade física do problema em questão. De fato, tal
argumentação sugere que, na presença de choques, a solução procurada deve ser descontínua.
Para elucidar um pouco mais este fato, considere o seguinte resultado.
Proposição 3.3 Sejam
t
S e
x
S as derivadas descritas nas equações (3.11) e (3.12),
respectivamente, e suponha que
0
f
. Assim, tomando apenas
0
t
:
(i) Se 0
0
S para todo
x
, então as funções
t
S e
x
S se mantém limitadas,
independente de
0
t
.
(ii) Se 0
0
S para algum
x
, então
t
S e
x
S tendem para o infinito quando a expressão
tfxS
0
1
tende para zero.
Demonstração: Para o primeiro caso, tem-se que
,
1
0
0
tSfS
SSf
S
t
(3.22)
.
1
0
0
tSfS
S
S
x
(3.23)
Assim, se 0
0
S então 0
xt
SS e a limitação das derivadas é óbvia. Se 0
0
S para todo
x
, das equações (3.22) e (3.23), nota-se que
t
S e
x
S só assumem valores arbitrariamente
55
grandes na vizinhança do tempo negativo
SfSt
0
/1 . Logo, como apenas são
considerados tempos positivos, neste caso, as derivadas
t
S e
x
S sempre se mantém limitadas,
independentemente do valor de
0
t
. Já para o segundo item, tem-se:
,
1
0
0
tSfS
SSf
S
t
(3.24)
.
1
0
0
tSfS
S
S
x
(3.25)
Desde que
0
0
xS
, faz sentido considerar
01
0
tfxS
, pois tal limite ocorre no
tempo positivo, então
,
1
0
fxS
t
(3.26)
onde as derivadas parciais
t
S
e
x
S
se tornam ilimitadas, o que conclui a demonstração.
A proposição acima mostra alguns possíveis aspectos da solução do problema descrito
em (3.6). Supõe-se que 0
f . Com isso, inicialmente, nota-se que um sinal negativo de
0
S
(a derivada da condição inicial) determina a formação de choque. De fato, se existir
1
tal
que
0
10
S , então
2010
SS para algum
12
. Como f
é por hipótese
crescente, segue que
1020
SfSf
, o que satisfaz a condição de choque descrita na
equação (3.19). Caso contrário, se
0
10
S , para todo
x
, então (desde que 0
f ) choques
não ocorrerão. Assim, considerando 0
f , o segundo item da Proposição 3.3 mostra que na
presença de choque as funções
t
S e
x
S , as derivadas parciais de primeira ordem da solução
txSS ,
, apresentam singularidades. Nesta situação, a EDP
0
SfSS
xt
certamente
não pode estar modelando corretamente o problema físico considerado, para todo
0
t
, pois
ela exige regularidade sobre
S
as quais claramente não mais existem. Logo, a princípio,
parece aceitável concentrar-se em soluções clássicas (ou regulares) que existem apenas para
56
um certo intervalo de tempo, definido para valores de
t
menores do que o instante de tempo
indicado na equação (3.7), onde
t
S
e
x
S
assumem valores arbitrariamente grandes,
caracterizando a existência de choques. Entretanto, como já enfatizado, leis de conservação
surgem naturalmente na modelagem matemática de vários fenômenos físicos, logo a análise e
obtenção de soluções, mesmo que descontínuas, fazem parte do estudo de tais fenômenos,
para qualquer instante de tempo. Na próxima seção, será estudada a maneira de tratar este
problema.
Na presença das hipóteses consideradas na Proposição 3.3, percebe-se que é possível
determinar o instante do primeiro choque ou “breaking point”, denotado aqui por
*t
. Para
isso, considera-se o par
**,tx
indicado na Figura 3.5.
Figura 3.5 – Instante do primeiro choque
Certamente, o par
**,tx
está sobre uma reta característica, que emana de um ponto
,
min
, descrita por
minmin0
tSfx
(3.27)
e, além disso, satisfaz a equação (3.26), ou seja,
57
min00
1
SfxS
t
(3.28)
para
*tt
e
*xx
. Com a hipótese 0
f e considerando que exista
x
tal que
0
0
xS , a
Proposição 3.3 indica que o instante de um choque satisfaz à equação (3.7). Assim, para que
*t
(satisfazendo a equação 3.28) seja de fato o instante do primeiro choque, é necessário que
os parâmetros desconhecidos,
*t
e
min
, minimizem a função
00
* SfxS
, o
denominador da equação (3.28), para todo
,
e
0
t
. Consequentemente, de
modo a se obter
*t
e
min
, pode-se substituir a equação (3.27) em (3.28). Com isso, nota-se
que a função a ser minimizada é no caso geral, uma função de duas variáveis, dada por:
.,
000
SftSfSth
(3.29)
Este resultado é resumido adiante.
Proposição 3.4 Suponha que 0
f e existe
x
tal que
0
0
xS . Se o par
*,
min
t
é um
minimizador global da função descrita na equação (3.29), então o primeiro choque do
problema (3.6) ocorre no instante
.
*,
1
min
th
t
(3.30)
Para ilustrar o resultado desenvolvido acima, considera-se o problema de Cauchy
relacionado com a equação de Burgers. Neste caso, sabe-se que 1
f . Para esse problema,
nota-se que a função indicada na equação (3.29) é mais facilmente descrita em termos da
variável
x
, ou seja,
xSxh
0
. Portanto, para a equação de Burgers, o “breaking point” é
dado pela relação simples:
,
*
1
*
0
xS
t
(3.31)
58
onde
*x
minimiza a função
xSS
00
, a derivada da condição inicial.
Deve-se notar que a expressão na equação (3.30) não pode ser utilizada para
determinar o primeiro instante do choque de um problema de Cauchy com condições iniciais
na forma de uma função em escada, como considerada no problema de Buckley-Leverett. De
fato, neste caso, não existe
x
tal que
0
0
xS .
3.3 Soluções descontínuas e a condição de Rankine-Hugoniot
Esta seção será iniciada com um exemplo que procura ilustrar a existência de uma
ampla classe de soluções para problemas de Cauchy (equação (3.6)).
Exemplo 3.1: É considerado o problema de valor inicial referente à equação de Burgers,
xStxS
txSSS
xt
0
0,
0,0
(3.32)
onde a condição inicial é
10
101
01
0
x
xx
x
xS
(3.33)
A derivada da função (3.33) é dada por:
10
101
00
0
x
x
x
xS
(3.34)
O esboço das funções
xS
0
e
xS
0
são mostrados abaixo:
59
Figura 3.6 – Os gráficos das funções
xS
0
e
xS
0
.
Desde que
0
0
xS para todo
1,0x
, então este problema deve apresentar choques. Como
mostrado na Figura 3.7, de qualquer ponto
emana uma reta característica definida por
uma das seguintes relações:
1,
10,1
0,
t
t
x
(3.35)
Figura 3.7 – Características definidas pela equação (3.35)
60
Das equações (3.31) e (3.34) nota-se que o “breaking point” deste problema é
1*
t
, onde
neste caso
1*
0
xS
. Como uma abordagem inicial, procura-se soluções apenas para
tempos que satisfazem
1*
tt
. Para isto, o subdomínio
1e;,
2
txtxD
(3.36)
é dividido em três conjuntos disjuntos:
1;,
2
1
txtxD
(3.37)
1;,
2
2
xttxD
(3.38)
xttxD 1;,
2
3
(3.39)
de tal maneira que
321
DDDD , como mostra a Figura 3.8.
Figura 3.8 – Divisão do domínio de interesse em três conjuntos disjuntos
Em seguida, nota-se que se
1
, Dtx
, como
t
x
, então existe
0
tal que
t
x
. Em outras palavras, qualquer ponto nesse conjunto encontra-se sobre uma reta
característica da forma
t
x
, para algum
0
. Visto que as soluções se mantêm
61
constantes sobre as características, pela equação (3.33) tem-se que
1,
0
StxS , para
todo
1
, Dtx
. Por outro lado, desde que
1
xt
, então
11/0 ttx
, para todo
2
, Dtx
. Neste caso, dado
2
, Dtx
, pode-se escrever
ttx 1/
para algum
1,0
. Isto equivale a dizer que qualquer ponto
tx,
de
2
D
encontra-se sobre uma reta
característica da forma
tx 1
, para um dado
1,0
. Logo, a solução é do tipo
ttxStxS 1/11,
0
, ou seja,
txtxS 1/1,
, para todo
2
, Dtx
.
Finalmente, se
3
, Dtx , tem-se que
1
x
. Assim, pela equação (3.35), nota-se que as
curvas características passando por qualquer ponto desse conjunto são retas verticais e, das
relações em (3.33), tem-se que
0,
0
xStxS , para todo
3
, Dtx . Em resumo, verifica-
se que no subdomínio
D
(veja equação 3.36) a solução é a função dada por
3
2
1
,0
,
1
1
,1
,
Dtxse
Dtxse
t
x
Dtxse
txS
(3.40)
Veja uma ilustração dessa solução na Figura 3.9.
Figura 3.9 – Solução em cada subconjunto considerado
62
Apesar da expressão na equação (3.40) representar uma função contínua no
subdomínio
D
, como será visto, para o problema do Exemplo 3.1 não se tem uma solução
contínua para todo
1
t
. Na realidade, será obtida mais adiante uma solução descontínua em
toda a curva no semiplano superior aberto
2
. Para construir essa solução necessita-se
conhecer a curva de descontinuidade, ou seja, conhecer o caminho onde a descontinuidade,
iniciada no ponto
1,1
, se propaga para todo
1
t
.
Tendo como motivação o exemplo acima, serão consideradas soluções que satisfazem
a EDP 0
xt
SSS , no sentido clássico, em cada lado de uma ou mais curvas diferenciáveis
do tipo
tx
, mas que são descontínuas ao longo dessa reta.
Definição 3.3 Diz-se que uma função
2
:S
pertence à classe
K
no semiplano
superior aberto
2
, se
txS ,
é sempre continuamente diferenciável, exceto possivelmente
para um número finito de curvas diferenciáveis em
2
, nos pontos das quais
txS ,
pode ter
descontinuidades. Uma função
txS ,
será chamada de solução fraca do problema de Cauchy
para uma lei de conservação, se
txS ,
pertencente a classe
K
satisfaz a lei de conservação
na sua forma integral, dada pela equação (3.2), e a condição inicial
xSxS
0
0, .
Os seguintes resultados são de demonstrações imediatas.
Proposição 3.5 Suponha que o fluxo f é uma função suficiente diferenciável.
(i) Então toda a solução clássica do problema de Cauchy (3.6) é também uma solução
fraca.
(ii) Se a função
txSS ,
é uma solução fraca do problema de Cauchy e, além disso,
possui derivadas parciais contínuas nas variáveis
x
e
t
, então
txS ,
é uma solução
clássica.
(iii) Se
txSS ,
é uma solução fraca do problema de Cauchy, então
txSS ,
satisfaz o problema descrito na equação (3.6), exceto possivelmente para um número
finito de curvas diferenciáveis em
2
.
Observação 3.1 Como todo conjunto limitado e conexo da reta real é necessariamente um
intervalo fechado, então para problemas unidimensionais, a lei de conservação integral
descrita na equação (3.2) toma a seguinte forma:
63
.,,,,,
batbSftaSfdxtxS
dt
d
b
a
(3.41)
Assim, no presente caso (onde tratamos de problemas em 1D), uma solução fraca
txS ,
do
problema de Cauchy, além de pertencer à classe
K
no semiplano superior aberto
2
,
satisfaz também a relação descrita na equação (3.41).
No que segue,
KtxS
,
representa uma solução fraca do problema descrito na
equação (3.6) e
tx
é uma curva no semiplano superior aberto
2
, nos pontos da qual
txS ,
possui descontinuidades. Seja um intervalo
ba,
do eixo
x
, onde para um certo tempo
t
tem-se
bta
, veja Figura 3.10. Denota-se por
L
S e
R
S os respectivos limites laterais
de
txS ,
, quando
tx
pela esquerda e pela a direita desta curva diferenciável.
Figura 3.10 – Curva de descontinuidade
Tomando a integral de
txS ,
ao longo do intervalo
ba,
tem-se a expressão
b
t
R
t
a
L
b
a
dxtxSdxtxSdxtxStW
,,,
(3.42)
64
a qual está definida para todo
t
com
bat ,
. Aplicando a regra da cadeia, obtemos
.,,,,
dt
d
txSdxtxS
dt
d
txSdxtxS
dt
tdW
R
t
b
t
RL
t
a
L
t
(3.43)
Desde que
x
LL
t
SfS
e
x
RR
t
SfS
, então a equação (3.43) pode ser escrita como:
.,,
b
t
x
R
t
a
x
LR
t
L
dxSfdxSftxStxS
dt
d
dt
tdW
(3.44)
Aplicando o teorema fundamental do cálculo na expressão dentro do parêntese, do lado
direito da equação (3.44), obtém-se
.,,
,,,,
tbSftaSf
txSftxSftxStxS
dt
d
dt
tdW
RLRL
(3.45)
Por outro lado da equação (3.41), pode-se escrever
.,, tbSftaSf
dt
tdW
(3.46)
Substituindo a equação (3.46) em (3.45), obtém-se a seguinte relação:
,f
dt
d
S
(3.47)
onde
LR
SSS
(3.48)
65
e
LR
SfSff
(3.49)
são os respectivos saltos que a solução
S
e o fluxo f apresentam ao cruzarem a curva de
descontinuidade
tx
. A expressão descrita na equação (3.47) é denominada de condição
de salto ou condição de Rankine-Hugoniot. Com essa relação é estabelecido o seguinte
resultado.
Proposição 3.6 Se a solução do problema de Cauchy (3.6) apresenta descontinuidade ao
longo de uma curva diferenciável
tx
que se encontra no semiplano superior aberto
2
,
então ela não pode ser de qualquer tipo. De fato, tal solução fraca deve necessariamente
relacionar os saltos
S
e
f
com a chamada velocidade de propagação da descontinuidade,
dada por
dttdv
s
/
*
. Além disso, temos a seguinte relação, conhecida como condição de
Rankine-Hugoniot:
.
*
LR
LR
Sw
SS
SfSf
S
f
dt
d
v
(3.50)
A classe de soluções fracas que obedecem à relação (3.50) é ainda relativamente grande, de
forma que em geral não há unicidade. Para garantir a unicidade, é necessário impor mais uma
condição, conhecida como princípio da entropia, que será vista mais adiante.
De volta ao Exemplo 3.1, pode-se agora aplicar a condição de Rankine-Hugoniot para
obter uma solução descontínua da equação de Burgers, com a condição inicial mostrada na
equação (3.33). Como 1
S , para todo 1
tx , e 0
S , para todo xt
1 , então é natural
procurar uma solução onde 1
L
S e 0
R
S , para todo
1
t
. Feito isso, tem-se
1S
,
2/1f
e consequentemente da equação (3.47) obtém-se
2/1/ dttd
. Em vista disso, a
equação da curva de descontinuidades pode ser calculada por integração, sendo da forma
ctetx 2/1
. Como esta curva começa no ponto
1,1
, a constante de integração é
imediatamente determinada, obtendo-se a reta
2/1 tx
, a qual é a curva de
descontinuidades deste problema. Assim, para
1
t
, chega-se à solução.
66
2/1,0
2/1,1
,
tx
tx
txS
(3.51)
Finalmente, a solução completa do problema descrito para o Exemplo 3.1 é a união
das expressões definidas na equação (3.40), para
1
t
, com as relações mostradas na equação
(3.51), para
1
t
, veja a Figura 3.11. Como se pode observar, ao contrário da Figura 3.7, a
ilustração na Figura 3.11 enfatiza a interceptação das características apenas na curva de
descontinuidade.
Figura 3.11 – Esboço da solução fraca do problema do Exemplo 3.1
Definição 3.4 Denota-se por
21
C o conjunto das funções
2
:
continuamente
diferenciáveis no semiplano aberto
2
. O suporte de
, indicado por
sup
, é o conjunto
dos pontos
2
,tx para os quais
0, tx
. Diz-se que a função
tem suporte compacto
em
2
, se existem
,bea
e
0
T
tal que
2
,0,sup Tba
, mas
0,,,
tbtaTx
, para todo
x
e
0
t
.
Assim, como indicado na Figura 3.12, deve-se notar que o
sup
está estritamente
contido no conjunto compacto
TbaD ,0,
. Denota-se por
21
0
C
o conjunto de todas
as funções
21
C
que possuem suporte compacto em
2
. Uma tal função
continuamente diferenciável de suporte compacto é denominada de função teste.
67
Figura 3.12 – Suporte da função teste
Proposição 3.7 Suponha que o fluxo f é suficientemente diferenciável e a função
S
pertence à classe
K
em
2
. Se
xxSxS ,0,
0
, então as seguintes
afirmações são equivalentes:
(i)
txSS ,
é uma solução fraca do problema de Cauchy (3.6);
(ii)
0,
dttxSfSdx
, para toda curva fechada
2
e seccionalmente
diferenciável que possui no máximo um número finito de pontos de descontinuidades
de
txSS ,
; (3.52)
(iii)
00,,,,,
0
0
dxxxSdxdttxtxSftxtxS
xt
, para toda função
teste
21
0
C
.
Demonstração: Inicialmente, será mostrado que )()( iii
. Para isso, seja
txSS ,
uma
solução fraca do problema de Cauchy e uma curva fechada no semiplano superior aberto
2
,
como mostra a Figura 3.13. Considera-se que
intR
é o subconjunto de
2
constituído da união de
com o seu interior, denotado por
int
. Desde que a curva
tenha no máximo um número finito de pontos de descontinuidades de
txS ,
, então
68
certamente
não é uma curva característica, o que faz sentido considerar a integral de
txS ,
ao longo de
. Por outro lado, a região
R
pode possuir um número finito de curvas
características, então a função definida pela relação
txStxSftxStx
xt
,,,,
é nula
em
R
, a menos, possivelmente, de um conjunto de medida nula onde ela é descontínua. Logo,
a integral dupla
R
xt
dxdttxStxSftxS ,,,
existe, no sentido de Lebesgue, e é também
nula. Em vista disso, considerando o campo vetorial definido pela a relação
txSftxStx ,,,
, e uma vez que a curva fechada
é seccionalmente diferenciável,
pode-se lançar mão do conhecido Teorema de Green no plano para obter a seguinte
identidade:
,,,,,,
R
tx
dtdxtxStxStxSfdttxSfdxtxS
(3.53)
a qual garante que a integral de linha do segundo tipo do campo vetorial
txSftxS ,,,
é
nula sobre a curva
, ou seja,
0
fdtSdx
. Assim, é provado que )()( iii
. Para mostrar
a implicação )()( iiiii
, supõe-se que
0
fdtSdx
, para toda curva fechada
com as
propriedades descritas no item )(ii deste teorema. Então, pela equação (3.53),
0,,,
dxdttxStxStxSf
R
tx
. Desde que o integrando dessa integral dupla seja
contínuo (a menos de um conjunto de medida nula) e a região
R
possa ser construída de
forma arbitrária, segue-se que
0
xt
SfS (a menos de um conjunto de medida nula).
Então, para toda função teste
21
0
C
, temos
dtdxtxSfSdtdxSfS
T b
a
xtxt
,0
00
Tba ,0,sup
dtdxtxSfdxdttxS
T b
a
x
b
a
T
t
,,
00
dtdxSfSfdxdtSS
T b
a
x
bx
ax
b
a
T
t
Tt
t
00
0
(integração por partes)
69
T b
a
x
b
a
T
t
b
a
dtdxSfdtdxSdxxxS
00
0,0,
0,,,
tbtaTx
dxxxSdtdxSfS
xt
0,
0
0
.
Em outras palavras, obtém-se a equação (3.52), provando que
)()( iiiii
. Resta mostrar que
)()( iiii
. Para isso, dado
0
t
, considera-se
Tba ,0,
. Em seguida, toma-se na
equação (3.52)
, onde por definição
1
, para os pontos pertencentes a
, e
0
,
para pontos que não estão em
. Com essa escolha, da demonstração da implicação anterior,
obtém-se
t b
a
xt
dtdxSfS
0
0
0
0
dxxSdtdxSfS de onde segue que
dtaSfbSfSdx
tb
a
0
, ou seja,
tbSftaSfdxtxS
dt
d
b
a
,,,
. O que
aparentemente conclui a demonstração. Infelizmente
, assim definida, não pode ser aceita
rigorosamente como uma função teste, pois não é diferenciável. Entretanto, de acordo com
Chorin e Marsden (2000), pode-se mostrar que existe uma sequencia de funções teste
n
tal
que
t b
a
xt
t b
a
nxt
dtdxSfSdtdxSfS
00
, quando
n
, o que efetivamente
conclui a demonstração.
Figura 3.13 – Curva fechada em
2
70
3.4 Perda de unicidade e a condição de entropia
O conceito de solução fraca juntamente com a condição de Rankine-Hugoniot
tornaram possível resolver o problema de Cauchy para uma equação hiperbólica não linear,
para todo
0
t
, o que não era sempre viável usando-se apenas a classe de soluções clássicas
ou regulares. No entanto, como será visto adiante, o alargamento do conjunto de soluções do
problema (3.6), com a introdução da classe
K
, muitas vezes permite a existência de mais de
uma solução fraca, todas satisfazendo à mesma condição inicial. Este fato pode ser notado no
seguinte exemplo.
Exemplo 3.2: Considere o problema de valor inicial para a equação de Burgers
xxStxS
txSSS
xt
0
0,
0,0
(3.54)
onde a condição inicial é dada por
0,1
0,0
0
x
x
xS
(3.55)
As características deste problema podem ser vistas na Figura 3.14.
0,
0,
set
se
tx
(3.56)
Na região do plano
t
x
que possui curvas características, a solução deste problema é dada
pela equação (3.57).
xt
tx
txS
0,1
0,0
,
(3.57)
71
Figura 3.14 – Retas características do problema do exemplo 3.2
Além disso, para
tx
0
, onde o problema não possui curvas características, é fácil
ver que a função
txtxS /,
satisfaz a EDP 0
xt
SSS . Logo, a seguinte função contínua
xt
txtx
tx
txS
0,1
0,/
0,0
,
(3.58)
é uma solução deste problema de Cauchy, para todo
0
t
e
x
, como mostra a
Figura 3.15.
Figura 3.15 – Uma solução contínua para o problema do exemplo 3.2
72
Por outro lado, partindo novamente da expressão mostrada em (3.56), pode-se
preencher a lacuna de informações, existentes na região sem características, construindo uma
solução descontínua ao longo de uma curva diferencial que inicia no ponto
0,0
e possui
inclinação 2/1/
dtd
. Assim, essa curva é dada por
2/tt
e a solução procurada, a
qual satisfaz a condição de Rankine-Hugoniot
Sfdtd /2/1/
, é claramente
2/,1
2/,0
,
tx
tx
txS
(3.59)
A Figura 3.16 mostra os valores da solução fraca (3.59) no plano
t
x
, observa-se que
agora algumas características partem de um ponto comum na curva de descontinuidade.
Figura 3.16 – Uma solução descontínua para o problema do exemplo 3.2
O Exemplo 3.2 mostra a necessidade de se impor mais condições às soluções fracas
para se obter uma única solução, que seja fisicamente correta. Na teoria das leis de
conservação, a condição que garante a unicidade de uma solução fraca é denominada de
condição de entropia. A condição de entropia pode ser vista como uma condição de
causalidade: o choque é determinado pelo dado inicial e não por eventos futuros. Existem
mais de uma formulação matemática para a condição de entropia; tais formulações são, de
certo modo, equivalentes. A princípio será utilizada uma abordagem puramente intuitiva.
Observando as três soluções mostradas nas Figuras 3.11, 3.15 e 3.16, nota-se que para cada
ponto no caminho do choque é possível traçar duas características, uma de cada lado do
choque, como nas Figuras 3.11 e 3.16, ou ambas as características podem ser traçadas para o
73
passado, ou seja, para a linha inicial, como na solução mostrada na Figura 3.15. Em vista
disso, seguindo o enfoque de Lax (1957), diz-se que uma solução fraca é admissível, isto é,
satisfaz à condição de entropia, se duas características (mesmo aquelas que se interceptam no
caminho do choque) podem ser traçadas até a linha inicial. Assim, a solução fraca do
Exemplo 3.1, vide a Figura 3.11, é admissível. No entanto, a solução fraca do Exemplo 3.2
esboçada na Figura 3.16 não satisfaz o princípio da condição de entropia. A solução
admissível desse exemplo é o leque de rarefação mostrado na Figura 3.15. Uma solução
fisicamente incorreta, como aquela exibida na Figura 3.16, é algumas vezes denominada de
choque de rarefação. Nota-se que uma tal solução não depende continuamente do dado inicial.
De fato, em uma vizinhança de um choque de rarefação é possível especificar
S
à direita e à
esquerda da curva
tx
, independente do valor
0
0, SxS , a única exigência é o
comprimento da condição de Rankine-Hugoniot.
No que segue, serão fornecidas expressões que caracterizam matematicamente a
condição de entropia. Para isso, seja
*
/
S
vdtdx
a velocidade de uma curva de
descontinuidade da solução do problema de Cauchy, onde, de acordo com a Proposição 3.6,
tem-se
LRLR
S
SSSfSfv /
*
. Supõe-se que a solução descontínua seja admissível,
no sentido da condição de entropia de Lax. Então, na vizinhança do choque deve-se observar
o aspecto geométrico esboçado na Figura 3.17, onde em cada lado do choque, as
características penetram no caminho do choque, de modo que à esquerda do choque ocorre
*
/
S
L
vSfdtdx
e a direita tem-se
*
/
S
R
vSfdtdx
. Isto permite formular a condição
de entropia como segue:
.
* R
S
L
SfvSf
(3.60)
74
Figura 3.17 – Uma típica descontinuidade admissível
A condição descrita na equação (3.60) afirma que, em uma descontinuidade
admissível, a velocidade da curva de descontinuidade é remediada pelas velocidades das
soluções à direita e à esquerda do caminho do choque. Muitas vezes resume-se a equação
(3.60) na forma
RL
SfSf
, a qual é denominada de desigualdade (de entropia) de Lax.
Proposição 3.8 Seja
S
uma solução descontínua do problema de Cauchy. Se f é uma
função convexa, isto é, se 0
f , então as seguintes afirmações são equivalentes:
(i)
R
S
L
SfvSf
*
;
(ii) Existe uma constante 0
c tal que
t
c
a
txStaxS
,,
para todo 0
a , 0
t e
,x
.
Demonstração: )()( iii
. Se ocorrer a desigualdade de entropia de Lax,
RL
SfSf
, e
se 0
f , então necessariamente temos que
RL
SS
(3.61)
75
ou seja, ao se cruzar a curva de descontinuidade da esquerda para direita, experimenta-se um
salto para baixo. Consequentemente, dado
0
t
, fazendo-se
x
variar de
até
, para
todo
0
a
, existe
0
c
tal que o quociente
atxStaxS /,,
nunca superará o valor
tc/
.
)()( iii
. Reciprocamente, supõe-se que existe
0
c
tal que a desigualdade do segundo item
da Proposição 3.8 ocorre, para todo
0
a
,
0
t
e
,x
. Se
S
é uma solução
descontínua do problema de Cauchy, então a desigualdade do segundo item dessa Proposição
implica que
S
experimentará um salto para baixo ao cruzar a curva de descontinuidade, isto
é,
RL
SS . Seja
LRLR
S
SSSfSfv /
*
. Então pelo teorema do valor médio, existe
S
, com
RL
SSS
, tal que
SfSSSfSfv
LRLR
S
/
*
. Como 0
f , então
f
é estritamente crescente e, portanto, tem-se
LR
SSfSf
, ou seja,
L
S
R
SfvSf
*
.
A desigualdade exposta no segundo item da Proposição 3.8 é também denominada de
desigualdade de entropia, sendo uma formulação mais apropriada para demonstrações. Note
que o salto para baixo, descrito na equação (3.61), efetivamente ocorre na solução
descontínua e admissível do Exemplo 3.1, veja Figura 3.11.
Definição 3.5 Diz-se que
S
pertencente à classe
K
em
2
é uma solução entrópica do
problema de Cauchy (3.6), se:
(i)
00,,,,,
0
0
dxxxSdxdttxtxSftxtxS
xt
, para toda função
teste
21
0
C
, e
(ii) existe
0
c
tal que
t
c
a
txStaxS
,,
, para todo
0
a
,
0
t
e
,x
.
A unicidade de uma solução fraca que satisfaz à condição de entropia de Lax é
formulada no seguinte teorema.
Teorema 3.1 (Unicidade com a condição de Lax) Se f é convexa, ou seja, 0
f , então, a
menos de um conjunto de medida nula, uma solução entrópica do problema de Cauchy (3.6) é
única.
76
Demonstração: Veja Smoller (1994)
Se o fluxo não é convexo, então a condição de entropia de Lax não garante a unicidade
de uma solução fraca do problema de Cauchy. Entretanto, uma condição adicional foi
descoberta pela matemática russa Olga A. Oleinik, em 1959. Existe uma tradução em língua
inglesa desse trabalho clássico, que inclui a referida condição, veja Oleinik (1963).
Definição 3.6 Diz-se que
S
pertencente à classe
K
em
2
é uma solução entrópica do
problema de Cauchy (3.6) no sentido de Oleinik, se:
(i)
0,
dttxSfdxS
, para toda curva fechada
2
e seccionalmente diferenciável
que possui no máximo um número finito de pontos de descontinuidade de
txSS ,
.
(ii)
LR
LR
R
R
SS
SfSf
SS
SfSf
, para todo
S
entre
R
S
e
L
S
.
A desigualdade descrita no item (ii) da Definição 3.6 é denominada condição
de entropia de Oleinik. A unicidade de uma solução fraca satisfazendo essa condição de
entropia é formulada no teorema adiante.
Teorema 3.2 (Unicidade segundo Oleinik) A menos de um conjunto de medida nula, uma
solução entrópica do problema de Cauchy (3.6) no sentido de Oleinik é única.
Demonstração: Veja Oleinik (1963).
Observação 3.2 A condição de entropia de Oleinik permite uma interpretação geométrica. A
corda
Sl
RL,
, unindo os pontos
LL
SfS , e
RR
SfS , , dada por
R
LR
LR
R
RL
SS
S
S
SfSf
SfSl
,
, satisfaz às seguintes condições:
(i)
SfSl
RL
,
, se, e somente se,
RL
SS ,
(ii)
SfSl
RL
,
, se, e somente se,
RL
SS
.
77
Em outras palavras, uma descontinuidade admissível ocorre na presença de duas
situações geométricas típicas, que dependem do sinal de
RL
SS
:
(i)
RL
SS
; neste caso, uma descontinuidade é admissível se, e somente se, o gráfico
de f , restrito ao intervalo ],[
LR
SS , encontra-se situado abaixo da corda
Sl
RL,
.
(ii)
RL
SS ; neste caso, uma descontinuidade é admissível se, e somente se, o gráfico
de f , restrito ao intervalo ],[
LR
SS , encontra-se situado acima da corda
Sl
RL,
.
A Figura 3.18 exibe possíveis descontinuidades admissíveis. Descontinuidades não
admissíveis são mostradas na Figura 3.19.
Figura 3.18 – Descontinuidades admissíveis
Figura 3.19 – Descontinuidades não admissíveis
Finalmente, observa-se que, se f é convexa, isto é, se 0
f , então as desigualdades
de entropia de Lax e de Oleinik são equivalentes. Este é o resultado da seguinte proposição,
de fácil demonstração.
78
Proposição 3.8 Se 0
f , então as seguintes desigualdades são equivalentes:
(i)
R
S
L
SfvSf
*
.
(ii)
LR
LR
R
R
SS
SfSf
SS
SfSf
, para todo
S
entre
R
S
e
L
S
.
3.5 Solução analítica para o problema de Burgers
Nesta seção, veremos algumas soluções analíticas para dados iniciais específicos para
o problema de Burgers, com o objetivo de comparar com as soluções numéricas oriundas dos
métodos considerados aqui.
Para isso vamos obter a equação das retas características, sabendo que elas são
constantes quando
0
SS
w
. Para obter uma equação para as características, avaliamos
x
em
função de
t
. Da equação (2.74) e chamando
www
SaSf
obtemos
.
0
Sa
dt
dx
(3.62)
Logo, a equação das retas características pode ser dada por:
.
000
xttSax
(3.63)
Para a equação de Burgers, temos a seguinte equação para as características:
000
0
0
xttSx
SS
S
dt
dx
e
dt
dS
w
w
w
(3.64)
3.5.1 Dado inicial 1
Para o primeiro problema utilizaremos o seguinte dado inicial:
79
.0
00
0
2
1
0,
t
x
x
xS
w
(3.65)
Figura 3.20 – Função de fluxo de Burgers para o dado inicial 1
As retas características são dadas por
)(
)(0
00
00
0
teixodopartindottSx
xeixodopartindoxtSx
SS
S
dt
dx
e
dt
dS
w
w
w
(3.66)
onde
0,000
0,02/1
2
1
txparaS
xtparaS
dt
dx
w
w
(3.67)
Integrando a equação (3.67), obtemos
0,0
0,0
2
1
0
0
txondex
xtondett
dt
dx
(3.68)
80
Para calcularmos a velocidade da descontinuidade no ponto de choque, utilizaremos a relação
de salto de Rankine-Hugoniot, como mostra a equação abaixo:
.
4
1
0
2
1
0
2
1
2
1
2
1
2
2
*
w
w
s
S
f
v
(3.69)
Logo, 1/4 será a velocidade no choque, e a solução analítica para esse problema é dada por:
txtem
txtem
txS
w
4/1,00
4/1,0
2
1
,
(3.70)
cujo as retas características e o esboço da solução analítica podem ser vistos nas Figuras 3.21
e 3.22.
Figura 3.21 – Retas características para o problema de Burgers para o dado inicial 1.
81
Figura 3.22 – Gráfico da solução analítica do problema de Burgers para condição inicial 1
3.5.2 Dado inicial 2
Utilizaremos também o seguinte dado inicial para o segundo problema de Burgers:
00
01
2
1
10
0,
x
x
x
xS
w
(3.71)
As retas características são dadas por:
00
01
2
1
10
,0
0
0
0
0
x
x
x
S
dt
dx
dt
dS
w
(3.72)
isto é
0 1 2 3 4 5 6 7 8
0
0.1
0.2
0.3
0.4
0.5
x
Sw
Solução Analítica para o problema de Burgers
t=8
t=16
t=24
82
00para00
001para2/1
2
1
01para00
txS
txS
txS
dt
dx
w
w
w
(3.73)
Integrando a equação (3.73), obtemos
00onde
001onde
2
1
01onde
0
0
0
txx
txtt
txx
x
dt
dx
(3.74)
A velocidade no ponto de choque pode ser dada por;
.
4
1
2
1
0
2
1
2
1
0
2
1
2
2
*
w
w
s
S
f
v
(3.75)
Ligando os pontos
0,0 e
2/1,0 aparece à linha reta que representa o choque. Observando
a Figura 3.20 podemos ver que a curva que tem segmento incialmente em 1/2 indo para a
origem abaixo da reta que direciona o choque, representará um leque ou rarefação. O choque
inicial para esse problema está em
0
x
com velocidade igual a
4/1
. Logo, podemos obter a
equação da reta que representa o choque, como mostrado abaixo:
.
4
1
4
1
0 ttx
s
(3.76)
Para satisfazer a condição de entropia temos que preencher o vazio com as retas
características da rarefação, como é dado em (Sweby, 2002).
.2/1
~
0,
~
1
1~
AtAx
t
x
A
dt
dx
(3.77)
83
Parte da expansão (em
2/1
~
A
)atinge o choque quando
.
4
1
2
1
1 tt
(3.78)
Do lado esquerdo da equação (3.78) obtemos o leque da rarefação, que é devido ao inicio da
expansão da onda com 2/1
w
S a 1
x antes da ocorrência dos choques. Do lado direito
foi obtido o choque. Resolvendo a equação (3.78) obtemos
.41
2
1
4
1
tt
(3.79)
Entretanto, quando
4
t
temos 0
w
S para ambos os casos quando
1
x
e para
0
x
a
direita do choque,
2/1
w
S
a esquerda do choque enquanto na região tx
1 temos
txtxS
w
/1, . Logo, para
4
t
, a solução analítica para o referido dado inicial é:
xt
tx
t
t
x
t
x
x
txS
w
4
1
0
4
1
1
22
1
1
2
1
1
10
,
(3.80)
Figura 3.23 – Retas características para o problema de Burgers para o referido dado inicial 2.
84
Figura 3.24 – Gráfico da solução analítica do problema de Burgers para o dado inicial 2
3.5.3 Dado inicial 3
O terceiro problema possui o seguinte dado inicial:
01
01
0,
x
x
xS
w
(3.81)
A velocidade no ponto do choque ou na descontinuidade é:
.0
11
2
1
2
1
*
w
w
s
S
f
v
(3.82)
Portanto, como foi visto, por exemplo, em 3.5.1 e 3.5.2, a solução analítica para esse dado
inicial compreende uma onda de rarefação:
-2 -1.5 -1 -0.5 0 0.5 1 1.5 2
0
0.1
0.2
0.3
0.4
0.5
x
Sw
Solução Analítica para o problema de Burgers- Dado inicial (2)
t=3
85
xt
txttx
tx
txS
w
,1
,/
,1
,
(3.83)
Como já foi mencionado, os problemas considerados nesta seção foram introduzidos
no presente trabalho para verificar o desempenho dos métodos numéricos em estudo, os quais
são apropriados para equações diferenciais parciais hiperbólicas. A figura abaixo mostra a
solução analítica (3.83) para t=1/4.
Figura 3.25 – Gráfico da solução analítica do problema de Burgers para o dado inicial 3
3.6 Análise do problema de Riemann
Denomina-se de problema de Riemann aquele que trata de encontrar uma solução
fraca do problema de Cauchy, para uma lei de conservação, com condições iniciais
descontínuas da forma:
0se;
0se;
0,
0,0
xS
xS
txS
txSfS
R
L
x
t
(3.84)
-0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1
-1
-0.8
-0.6
-0.4
-0.2
0
0.2
0.4
0.6
0.8
1
x
Sw
Solução Analítica para o problema de Burgers Dado inicial (3)
t=1/4
86
onde
L
S e
R
S são constantes.
Para analisar o problema de Riemann, supõe-se inicialmente o caso particular onde
0
f . Neste caso, f
é uma função estritamente crescente. Assim, se
RL
SS , então
RL
SfSf e a informação à esquerda se propaga mais lentamente do que a informação à
direita. A condição de entropia de Lax, a qual exige que as características possam ser traçadas
para o passado, garante que a solução não é uma onda de choque. De fato, nesse caso, a
condição de entropia dita que a única solução fisicamente correta do problema de Riemann é
uma onda de rarefação, veja a Figura 3.26. Se
RL
SS
, então
RL
SfSf
e a informação
à esquerda se propaga mais rapidamente do que a informação à direita. A condição de
entropia garante que a solução é uma onda de choque, cuja descontinuidade viaja com
velocidade
LRLR
S
SSSfSfv /
*
, veja a Figura 3.27.
Figura 3.26 – Solução do problema de Riemann, f convexa e
RL
SS
Figura 3.27 – Solução do problema de Riemann, f convexa e
RL
SS
87
Nos casos mais gerais, onde f não é estritamente convexa, a condição de entropia de
Oleinik determina o tipo de solução admissível do problema de Riemann. Considera-se, por
exemplo, que o fluxo f possui uma forma sigmoidal com um único ponto de inflexão. Neste
caso, a onda conectando os estados
L
S e
R
S tanto pode ser um choque como uma rarefação,
ou uma combinação das duas. De fato, se a corda unindo
LL
SfS , e
RR
SfS , não
intercepta totalmente a curva f , então, entre estes dois pontos, a referida curva possui
concavidade voltada para cima ou para baixo. Se
RL
SfSf , então a condição de
entropia de Oleinik garante que a solução é uma onda de choque. Observa-se duas possíveis
ocorrências esboçadas na Figura 3.28. Por outro lado, se
RL
SfSf
, então a solução
conectando os estados
L
S
e
R
S
é uma onda de rarefação, como evidenciado nos dois casos
mostrados na Figura 3.29.
Figura 3.28 – Choque como solução do problema de Riemann, f sigmoidal
88
Figura 3.29 – Rarefação como a solução do problema de Riemann,
f
sigmoidal
Supondo-se que a corda unindo os pontos
LL
SfS , e
RR
SfS , intercepta
totalmente o gráfico de f , existe um ponto
*S
, entre
R
S e
L
S , tal que a referida corda é
tangente a f em
**
, SfS . A solução conectando os estados
R
S
e
*S
é uma onda de
choque e a solução conectando os estados
*S
e
L
S
é uma onda de rarefação. A Figura 3.30
ilustra uma rarefação seguida por uma onda de choque, para o caso
LR
SS
. Se
RL
SS
,
tem-se, de acordo com a Figura 3.31, um choque seguido por uma onda de rarefação.
Figura 3.30 – Rarefação seguida de um choque, f sigmoidal
89
Figura 3.31 – Um choque seguido de rarefação, f sigmoidal
Observação 3.3 É evidente que o presente estudo engloba a análise do problema de Buckley-
Leverett, definido pela a equação (2.70), com o modelo de Stone para a função de fluxo,
equação (2.71). Esse é um caso particular do problema de Riemann, com
LR
SS
. Para
resumir, observa-se que (em vista do exposto acima) a solução do problema de Buckley-
Leverett pode apresentar três formas diferentes: uma onda de choque, evidenciada por uma
das possibilidades indicada na Figura 3.28, uma onda de rarefação, de acordo com um dos
casos mostrados na Figura 3.29, ou uma rarefação seguida de um choque, do tipo exibido na
Figura 3.30.
3.7 Teorema de existência e unicidade de Kruzkhov
O estudo de modelos matemáticos para problemas físicos mostra que leis de
conservação são, muitas vezes, simplificações de uma realidade mais complexa, descrita (por
exemplo) por equações diferenciais da forma
xxStxS
txSSfS
xxxt
0
0,
0,
(3.85)
onde
0
é o coeficiente de difusão, o qual é muito pequeno quando comparado com a
unidade.
90
Proposição 3.9 Supõe-se que exista uma solução
S do problema (3.85), para todo
10
.
Se, a menos de um subconjunto de medida nula de SS
setem,
2
, quando
0
,
então
S
é uma solução fraca do problema de Cauchy (3.6).
Demonstração: Veja Serre (1996)
Considerando a lei de conservação
0
x
t
SfS e integrando sobre um intervalo
xx ,
, obtém-se
.0
x
x
x
x
x
t
dxSfdxS
(3.86)
Supõe-se que
00 f
. Então, se
0S
suficientemente rápido, quando
x
, em vista
da equação (3.86) pode-se afirmar que
.0
dxS
t
(3.87)
Supondo que
S
é diferenciável, a equação (3.87) possui a seguinte forma equivalente
0
Sdx
dt
d
(3.88)
A equação (3.88) mostra que a propriedade
S
é conservada, isto é,
.0,,
dxxSdxtxS
(3.89)
A relação obtida acima leva ao seguinte questionamento: em que condições a equação
(3.89) implica na existência de uma lei de conservação adicional? Isto é, existe uma equação
da forma
91
,0
xt
(3.90)
onde
S
e
S
são funções reais?
Definição 3.7 Quando (3.90) ocorre, diz-se que
S
é a função entropia (ou
simplesmente entropia),
S
é a função fluxo de entropia (ou apenas fluxo de entropia),
e 0
xt
é uma lei de conservação entrópica associada à equação
0
x
t
SfS .
Proposição 3.10 Considere que a equação
0
x
t
SfS possui uma lei de conservação
entrópica associada, 0
xt
, onde
E
é uma função convexa e
0
. Suponha que
S
é
uma solução fraca do problema de Cauchy (3.6), obtido conforme descrito na Proposição 3.9,
onde a convergência
SS
é de tal forma que
S
e a sua derivada com relação à variável
x
convergem limitadamente. Então
S
satisfaz fracamente a seguinte desigualdade:
.0
xt
SS
(3.91)
O que significa que para toda a função teste não negativa
tem-se
.00,0,
0
dxxSxdxdtSS
xt
(3.92)
Demonstração: Veja Smoller (1994) ou Serre (1996).
Observação 3.4 A equação (3.92) pode ser escrita de uma maneira equivalente que facilita a
sua identificação com a forma indicada em (3.91):
.00,0,
0
dxxSxdxdtSS
xt
(3.93)
A desigualdade (3.91), ou mais precisamente a desigualdade (3.93), é denominada de
desigualdade de entropia de Kruzkhov (1969).
92
Definição 3.8 Diz-se que uma solução fraca do problema de Cauchy (3.6),
S
, é uma solução
entrópica no sentido de Kruzkhov, se
S
satisfaz à desigualdade de entropia de Kruzkhov
(3.93) para toda função entropia
contínua, convexa e de fluxo
.
Teorema 3.3 (Kruzkhov 1969) Para toda função
0
S limitada e integrável (a menos de um
conjunto de medida nula) sobre
existe uma, e somente uma, solução entrópica no sentido
de Kruzkhov para o problema de Cauchy (3.6).
Demonstração: Veja Serre (1996).
93
4. ABORDAGENS NUMÉRICAS PARA LEIS DE
CONSERVAÇÃO
4.1 Esquemas conservativos
Considera-se a discretização do domínio espacial em célula
2/12/1
,
iii
xx , de
comprimento uniforme
2/12/1
ii
xxx . O ponto médio
2/
2/12/1
iii
xxx é
denominado nó da célula, enquanto que os pontos
2/1i
x são os bordos de
i
, pois estão
localizados na fronteira. Assim, obtém-se uma grade de células centradas, como mostra a
Figura 4.1.
Figura 4.1 Grade de células centradas
O domínio temporal é discretizado utilizando-se passos de tempo constantes, de
valores
t
, de modo que
tnt
n
representa o nível de tempo
n
, para algum inteiro
0
n
. O
valor médio de
txS ,
na célula
2/12/1
,
iii
xx , em um certo intervalo de tempo
n
tt , é
dado por:
.,
1
,
1
2/1
2/1
i
i
i
dxtxS
x
dxtxS
x
U
n
x
x
n
i
(4.1)
A formulação integral para uma lei de conservação escalar (equação (3.41)) restrita ao
intervalo
2/12/1
,
iii
xx
toma a seguinte forma,
94
.,,,
2/12/1
txSftxSfdxtxS
dt
d
ii
i
(4.2)
Integrando a equação (4.2) ao longo do intervalo de tempo
1
,
nn
tt , obtém-se,
11
,,,,
2/12/11
n
n
n
nii
t
t
i
t
t
inn
dttxSfdttxSfdxtxSdxtxS
(4.3)
dividindo a equação (4.3) por
x
e usando a definição dada em (4.1), pode-se escrever,
.,,
1
11
2/12/1
1
n
n
n
n
t
t
i
t
t
i
n
i
n
i
dttxSfdttxSf
x
UU
(4.4)
Em seguida, define-se o fluxo médio no intervalo de tempo
1
,
nn
tt
, para um dado bordo
2/1i
x , como sendo,
1
,
1
2/12/1
n
n
t
t
i
n
i
dttxSf
t
(4.5)
Combinando as equações (4.4) e (4.5), obtém-se a forma geral do esquema explícito
pretendido,
n
i
n
i
n
i
n
i
x
t
UU
2/12/1
1
(4.6)
a qual exige um dado inicial da forma
0
i
U
, para todo
Zi
.
O número de pontos nos quais os fluxos numéricos
n
i 2/1
são avaliados dependem da
qualidade da solução numérica pretendida. Um esquema da forma mostrada na equação (4.6)
é dito conservativo, se existem inteiros
p
e
q
de modo que ocorrem as seguintes
95
dependências funcionais:
n
qi
n
pi
n
i
UKU
,,
2/1
e
n
qi
n
pi
n
i
UKU
112/1
,,
. Por exemplo,
se 0
p e 1
q , então se tem
n
i
n
i
n
i
UU
12/1
,
e
n
i
n
i
n
i
UU ,
12/1
. Neste caso, o
esquema conservativo é dito de três pontos.
Dependendo da escolha de
n
i 2/1
(os fluxos numéricos nos bordos das células) vários
esquemas conservativos podem ser construído a partir da equação (4.6), onde a solução
procurada
1n
i
U
é o valor médio da função
S
na célula
2/12/1
,
iii
xx , no nível de tempo
futuro
1
n
.
Esses esquemas conservativos são na verdade esquemas numéricos desenvolvidos para
resolver leis de conservação do tipo hiperbólico. Para isto, às vezes, se utilizam os métodos
numéricos de primeira e segunda ordens, como os esquemas de Lax-Friedrichs e Lax-
Wendroff, respectivamente, os quais serão vistos com detalhes mais adiante. Os métodos de
primeira ordem têm a vantagem de aproximar a solução de leis de conservação sem oscilações
ao redor da descontinuidade. Entretanto necessitam de uma malha muito refinada para
produzir uma solução aproximada que seja satisfatória. Por outro lado, os esquemas de
segunda ordem convergem para a solução buscada em uma malha menos refinada, porém pelo
seu caráter oscilatório faz com que as soluções obtidas apresentam severas oscilações ao redor
das ondas de choque. Os métodos de alta resolução são esquemas numéricos que unificam em
um só algoritmo os benefícios dos métodos de primeira e segunda ordem, de tal maneira que
aproximam melhor a solução nas aproximidades das descontinuidades, enquanto que nas
regiões suaves a solução se comporta como a de um método de segunda ordem.
4.2 A condição CFL
O domínio de dependência de um esquema numérico é, por definição, o conjunto de
pontos onde a condição inicial pode afetar a solução numérica em um ponto
ni
tx , . A Figura
4.2 ilustra esse fato para um esquema de três pontos, mostrando que a solução em um
determinado ponto
2
,tx
i
depende das informações contidas nos pontos iniciais
02
,tx
i
,
01
,tx
i
,
0
,tx
i
,
01
,tx
i
e
02
,tx
i
. Em outras palavras, o domínio de dependência
numérico relativo ao ponto
2
,tx
i
é o intervalo
22
,
ii
xx . De forma análoga, pode-se notar
que o domínio de dependência numérico para um ponto genérico
1
,
ni
tx é o intervalo
11
,
nini
xx .
96
Figura 4.2 – Domínio de dependência para um esquema de três pontos
O domínio de dependência teórico é constituído pelos pontos através dos quais a
condição inicial pode afetar a solução
txS ,
da equação diferencial. Sabe-se que para uma
lei de conservação escalar a solução em um ponto
tx,
é determinada pela condição inicial
definida no “pé” da curva característica passando por
tx,
,
txSftx
tx 0,
. Como
S
se mantém constante ao longo de uma característica, então para o referido ponto
tx,
(para uma solução suficientemente diferenciável) tem-se
xStxS
0
, . Logo, de forma
equivalente (para uma solução suave) pode-se reescrever a expressão que define
tx,
como
ttxSftx
tx
,
,
.
Diz-se que a condição CFL (Courant-Friedrichs-Lewy) é satisfeita se o domínio de
dependência teórico estiver contido no domínio de dependência numérico. Assim, como
esboçado na Figura 4.3, considerando o “pé” da característica que passa pelo ponto
1
,
ni
tx
como sendo
111,
,
nniini
ttxSfx
, a condição CFL é satisfeita se,
111,
,
ninini
xx
, para todo
i
e
0
n
. (4.7)
97
Figura 4.3 – Representação geométrica da condição CFL.
Usando a notação
1
1
,
ni
n
i
txSff , pode-se obter uma forma mais apropriada
para a condição CFL. Para isso, escreve-se a relação em (4.7) como segue:
0e111
1
nixnitnfxixni
n
i
(4.8)
ou seja,
0e11
1
nif
x
t
n
i
(4.9)
isto é,
0e1
1
nif
x
t
n
i
(4.10)
Na prática, considera-se a condição CFL (equação (4.10)) na seguinte forma equivalente:
,0e;1max
,
nif
x
t
n
i
ni
(4.11)
98
onde se costuma-se utilizar a aproximação:
n
i
n
i
n
i
n
i
n
i
UUfff
11
/ , se
0
1
n
i
n
i
UU
.
A condição descrita na equação (4.11) é empregada em métodos numéricos como uma
condição necessária, mas não suficiente, para a estabilidade de esquemas explícitos destinados
à resolução de equações hiperbólicas, veja Thomas (1999) e LeVeque (2000), por exemplo.
Em geral, escolhido o valor
x
, a condição CFL é usada para estimar o tamanho do passo de
tempo
t
, de acordo com a relação (4.11).
4.3 Método de Lax-Friedrichs
O esquema numérico de Lax-Friedrichs é de primeira ordem no tempo e no espaço, e
sua região de estabilidade deve satisfazer a relação (4.11). Para esse esquema, o fluxo
numérico é definido pela a expressão:
,
2
2
1
112/1
n
i
n
iii
n
i
UU
t
x
ff
(4.12)
onde
n
i
n
i
Uff
2/12/1
é a função de fluxo. Substituindo a função de fluxo de Burgers,
(equação (2.94)) na equação (4.12) obtém-se
n
i
n
i
nn
i
n
i
UU
t
x
UU
i
1
22
12/1
22
1
2
1
2
1
(4.13)
e
.
22
1
2
1
2
1
1
22
2/1
1
n
i
n
i
nn
i
n
i
UU
t
x
UU
i
(4.14)
Substituindo as equações (4.13) e (4.14) na equação (4.6), obtém-se o esquema conservativo
de diferenças finitas de Lax-Friedrichs para a solução da equação de Burgers, como mostra a
seguir:
99
.
22
1
2
1
2
1
22
1
2
1
2
1
1
22
1
1
22
1
1
n
i
n
i
n
i
n
i
n
i
n
i
n
i
n
i
n
i
n
i
UU
t
x
UU
UU
t
x
UU
x
t
UU
(4.15)
A equação (4.15) pode ser simplificada, ficando da seguinte forma:
.
2
1
4
11
2
1
2
1
1 n
i
n
i
n
i
n
i
n
i
n
i
UUUU
x
t
UU
(4.16)
De forma análoga, substituindo o modelo da função de fluxo (equação (2.71)) do
problema de Buckley-Leverett na equação (4.12) obtém-se
n
i
n
i
n
i
o
w
n
i
n
i
n
i
o
w
n
i
n
i
n
i
UU
t
x
UU
U
UU
U
1
22
2
2
1
2
1
2
1
2/1
2
11
2
1
(4.17)
e
.
2
11
2
1
1
2
1
2
1
2
1
22
2
2/1
n
i
n
i
n
i
o
w
n
i
n
i
n
i
o
w
n
i
n
i
n
i
UU
t
x
UU
U
UU
U
(4.18)
Substituindo as equações (4.17) e (4.18) na equação (4.6), obtém-se o esquema conservativo
de diferenças finitas de Lax-Friedrichs para a solução da equação de Buckley-Leverett
100
n
i
n
i
n
i
o
w
n
i
n
i
n
i
o
w
n
i
n
i
n
i
n
i
n
i
o
w
n
i
n
i
n
i
o
w
n
i
n
i
n
i
n
i
UU
t
x
UU
U
UU
U
UU
t
x
UU
U
UU
U
x
t
UU
1
2
1
2
1
2
1
22
2
1
22
2
2
1
2
1
2
1
1
2
11
2
1
2
11
2
1
(4.19)
A equação (4.19) pode ser simplificada, ficando da seguinte forma:
.
2
1
11
2
11
2
1
2
1
2
1
2
1
2
1
2
1
1
1
n
i
n
i
n
i
o
w
n
i
n
i
n
i
o
w
n
i
n
i
n
i
n
i
UU
UU
U
UU
U
x
t
UU
(4.20)
Em outras palavras, o esquema Lax-Friedrichs pode ser escrito como
.
2
2
1
1111
1 n
i
n
i
n
i
n
i
n
i
n
i
ff
t
x
UUUU
(4.21)
A Figura 4.4 mostra o domínio de dependência para esse esquema. Se a condição
(4.11) concordar com a inclinação de AB então a condição CFL é satisfeita, entretanto a linha
AC é uma violação da condição de CFL, pois se encontra fora do domínio de dependência.
101
Figura 4.4 – Ilustração do domínio de dependência para o esquema Lax-Friedrichs
4.4 Método de Lax-Wendroff
O método de Lax-Wendroff é um esquema de diferenças finitas com acurácia de
segunda ordem no espaço. Como já dito antes, esse método apresenta oscilações espúrias na
vizinhança de um choque. Assim como o esquema Lax-Friedrichs, a região de estabilidade
deve satisfazer a relação (4.11) e o domínio de dependência é análogo do método de primeira
ordem. (Veja Figura 4.4).
O fluxo numérico deste método é escrito como
,
2
1
12/112/1
n
i
n
i
n
i
nn
i
n
i
ffff
i
(4.22)
onde
.
2
2/1
2/1
n
i
n
i
f
x
t
(4.23)
Ao substituir a função de fluxo de Burgers (equação (2.94)) na equação (4.22) obtém-se
22
12/1
22
12/1
2
1
2
1
2
1
2
1
2
1
nn
i
n
i
nn
i
n
i
ii
UUUU
(4.24)
102
e
22
2/1
22
2/1
11
2
1
2
1
2
1
2
1
2
1
nn
i
n
i
nn
i
n
i
ii
UUUU
(4.25)
De forma análoga, é introduzido o modelo da função de fluxo (equação (2.71)) do problema
de Buckley-Leverett na equação (4.22) obtendo-se
22
2
2
1
2
1
2
1
2/1
22
2
2
1
2
1
2
1
2/1
11
11
2
1
n
i
o
w
n
i
n
i
n
i
o
w
n
i
n
i
n
i
n
i
o
w
n
i
n
i
n
i
o
w
n
i
n
i
n
i
UU
U
UU
U
UU
U
UU
U
(4.26)
e
.
11
11
2
1
2
1
2
1
2
1
22
2
2/1
2
1
2
1
2
1
22
2
2/1
n
i
o
w
n
i
n
i
n
i
o
w
n
i
n
i
n
i
n
i
o
w
n
i
n
i
n
i
o
w
n
i
n
i
n
i
UU
U
UU
U
UU
U
UU
U
(4.27)
Agora pode-se escrever o esquema de Lax-Wendroff para a solução da equação de Burgers
.)25.4()24.4(
1
equaçãoequação
x
t
UU
n
i
n
i
(4.28)
103
e para o problema de Buckley-Leverett
.)27.4()26.4(
1
equaçãoequação
x
t
UU
n
i
n
i
(4.29)
Costuma-se empregar a aproximação
2/2/
12/1
n
i
n
i
n
i
UUfxt
. Assim se pode obter
uma formalização mais geral do esquema de Lax-Wendroff, como mostra a seguir,
.
2
22
1
2
1
1
1
1
2
11
1
n
i
n
i
n
i
n
i
n
i
n
i
n
i
n
i
n
i
n
i
n
i
n
i
ff
UU
f
ff
UU
f
x
t
ff
x
t
UU
(4.30)
4.5 Métodos TVDs com limitadores de fluxo
Os esquemas do tipo TVDs (Total Variation Diminishing) são métodos numéricos de
alta resolução baseados no método das diferenças finitas para solução da equação da
conservação. Tais métodos foram desenvolvidos com o objetivo de evitar as oscilações
espúrias tipicamente encontradas nas soluções numéricas fornecidas pelos esquemas
clássicos, como aqueles descritos, por exemplo, nos textos de Aziz e Settari (1979) e
Peaceman (1977). Um método é dito TVD se:
,
1 nn
TVTV
(4.31)
onde
.
1
i
n
i
i
n
i
n
i
n
UUUTV
(4.32)
A equação (4.32) define o conceito de variação total, introduzido por Harten (1983).
104
Segundo Sweby (1984), um esquema TVD possui a notável propriedade de produzir
soluções que obedecem ao princípio de entropia, e que, portanto, são as soluções corretas do
problema.
A discretização para esse tipo de método é da seguinte forma:
.0
2/12/1
1
x
UfUf
t
UU
n
i
n
i
n
i
n
i
(4.33)
Nota-se que a discretização do termo convectivo
xf /
foi obtida através de uma diferença
centrada de segunda ordem. As saturações nos pontos médios
n
i
U
2/1
serão definidas por,
,
2
1
2/1
2/1
n
i
n
i
n
i
n
i
n
i
UU
r
UU
(4.34)
onde a função
n
i
r
2/1
é denominada de limitador de fluxo, as quais dependem do parâmetro
n
i
r
2/1
. Esse parâmetro mede a razão de inclinação da saturação nos bordos da célula
2/12/1
,
ii
xx
, como descreve a equação seguinte.
.,0max
1
1
2/1
n
i
n
i
n
i
n
i
n
i
UU
UU
r
(4.35)
Se
U
é suave, espera-se que o parâmetro
n
i
r
2/1
se aproxime da unidade, por outro lado, na
presença de uma descontinuidade, essa razão deve se afastar do valor unitário. Assim,
utilizando a notação
n
i
n
i
n
i
n
i
n
i
n
i
UU
r
UfUff
1
2/1
2/12/1
2
(4.36)
e a partir das equações (4.33) e (4.34), obtém-se, para todo 1,...,2
Ni a seguinte equação
discretizada
105
.0
2/12/1
1
x
ff
t
UU
n
i
n
i
n
i
n
i
(4.37)
A função de fluxo
n
i
r
2/1
define uma família de esquemas do tipo limitadores de
fluxo. Um membro dessa família difere de outro pela escolha da função
r
. A Tabela 1
mostra alguns limitadores de fluxo.
Os esquemas do tipo TVD’s vem ao longo do tempo obtendo excelentes resultados
quando aplicados as equações de conservação, por isso tais esquemas foram introduzidos
neste trabalho afim de obter uma comparação com outros métodos
Tabela 1 – Limitadores de Fluxo
Função
( )
r
com
0
r
Denominação
( ) 2 1
r r r
Limitador de Van Leer
( ) min 1,
r r
MINMOD
( ) min 2,
r r
Limitador de Osher
( ) min 2 , ( 1) 2, 2
r r r
MUSCL
( ) min 2 , (2 1) 3, 2
r r r
Limitador de Koren
( ) max min 1,2 , min 2,
r r r
SUPERBEE
Os métodos numéricos conservativos utilizados nesse trabalho são todos baseados nos
métodos de diferenças finitas. Outros exemplos de esquemas conservativos podem ser
encontrados nos livros de Thomas (1999) e LeVeque (2002), por exemplo.
4.6 Teorema de Lax-Wendroff
Na presente seção, um esquema numérico para lei de conservação é representado de
uma maneira genérica como segue,
.0
n
i
n
i
UL
(4.38)
106
Definição 4.1 O esquema (4.38) é dito consistente com a lei de conservação
0
xt
SfS ,
no ponto
tx,
, se
0, tnxiSL
n
i
, quando 0,
tx e
txtnxi ,,
.
Pode-se mostrar que o método de Lax-Friedrichs é, por exemplo, consistente com a lei
de conservação escalar
0
xt
SfS .
Teorema 4.1 (Lax-Wendroff) Seja
n
i
U
uma solução numérica do problema de Cauchy (3.6)
obtida através de um esquema conservativo e consistente com a equação
0
xt
SfS . Se
SU
n
i
, quando
0,
tx
, então
txSS ,
é uma solução fraca do problema de Cauchy.
Demonstração: Veja Thomas (1999).
Essencialmente, se a solução de um esquema conservativo converge, então o teorema
de Lax-Wendroff mostra que ela deve convergir para uma solução fraca do problema de
Cauchy para uma lei de conservação. No entanto, observa-se que esse teorema não garante
que
S
é uma solução fraca fisicamente admissível, ou mais precisamente, não garante que
S
seja exatamente a única solução entrópica no sentido do teorema de unicidade de Kruzkhov.
4.7 Desigualdade de entropia discreta
De acordo com o exposto na Seção 3.6, aqui
SEE
denota uma função entrópica e
S
é a função de fluxo de entropia, as quais satisfazem a lei de conservação entrópica
0
xt
associada à lei de conservação
0
xt
SfS
.
No que segue, considera-se a existência de uma função denominada de fluxo numérico
entrópico
n
i 2/1
, obtida a partir dos valores de uma função escalar
definida nos bordos de
cada célula da grade de discretização. Para um esquema conservativo, supõe-se a existência
de inteiros
p
e
q
de modo que,
.,,...
2/1
n
qi
n
pi
n
i
UU
(4.39)
107
Em particular, para um esquema de três pontos, a função
n
i 2/1
é da forma,
.,
12/1
n
i
n
i
n
i
UU
(4.40)
Definição 4.2 Diz-se que a solução
n
i
U
de um esquema numérico satisfaz a desigualdade de
entropia discreta, se existe um fluxo numérico entrópico
n
i 2/1
tal que,
.
2/12/1
1 n
i
n
i
n
i
n
i
x
t
UU
(4.41)
Definição 4.3 Um esquema conservativo (equação 4.6) é dito um E-esquema, se:
112/1
se,
iiiii
UUUUSSf
(4.42)
iiiii
UUUUSSf
112/1
se,
(4.43)
Pode-se mostrar que o método de Lax-Friedrichs é um E-esquema.
Proposição 4.1 Considere a lei de conservação
0
xt
SfS
juntamente com uma função
entropia associada,
, convexa de fluxo
. Suponha que
n
i
U
é uma solução numérica para o
problema de Cauchy (3.6) obtida a partir de um E-esquema. Se
n
i
U
satisfaz a condição
,
12/112/1
n
i
n
i
n
i
n
i
n
i
n
i
UUff
x
t
(4.44)
então
n
i
U
satisfaz também a desigualdade de entropia discreta (equação (4.41)), para algum
fluxo numérico entrópico
n
i 2/1
.
Demonstração: Veja Tadmor (1984).
108
Definição 4.4 Um fluxo numérico entrópico de um esquema conservativo
n
qi
n
pi
n
i
UKU
,,
2/1
é dito consistente com o fluxo de entropia
, se
UUKU
,,
.
Teorema 4.2 Sejam
e
, a função entropia e fluxo de entropia de uma lei de conservação
entrópica associada à lei de conservação
0
xt
SfS
. Seja
n
i
U uma solução numérica do
problema de Cauchy (3.6) obtida a partir de um esquema conservativo. Suponha que exista
um fluxo numérico entrópico
n
i 2/1
consistente com o fluxo de entropia
tal que
n
i
U
satisfaz a desigualdade de entropia discreta (equação (4.41)). Se
SU
n
i
, quando
0,
tx
, então
txS ,
satisfaz a desigualdade de entropia de Kruzkhov.
Demonstração: Veja Thomas (1999).
Nota-se que a Proposição 4.1 e o Teorema 4.2, juntamente com o teorema de Lax-
Wendroff, indicam que a solução de um E-esquema converge para a única solução entrópica
(segundo Kruzkhov) do problema de Cauchy (3.6). Neste sentido, pode-se afirmar que o
método de Lax-Friedrichs, por exemplo, gera (para
x
e
t
suficientemente pequenos)
soluções para o problema (3.6) as quais são compatíveis com o principio de entropia de
Kruzkhov. Infelizmente, o próximo resultado mostra que a classe dos E-esquemas é menor do
que desejado.
Proposição 4.2 Um E-esquema possui no máximo acurácia de primeira ordem.
Demonstração: Veja Osher (1984).
4.8 Esquema composto
Assim como os métodos TVD’s, o esquema composto é um esquema desenvolvido
por Liska e Wendroff (1998) com o objetivo de evitar as oscilações espúrias tipicamente
encontradas nas soluções numéricas fornecidas pelos esquemas de segunda ordem, obtendo
uma solução mais exata na captura do choque, para isso é feita uma composição entre os
métodos de primeira e segunda ordens, difusivos e oscilatórios respectivamente. Os métodos
utilizados nessa composição foram o Lax-Wendroff e Lax-Friedrichs.
109
Para o método de Lax-Wendroff será proposta uma modificação devido à dificuldade
inicial com a exigência do cálculo da derivada da função de fluxo em cada passo de tempo. A
tal modificação é conhecido como método de Lax-Wendroff em duas etapas, esse método
modificado evita essa exigência da derivada. A primeira etapa desse método modificado é de
predição, sendo dada por:
.
2
2
1
11
2/1
2/1
n
i
n
i
n
i
n
i
n
i
UfUf
x
t
UUU
(4.45)
A segunda etapa é de correção e tem a forma do método das diferenças centradas, como se vê
na seguinte equação:
.
2/1
2/1
2/1
2/1
1
n
i
n
i
n
i
n
i
UfUf
x
t
UU
(4.46)
Nota-se que esse método modificado utiliza estratégia das grades deslocadas de forma que os
valores da saturação, representados computacionalmente por
U
, nas bordas das células nos
tempos médios (equação (4.45)) são cálculos intermediários, pois os valores efetivos da
saturação são obtidos em cada nó da grade no tempo 1
n .
A idéia para a estratégia do esquema composto é simples, serão empregados vários
passos do esquema de Lax-Wendroff em duas etapas seguidos de uma aplicação do método de
Lax-Friedrichs. Desta forma, o esquema de Lax-Friedrichs, cuja solução satisfaz o princípio
de entropia, funciona como um filtro para o esquema de segunda ordem aplicado antes. Nessa
composição o esquema de Lax-Friedrichs também trabalha em duas etapas. A primeira etapa
é a mesma da equação (4.45) e a segunda etapa é escrita na seguinte forma:
.
2
2
1
2/112/1
2/1
2/11
2/1
2/1
1
1
n
i
n
i
n
i
n
i
n
i
UfUf
x
t
UUU
(4.47)
Será denotado por
W
L o método em duas etapas (equações ((4.45) e (4.46)) e por
F
L
o
operador correspondente às equações (4.45) e (4.47).
k
o operador definido por
1k
110
etapas de
W
L seguidas por uma aplicação de
F
L
. Logo, a estratégia do esquema composto é
dada como segue:
,...
1
vezesk
WWFk
oLoLoL
(4.48)
onde a solução final é dada por
.
1 n
ik
n
i
UU
(4.49)
O esquema composto definido nas equações (4.48) e (4.49) será denotado por
kLWLF
. Liska e Wendroff (1998) sugerem o valor
4
k
devido aos melhores resultados
encontrados. O esquema composto
kLWLF
é de primeira ordem. Liska e Wendroff
propõem também um esquema composto de segunda ordem para problemas unidimensionais.
Para isso, os autores substituiem na composição (4.49) o operador
F
L
por um filtro de
segunda ordem do tipo WENO, para mais detalhes veja Liska e Wendroff (1998). No presente
trabalho será utilizado o esquema composto
4
LWLF
.
111
5. RESULTADOS
5.1 Aspectos gerais
Neste capítulo serão apresentados os resultados numéricos para as equações de
Burgers e Buckley-Leverett. Ambas as equações estão baseadas no problema de Cauchy para
lei de conservação, como segue:
0,
0,
0,
0,0
xseS
xseS
xS
tx
x
Sf
t
S
R
L
(5.1)
Será adotado o modelo de Stone (equação (2.71)) para a função de fluxo fracionária para o
problema de Buckley-Leverett, o qual depende do valor da razão de viscosidade dos fluidos
presentes.
Em todas as simulações foi tomado
t
suficientemente pequeno, procurando sempre
satisfazer à condição CFL (equação (4.11)), ou seja:
.
max
i
i
f
x
t
(5.2)
Assim verifica-se que a condição CFL depende do valor máximo da derivada da função de
fluxo. Para o problema de Burgers tem-se,
,SSf
(5.3)
112
Figura 5.1 – Gráfico da Função
Sf
para o problema de Burgers
podendo ser visto na Figura 5.1 que
1max
Sf
, para o problema de Burgers.
No problema de Buckley-Leverett, teremos pontos de máximos distintos para cada valor da
razão da viscosidade
óleoágua
c
/
,
.
1
12
2
2
2
ScS
ScS
Sf
(5.4)
Assim, nota-se conforme a Figura 5.2 que
97.2max
Sf
, se
1.0
c
. Na Figura 5.3 tem-se
que
08.2max
Sf para
5.0
c
. Da mesma forma, de acordo com a Figura 5.4,
2max
Sf se
1
c
. Assim, por exemplo, em um reservatório de 10000 cm com
cmx 100
, para
1
c
, o passo de tempo
t
para satisfazer a relação CFL tem que cumprir
50
t
.
0 0.2 0.4 0.6 0.8 1
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
S
Derivada do Fluxo Fracionário (Burgers) x Saturação da água
f'(S)
113
Figura 5.2 – Gráfico da Função
Sf
para o problema de Buckley-Leverett, c=0.1
Figura 5.3 – Gráfico da Função
Sf
para o problema de Buckley-Leverett, c=0.5
0 0.2 0.4 0.6 0.8 1
0
0.5
1
1.5
2
2.5
3
S
Derivada do Fluxo Fracionário (BuckleyLeverett)xSaturação da água
f'(S)
max|f'(S)|=2.97
c=0.1
0 0.2 0.4 0.6 0.8 1
0
0.5
1
1.5
2
S
Derivada do Fluxo Fracionário(Buckley-Leverett)xSaturação da água
f'(S)
maxf'(S)|=2.08
c=0.5
114
Figura 5.4 – Gráfico da Função
Sf
para o problema de Buckley-Leverett, c=1
Por questão de brevidade, será denotado o método de Lax-Wendroff por
LW
, o
método de Lax-Friedricks por
LF
, o método composto será apresentado por
4
LWLF
e o
método TVD pelo seu limitador de fluxo.
5.2 Condições de injeção e extração
No desenvolvimento do presente trabalho, como já dito antes, os fuidos são
incompressíveis e a saturação inicial é dada pela condição de Cauchy no tempo zero.
Para simular um reservatório de petróleo, considera-se que o fluido é introduzido e
retirado da região computacional em um domínio unidimensional. Nos métodos que serão
testados nesse trabalho, o fluxo é introduzido por injeção lateral através da fronteira esquerda
e retirado pela fronteira direita, técnica esta mais adequada para leis de conservação.
5.3 Comparações dos métodos LF, LW e LWLF-4
Para as comparações dos métodos consideraremos uma grade uniforme com 100
células, com os valores iniciais 5.0
L
S e 0
R
S em um domínio unidimensional de 10000
cm. Nas Figuras 5.5 e 5.6 são mostradas comparações dos métodos LW e
LF
, para um
0 0.2 0.4 0.6 0.8 1
0
0.5
1
1.5
2
S
Deriva do Fluxo Fraionário(Buckley-Leverett)xSauração da água
f'(S)
max|f'(S)|=2
c=1
115
tempo
st 900
com distintos valores para a razão de viscosidades e
t
, para os problemas de
Buckley-Leverett e Burgers respectivamente.
Figura 5.5 – Comparação dos esquemas
LW
e
LF
para o problema de Buckley-Leverett com
o passo de tempo
st 5
e razão de viscosidade
5.0
c
.
Figura 5.6 – Comparação dos esquemas
LW
e
LF
para o problema de Burgers para o passo
de tempo
st 20
e razão de viscosidade
1
c
.
1000 2000 3000 4000 5000 6000 7000 8000 9000 10000
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
Comprimento do Reservatório [cm]
Saturação da Água
Problema de Buckley-Leverett - c=0.5
LW
LF
1000 2000 3000 4000 5000 6000 7000 8000 9000 10000
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
Comprimento do Reservatório [cm]
Saturação da Água
Problema de Burgers - c=1
LW
LF
116
Como já foi mencionado antes, a solução do método
LF
apresenta uma considerável difusão
numérica. Já na solução do
LW
podem ser observadas fortes oscilações espúrias, de origem
numérica.
Nas Figuras 5.7 e 5.8 são comparadas as soluções calculadas com os três métodos
numéricos para o problema de Buckley-Leverett, com passos de tempo e valores da razão de
viscosidade distintos para os dois casos. Como era de se esperar, o método
LW
apresenta
oscilações espúrias e o
LF
se mostra bastante difusivo com respeito à solução numérica do
esquema composto
4
LWLF
, o qual se encontra mais próximo da solução descontínua. Para
um maior destaque, o domínio é descrito em termos percentuais do comprimento do
reservatório.
Figura 5.7 – Comparação dos esquemas
LW
,
LF
e
4
LWLF
para o problema de Buckley-
Levertt com o passo de tempo
st 5
no instante
st 360
e razão de viscosidade
1
c
.
0 20 40 60 80 100
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
Porcentagem do Reservatório [%]
Saturação da Água
Saturação da Água x Porcentagem do Reservatório para c=1
LW
LF
LWLF-4
117
Figura 5.8 – Comparação dos esquemas
LW
,
LF
e
4
LWLF
para o problema de Buckley-
Levertt com o passo de tempo
st 5
no instante
st 900
e razão de viscosidade
5.0
c
.
5.4 Aplicação do método LWLF-4 ao problema de Burgers
Veremos o desempenho do esquema composto
4
LWLF
em relação às soluções
analíticas expostas na Seção (3.5), para o problema de Burgers em um domínio
unidimensional com distintos valores para o espaçamento
x
e um valor fixo para o passo de
tempo,
6
10.5
t
. O método mostra-se muito viável para a solução do referido problema.
5.4.1 Resultados para o dado inicial 1
As Figuras 5.9, 5.10 e 5.11 mostram resultados para o problema exposto na Seção 3.51
para um tempo fixo
8.1
t
. Observa-se que quanto mais fina a grade maior é a acurácia do
método composto.
0 20 40 60 80 100
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
Porcentagem do Reservatório [%]
Saturação da Água
Saturação da Água x Porcentagem do Reseratório para c=05
LWLF-4
LF
LW
118
Figura 5.9 – Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 150 células com
0067.0
x
.
Figura 5.10 – Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 350 células com
0029.0
x
.
0 0.2 0.4 0.6 0.8 1
0
0.1
0.2
0.3
0.4
0.5
x
Saturação
Problema de Burgers
LWLF-4
Solução Analítica
0 0.2 0.4 0.6 0.8 1
0
0.1
0.2
0.3
0.4
0.5
x
Saturação
Problema de Burgers
LWLF-4
Solução Analítica
119
Figura 5.11 – Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 450 células com
002.0
x
.
As soluções analíticas para os dados iniciais 2 e 3 compreendem uma onda de
rarefação, como se pode ver nas Seções 3.52 e 3.53.
5.4.2 Resultados para o dado inicial 2
As Figuras 5.12, 5.13, 5.14 e 5.15 mostram resultados para o problema exposto na
Seção 3.52 para um tempo fixo
3
t
. Assim como o problema anterior, observa-se quanto
mais fina a grade maior é a acurácia do método composto.
0 0.2 0.4 0.6 0.8 1
0
0.1
0.2
0.3
0.4
0.5
x
Saturação
Problema de Burgers
LWLF-4
Solução Analítica
120
Figura 5.12 – Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 400 células com
01.0
x
.
Figura 5.13 – Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 800 células com
005.0
x
.
-2 -1.5 -1 -0.5 0 0.5 1 1.5 2
0
0.1
0.2
0.3
0.4
0.5
x
Solução (S)
Problema de Burgers
LWLF-4
S.Analítica
-2 -1.5 -1 -0.5 0 0.5 1 1.5 2
0
0.1
0.2
0.3
0.4
0.5
x
Solução (S)
Problema de Burgers
LWLF-4
S.Analítica
121
Figura 5.14 – Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 1000 células com
004.0
x
.
Figura 5.15 – Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 1200 células com
0033.0
x
.
-2 -1.5 -1 -0.5 0 0.5 1 1.5 2
0
0.1
0.2
0.3
0.4
0.5
x
Solução (S)
Problema de Burgers
LWLF-4
S.Analítica
-2 -1.5 -1 -0.5 0 0.5 1 1.5 2
0
0.1
0.2
0.3
0.4
0.5
x
Solução (S)
Problema de Burgers
LWLF-4
S.Analítica
122
5.4.3 Resultados para o dado inicial 3
As Figuras 5.16, 5.17 e 5.18 mostram resultados para o problema exposto na Seção
3.53 para um tempo fixo
4/1
t
. Assim como os demais problemas, o método composto tem
um melhor desempenho com o refinamento da grade.
Figura 5.16 – Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 300 células com
0067.0
x
.
Figura 5.17 – Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 400 células com
005.0
x
.
-0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1
-1
-0.8
-0.6
-0.4
-0.2
0
0.2
0.4
0.6
0.8
1
x
Solução (S)
Problema de Burgers
LWLF-4
S.Analítica
-0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1
-1
-0.8
-0.6
-0.4
-0.2
0
0.2
0.4
0.6
0.8
1
x
Solução (S)
Problema de Burgers
LWLF-4
S.Analítica
123
Figura 5.18 – Comparação do esquema composto
4
LWLF
com a solução analítica do
Problema de Burgers para uma grade de 600 células com
0033.0
x
.
5.5 Aplicação do Método LWLF-4 no Problema de Buckley-Leverett
Os resultados numéricos encontrados foram comparados com a solução analítica de
Welge (1952) (equação (2.93)) para os valores de saturação de corte para os distintos valores
da razão de viscosidade e do tempo. Os mesmos mostraram-se satisfatórios, pois para todos os
tempos apresentaram poucos pontos são encontrados no salto, porém o valor da saturação no
choque ainda foi um pouco elevado, apesar do método apresentar solução pouca difusiva.
Com o objetivo de encontrar um valor numérico mais próximo do valor analítico, foi
feito um estudo de refinamento da malha. Esse estudo foi realizado primeiramente fixando
x
e variando o passo de tempo
t
, e depois a recíproca, ou seja, variando
x
, sempre
respeitando a condição CFL. Para uma melhor escolha da grade foi levado em conta nesta
análise a menor difusibilidade e a melhor aproximação do valor do choque, sendo
cmx 40
e
st 7
os melhores valores encontrados, em uma grade não muito fina com 250 células
para um domínio unidimensional de 10000 cm. Esses resultados podem ser vistos nas Figuras
5.21, 5.22 e 5.23 com distintos valores para a razão de viscosidade e do tempo, mostrando-se
realmente mais satisfatórios principalmente para a razão de viscosidade mais alta, onde a
solução numérica apresentou menos pontos no salto e melhor valor da saturação no choque.
-0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1
-1
-0.8
-0.6
-0.4
-0.2
0
0.2
0.4
0.6
0.8
1
x
Solução (S)
Problema de Burgers
LWLF-4
S.Analítica
124
Os resultados numéricos foram simulados com os seguintes dados iniciais: 1
L
S e
0
R
S .
Figura 5.19 – Aplicação do método
4
LWLF
para o problema de Buckley-Leverett com
150 células, cmx 67,66
,
st 10
e razão de viscosidade igual a
1
c
.
Figura 5.20 – Aplicação do método
4
LWLF
para o problema de Buckley-Leverett com
150 células, cmx 67,66
,
st 10
e razão de viscosidade igual a
5.0
c
.
0 2000 4000 6000 8000 10000
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
Comprimento do Reservatório [cm]
Esquema Composto LWLF-4
Saturação (S)
T=260s
T=230s
T=200s
T=170s
T=140s
T=110s
T=80s
T=50s
S*=0.707
0 2000 4000 6000 8000 10000
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
Comprimento do Reservatório [cm]
Esquema Composto LWLF-4
Saturação (S)
T=230s
T=200s
T=170s
T=140s
T=110s
T=80s
T=50s
T=20s
S*=0.577
125
Figura 5.21 – Aplicação do método
4
LWLF
para o problema de Buckley-Leverett com
250 células,
cmx 40
,
st 7
e razão de viscosidade igual a
1
c
.
Figura 5.22 – Aplicação do método
4
LWLF
para o problema de Buckley-Leverett com
250 células,
cmx 40
,
st 7
e razão de viscosidade igual a
5.0
c
.
0 1000 2000 3000 4000 5000 6000 7000 8000 9000 10000
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
Comprimento do Reservatório [cm]
Esquema Composto LWLF-4
Saturação (S)
T=260s
T=170s
T=80s
S*=0.707
0 1000 2000 3000 4000 5000 6000 7000 8000 9000 10000
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
Comprimento do Reservatório [cm]
Esquema Composto LWLF-4
Saturação (S)
T=260s
T=170s
T=80s
S*=0.577
126
Figura 5.23 – Aplicação do método
4
LWLF
para o problema de Buckley-Leverett com
250 células,
cmx 40
,
st 7
e razão de viscosidade igual a
1.0
c
.
5.6 Comparações do método LWLF-4 com os métodos tipo TVDs com limitadores de
fluxo na captura do choque
Foi comparado o método LWLF-4 com os métodos do tipo TVD com limitadores de
fluxo. Foi usada também a solução analítica de Welge (1952). Para todas as simulações foram
usadas uma grade uniforme com 350 células, para T=230 e
0028.0
x
com um passo de
tempo muito pequeno
0005.0
t
, a fim de satisfazer a condição CFL para o problema de
Buckley-Leverett. Foram testados três limitadores de fluxo, os mais utilizados em problemas
de petróleo: Superbee, Van Leer e o Minmod, todos obtendo resultados satisfatórios na
captura do choque. Entretanto, o método composto LWLF-4 mostrou mais capacidade para a
captura do choque em relação aos limitadores de fluxo Superbee e Van Leer, que obtiveram
uma pequena elevação para o valor da saturação no ponto de choque, para os diferentes
valores da razão da viscosidade. No entanto o limitador de fluxo Minmod obteve uma grande
compatibilidade com o método LWLF-4, sendo mais robusto que os outros limitadores
quando comparado com a solução analítica de Welge (1952), como mostram os gráficos e
tabelas adiante.
O cálculo para os erros inseridos nas tabelas foi dado pela seguinte relação:
0 2000 4000 6000 8000 10000
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
Comprimento do Reservatório [cm]
Esquema Composto LWLF-4
Saturação (S)
T=210s
T=140s
T=70s
S*=0,302
127
,
*
*
ERRO
Numérico
S
SS
(5.5)
onde
*S
é o valor analítico da saturação no ponto de choque e
Numérico
S o valor numérico.
Figura 5.24 - Comparação do método
4
LWLF
com o limitador de fluxo Superbee, para
1
c
.
Figura 5.25 - Comparação do método
4
LWLF
com o limitador de fluxo Superbee, para
5.0
c
.
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
x
Comparação dos métodos na captura do choque
Saturação
LWLF-4
SUPERBEE-TVD
S*=0.707
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
x
Comparação dos métodos na captura do choque
Saturação
LWLF-4
SUPERBEE-TVD
S*=0.577
128
Figura 5.26 - Comparação do método
4
LWLF
com o limitador de fluxo Superbee, para
1.0
c
.
Figura 5.27 - Comparação do método
4
LWLF
com o limitador de fluxo Van Leer, para
1
c
.
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
x
Comparação dos métodos na captura do choque
Saturação
LWLF-4
SUPERBEE-TVD
S*=0.302
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
x
Comparação dos métodos na captura do choque
Saturação
LWLF-4
Van Leer-TVD
S*=0.707
129
Figura 5.28 - Comparação do método
4
LWLF
com o limitador de fluxo Van Leer, para
1.0
c
.
Figura 5.29 - Comparação do método
4
LWLF
com o limitador de fluxo Minmod, para
1
c
.
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
x
Comparação dos métodos na captura do choque
Saturação
LWLF-4
Van Leer-TVD
S*=0.302
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
x
Comparação dos métodos para a captura do choque
Saturação
LWLF-4
MINMOD-TVD
S*=0.707
130
Figura 5.30 - Comparação do método
4
LWLF
com o limitador de fluxo Minmod, para
1.0
c
.
Tabela 2- Valores da saturação no ponto de choque da comparação dos métodos para a razão
de viscosidade c=1.
Tabela 3- Valores da saturação no ponto de choque da comparação dos métodos para a razão
de viscosidade c=0.1.
0 0.2 0.4 0.6 0.8 1
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
x
Comparação dos métodos para a captura do choque
Saturação
LWLF-4
MINMOD-TVD
S*=0.302
131
6. CONCLUSÕES E PERSPECTIVAS
6.1 Conclusões
Neste trabalho foram considerados os problemas de Burgers e Buckley-Leverett
visando ao escoamento imiscível bifásico água-óleo que ocorre na recuperação secundária de
petróleo.
Como tais problemas são modelados por equações diferenciais parciais hiperbólicas
não lineares, foram considerados os esquemas de diferenças finitas, incluindo um método
composto propriamente elaborado para esses tipos de equações.
Entretanto, o referido trabalho teve também o objetivo de incorporar conhecimentos
matemáticos mais avançados, para em seguida utilizá-lo no entendimento do problema
estudado. Abordou-se com profundidade razoável a teoria das leis de conservação em
domínios unidimensionais. Nesse ponto foram tratados conceitos como as condições de
Rankine-Hugoniot e estudaram-se as soluções fracas que, fisicamente, podem ser
interpretadas como ondas de choque ou rarefação. Para que se distinguissem soluções fracas
fisicamente corretas, incorporou-se ao presente trabalho o princípio de entropia, nas diversas
formas, desde as mais teóricas até aquelas com um perfil de aplicação prática, onde foram
incluídos os importantes princípios de entropia de Lax, Oleinik e Kruzkhov.
O esquema composto LWLF-4, recentemente introduzido por Liska e Wendroff,
combinou de forma apropriada o método clássico de Lax-Friedrichs com o esquema de Lax-
Wendroff, de modo que o primeiro funcionou como um filtro para o método de Lax-
Wendroff.
Os resultados com o método composto foram satisfatórios, uma vez que as soluções
obtidas não apresentaram as oscilações espúrias típicas do esquema Lax-Wendroff e nem tão
pouco excessivamente difusiva como o esquema de Lax-Friedrichs, obtendo assim o valor da
saturação no ponto de choque bem próximo ao valor encontrado analiticamente, sendo assim
um método viável para a captura do mesmo.
No problema de Burgers, as soluções numéricas foram se aproximando cada vez mais
da solução analítica quando a grade era refinada, obtendo-se resultados totalmente
satisfatórios em todas as condições iniciais propostas no trabalho.
132
O esquema composto mostrou também robustez quando comparados ao método TVD
com limitadores de fluxo. Os limitadores Superbee e Van Leer tiveram um valor de saturação
no ponto de choque um pouco superior ao valor analítico proposto por Welge (1952), sendo
ainda assim métodos viáveis para os problemas propostos. Entretanto o limitador de fluxo
Minmod assim como o método Composto obtiveram resultados altamente satisfatórios na
captura do choque para o problema de Buckley-Leverett.
Resoluções numéricas extensivas para as equações de Burgers e Buckley-Leverett
unidimensional indicam que o esquema composto LWLF-4 pode ser usado de forma bastante
satisfatória na simulação do processo de recuperação secundária de petróleo, o que recomenda
a continuidade do presente estudo para problemas bidimensionais, por exemplo.
6.2 Perspectivas
Como perspectiva futura, propomos o estudo numérico do escoamento imiscível
trifásico, por exemplo, água, óleo e gás, em meios porosos. Tal escoamento assim como o
bifásico, surge em muitos processos de importância prática em meios porosos, como na
produção primária de petróleo em condições próximas ao ponto de bolha do óleo, na injeção
de água em reservatórios contendo gás e óleo, na injeção de vapores, na injeção imiscível de
CO
2
, em alguns reservatórios de gases condensados (onde uma pequena variação na pressão
pode provocar o surgimento da fase óleo), na drenagem gravitacional de camadas de gás
contendo água e óleo, em processos de injeções alternadas de água e gás em reservatórios de
óleo pesado e na contaminação de lençóis freáticos através de misturas de componentes
químicos que em tais condições podem exibir duas fases, uma delas um gás e a outra um
líquido poluente.
Para ilustrar, descrevem-se a seguir as equações do modelo clássico para o escoamento
trifásico (água, óleo, gás) em meios porosos, na presença das seguintes hipóteses: (1)
problema unidimensional, (2) escoamento miscível, (3) fluidos incompressíveis, (4) efeito de
pressão capilar desconsiderado e (5) efeito de gravidade desconsiderado. De acordo com
Chavent e Jaffré (1986), se
[0, ]
x L
denota a posição e
0
t
o tempo, tais equações
governantes podem ser escritas na forma:
,0
Tx
u
(6.1)
133
,
0
0
g
w
xT
g
w
t
f
f
u
S
S
(6.2)
para todo
0
x L
,
0
t
.
Valores iniciais para as saturações da água e do gás devem ser conhecidas no instante
inicial. Juntamente com essas condições iniciais, devem-se ter condições de contorno
especificando diferentes valores de pressão na fronteira do domínio linear.
As equações descritas em (6.1) e (6.2) podem ser adimensionalizadas. Assim,
introduzindo-se as variáveis adimensionais
D
x
x
L
e
0
1
( )
t
D T
t d
L
v , o sistema (6.2) torna-
se
.
0
0
g
w
x
g
w
t
f
f
S
S
DD
(6.3)
Para simplificar a notação, o sistema mostrado na equação (6.3) será denotado por
,
0
0
g
w
x
g
w
t
f
f
S
S
(6.4)
ou simplesmente por
,0
~
fS
xt
(6.5)
onde
g
w
S
S
S
,
w
g
f
f
f
e
0
0
0
.
Uma vez que
0 1
s
, para todo
, ,
w g o
, e
1
s
, então o sistema (6.5)
encontra-se definido no triângulo de saturações:
.1,10,10;,
gwgwgws
SSSSSST
(6.6)
Experimentalmente, conhece-se que no escoamento em meios porosos uma fase
se
mantém imóvel se o seu valor de saturação é menor do que um certo valor crítico (ou
134
irredutível), denotado aqui por
c
S , para todo
, ,
w g o
. A partir do seu valor de saturação
crítica a fase pode escoar durante o fluxo multifásico em um material poroso. Para tratar dessa
questão de ordem prática, pode-se considerar que toda fase fluida encontra-se em quantidade
suficiente para escoar. Com essa hipótese, a equação (6.5) é escrita em termos das chamadas
saturações reduzidas das fases água e gás, dadas, respectivamente, por
cwwrw
SSS e
cggrg
SSS
, as quais, nessas condições, são quantidades não negativas, isto é, 10
rw
S e
10
rg
S
. Desse modo pode-se mostrar que
cgcwrgrw
SSSS
1
. Assim, tem-se o
chamado triângulo de saturações reduzidas, no qual efetivamente a equação (6.5) está definida
durante um deslocamento multifásico de interesse prático. Tal triângulo é dado por
.1,10,10;,
cgcwrgrwrgrwrgrws
SSSSSSSST
r
(6.7)
De acordo com a Figura 6.1, seja
rgrwgrgrwwrgrw
SSSSSS ,,,,
a aplicação
ssr
TT
:
(6.8)
definida por
cwrwrgrwww
SSSSS
,
e
cgrgrgrwgg
SSSSS
,
Figura 6.1 - Mudança de variáveis.
Nota-se que a aplicação
assim definida é uma bijeção de
sr
T
em
s
T
, logo é
inversível. Na realidade, essa aplicação é uma simples translação, sendo, portanto, uma
mudança de variáveis que preserva o sentido de qualquer vetor
r
srgrw
TSS ,
.
135
Com essa mudança de variáveis as permeabilidades relativas das fases e,
conseqüentemente, o vetor
f
do fluxo difusivo das fases água e gás podem ser escritos em
termos das saturações reduzidas. Assim, a equação (6.5) toma a forma
.0
fS
xrr
(6.9)
Na Eq. (6.9),
r
s
é o vetor das saturações reduzidas e
SfJSf
r
~
1
, onde
J
é a matriz
Jacobiana da aplicação
.
Nas últimas duas décadas, o estudo do sistema (6.9), que modela o escoamento
trifásico em meios porosos, tem despertado o interesse de muitos pesquisadores. Em
particular, a análise desse sistema com relação à perda de hiperbolicidade é uma área que tem
concentrado um grande esforço de pesquisa. Bell et al. (1986) observaram que os modelos de
permeabilidade de Stones (1970) proporcionam o surgimento de regiões elípticas dentro do
triângulo de saturações. Shearer (1988) e Holden (1990), entre outros autores, mostraram que
o aparecimento de tais regiões elípticas é na realidade uma característica comum decorrente
de muitos modelos clássicos para a permeabilidade relativas das fases.
Esses fatos têm (de certa forma) inibido a resolução da equação (6.9) através de
modernos métodos numéricos desenvolvidos especificamente para a resolução de sistemas
hiperbólicos não lineares, retardando, dessa maneira, o desenvolvimento de novos
simuladores para o escoamento trifásico em meios porosos.
Observando que a ocorrência de regiões elípticas contradiz o (esperado)
comportamento físico da solução de um sistema que modela o escoamento trifásico,
recentemente Juanes e Patzek (2004) desenvolveram e analisaram modelos para as
permeabilidades relativas que tornam o sistema (6.9) estritamente hiperbólico. Como um
exemplo, esses autores consideram as seguintes mobilidades relativas, escritas em termos das
saturações reduzidas,
,
1
2
rw
w
w
S
(6.10)
,1
1
2
rwgrgg
g
g
SS
(6.11)
136
,111
1
rgrwrgrw
o
o
SSSS
(6.12)
onde o coeficiente do lado direito da equação (6.11) cumpre a condição
g
g
o w
e
2
w o
. As Figuras 6.2 e 6.3 mostram o modelo para as funções de fluxo do fluido água e
gás.
A perspectiva de pesquisas futuras tem como objetivo utilizar os esquemas compostos
referidos neste trabalho (os quais foram desenvolvidos por Wendroff e colaboradores) na
resolução de equações e sistemas de equações hiperbólicas não-lineares representando leis de
conservação em uma, duas e três variáveis espaciais. Particular ênfase será dada ao estudo
numérico do sistema
2 2
resumido nas equações (6.9)-(6.12). Com isso, objetiva-se de uma
forma mais ampla a construção de simuladores para o escoamento imiscível trifásico em
reservatórios de petróleo.
Figura 6.2 – Função de fluxo fracionária da água (Modelo de Juanes e Patzek)
137
Figura 6.3 – Função de fluxo fracionária do gás (Modelo de Juanes e Patzek)
138
REFERÊNCIAS
AHMED, R. Numerical schemes applied to the burgers and Buckley-Leverett equations.
Dissertação (Mestrado) - University of Reading, Berkshire, 2004.
AZIZ, K.; SETTARI, A. Petroleum reservoirs simulation. London: Applied Science, 1979.
BEAR, J.; BACHMAT, Y. Introduction to modeling of transport phenomena in porous
media. Dordrecht : Kluwer, 1991.
BLUNT, M.; RUIN, B. Implicit flux limiting schemes for petroleum reservoir simulation.
Journal of Computation Physics, 102, 194, 1992.
BROOKS S.E.; COREY A.T. Hydraulic properties of porous media. Hydrol Pap., v. 3. Fort
Collins: Colorado State Univ., 1964.
BUCKLEY, S. E.; LEVERETT, M. C. Mechanism of fluid displacement in sands. Trans.
AIME, 146, 107, 1942.
CAPELEIRA PINTO, A. C. Esquemas de alta resolução para controle de dispersão
numérica em simulaçao de reservatórios. Dissertação (Mestrado) – Universidade Estadual de
Campinas, Campinas, 1991.
CHEN, Z.; HUAN, G.; M. A, Y. Computational methods for multiphase flow in porous
media. Texas: Southerm Methodist University, 2005.
CHORIN, A.; MARSDEN, J. E. A mathematical introduction to fluid mechanics. New
York: Springer, 2000.
COURANT, R.; ISSACSON, E.; REES, M. On the solution of nonlinear hyperbolic
differential equations by finite differences. Comm.Pure Appl. Math., 5, p. 243, 1952.
DARCY, H. Les fountaines publiques de la ville de Dijon. Paris: Dalmont, 1856.
GUEVARA-JÓRDAN, J.M.; CELLA, M. B. Métodos de alta resolución para la ecuación de
Buckley-Leverett. Universidad Central de Venezuela, v.7, n.2, p.105-131, 1999.
HARTEN, A.High resolution schemes for hyperbolic conservation laws. J. Comput. Phys., v.
49, p. 357-393, 1983.
HARTEN, A.; LAX, P.D.; VAN LEER, B. On upstream differencing and Godunov type
schemes for hyperbolic conservation laws. SIAM Review, v. 25, n.1, 1983.
KRUZKHOV, S. N. Generalized solutions of the Cauchy problem in the large for nonlinear
equations of first order. Dokl. Akad. Nauk. SSSR, v. 187, p. 29-32, 1969
LAKE, L.W. Enhanced oil recovery. New Jersey: Prentice-Hall, 1989.
139
LAX, P. Hyperbolic systems of conservation law II. Comm. Pure Applied Math., v. 10, p.
537-566, 1957.
LEVEQUE, R.J. Finite volume methods for hyperbolic problems. Cambridge: Cambridge
University Press, 2002. Cambridge Texts in Applied Mathematics.
LISKA, R.; WENDROFF, B. Composite schemes for conservation laws. Society for
Industrial and Applied Mathematics, v. 35, n. 6, 1998.
OLEINIK, O. A. Uniqueness and stability of the generalized solution of the Cauchy problem
for a quasilinear equations. Amer. Math. Soc. Transl., v. 33, p. 285-290, 1963.
OSHER, S. Riemann solvers, the entropy condition, and difference approximations. SIAM J.
Numer. Anal., v. 21, p. 217-235, 1984.
PEACEMAN, D.W. Fundamentals of numerical reservoir simulation. New York: Elsevier,
1977.
POSTON, S.W.; YSRAEL, S.C.; HASSAIN, A.K.; MONTGOMERY, E.F.; ROMEY,
H.J.JR. The effect of temperature on irreducible of uncomolidated sands. Trans.
SPE AIME, v. 249, p.178-180, 1970.
SERRE, D. Systèmes de lois de conservation I: hyperbolicité, entropies, ondes de choc.
Paris: Diderot, 1996.
SETHLAN, J. A.; CHORIN, A. J.; CONCUS, P. Numerical solution of the Buckley-Leverett
equations. Society of Petroleum Engeneers, Paper SPE-12254, 1983.
SMOLLER, J. Shock waves and reaction-diffusion equations. 2nd. ed. New York: Springer,
1994.
STONE, H. L. Probability model for estimating three-phase relative permeability. Trans. SPE
of AIME, v. 249, p. 214-218, 1970.
SWEBY, P. K. High resolution schemes using flux limiters for hyperbolic conservation laws.
SIAM Journal of Numerical Analysis, v. 21, n.5, 1984.
SWEBY, P. K. Theory of differential equations. Berkshire: University of Reading,
2002. Lecture notes.
TADMOR, E. Numerical viscosity and the entropy condition for conservative difference
schemes. Mathematics of Computation, v. 43, n. 369, 1984.
THOMAS, J. W. Numerical partial differential equations: conservation laws and
elliptic equations. New York: Springer, 1999. Texts in Applied Mathematics, 33.
TRANGENSTEIN, J. A. Numerical solution of hyperbolic partial differential equations. New
York: Cambridge University Press, 2007.
140
WELGE, H. J. A simplified method for computing oil recovery by gas or water drive. Trans.
AIME, n. 195, p. 91-99, 1952.
WENDLAND E.; FLENSBERG D. Numerical solution of two-phase flow for the advection
dominated and non-linear case. Advances in Water Resources, v. 28, p. 643-660, 2005.
Livros Grátis
( http://www.livrosgratis.com.br )
Milhares de Livros para Download:
Baixar livros de Administração
Baixar livros de Agronomia
Baixar livros de Arquitetura
Baixar livros de Artes
Baixar livros de Astronomia
Baixar livros de Biologia Geral
Baixar livros de Ciência da Computação
Baixar livros de Ciência da Informação
Baixar livros de Ciência Política
Baixar livros de Ciências da Saúde
Baixar livros de Comunicação
Baixar livros do Conselho Nacional de Educação - CNE
Baixar livros de Defesa civil
Baixar livros de Direito
Baixar livros de Direitos humanos
Baixar livros de Economia
Baixar livros de Economia Doméstica
Baixar livros de Educação
Baixar livros de Educação - Trânsito
Baixar livros de Educação Física
Baixar livros de Engenharia Aeroespacial
Baixar livros de Farmácia
Baixar livros de Filosofia
Baixar livros de Física
Baixar livros de Geociências
Baixar livros de Geografia
Baixar livros de História
Baixar livros de Línguas
Baixar livros de Literatura
Baixar livros de Literatura de Cordel
Baixar livros de Literatura Infantil
Baixar livros de Matemática
Baixar livros de Medicina
Baixar livros de Medicina Veterinária
Baixar livros de Meio Ambiente
Baixar livros de Meteorologia
Baixar Monografias e TCC
Baixar livros Multidisciplinar
Baixar livros de Música
Baixar livros de Psicologia
Baixar livros de Química
Baixar livros de Saúde Coletiva
Baixar livros de Serviço Social
Baixar livros de Sociologia
Baixar livros de Teologia
Baixar livros de Trabalho
Baixar livros de Turismo
