Transcript Advecção Difusão TEM

1ª Aula Advecção - Difusão
Objectivos deste capítulo e
Método dos volumes finitos.
Objectivos
• Este capítulo tem como objectivos apresentar métodos
de resolução da equação de Advecção – Difusão e fazer
uma aplicação num sistema unidimensional.
• Este capítulo dá continuidade ao problema de difusão
resolvido em Mecânica dos Fluidos Ambiental. Usa o
mesmo código desenvolvido em VBA, adicionando o
transporte pela velocidade e juntando alguma
complexidade às condições de fronteira num problema
com superfície livre.
• O trabalho desenvolvido dá suporte teórico para
Modelação Ambiental.
Programa deste capítulo
• Revisão do método do volume finito para quantificação do princípio
de conservação “a taxa de acumulação é igual ao que entra, menos
o que sai, mais o que se produz menos o que se consome”.
• Particularidade da advecção por necessitar dos valores sobre as
faces do volume finito. Método upwind e método do valor médio
(diferenças centrais). Outros métodos de resolução.
• A questão do tempo: métodos explícitos, implícitos e CrankNicholson (semi-implícitos).
• A questão da difusão numérica e da estabilidade. Relação entre as
propriedades dos métodos numéricos e os princípios físicos. Nº de
Courant e nº de Difusão.
• Dedução das equações algébricas a partir das equações diferenciais
e das séries de Taylor. Erro de truncatura e precisão do método.
Processos
• Taxa de acumulação:


TaxadeAcumulação   cdvol 
t  vol

• Fluxos:
– Advectivo:
– Difusivo:


cu .n dA


A
A

 
 c .n dA
– Porquê o sinal “-” antes dos integrais?
Localização das variáveis no volume de
controlo: Fluxos advectivo e difusivo
através das faces
cxt txVolxt tx  cxt xVolxt x
cxt tVolxt t  cxt Volxt

*
*
  cu.n dA  ucdA x  x / 2  cQ x x / 2
A
 
 
 *

c*x  c*x x




c
.
n
dA


dA
A
 x x / 2
x


*
cxt txVolxt tx  cxt xVolxt x
ucdA*xx / 2  cQ*xx / 2
 *

c
c
x  x / 2
dA
x


*
x  x
*
x
*
Aplicando o princípio de conservação
A taxa de acumulação é igual ao que entra menos o que
sai, mais o que se produz menos o que se destrói. Se não
existir produção nem destruição, fica só:
c
t  t
x

Volxt  t  c xt Volxt

t

c c

 Qcx  x / 2    Ax  x / 2  x x  x   
 x  


 c x  x  c x  
 Qcx  x / 2    Ax  x / 2 
 
 x  

Hipótese Upwind para a concentração
na face
c 
 c  se : u 
0
c *x x / 2  c *x se : u *x x / 2  0
*
*
*
c x  x / 2  c x se : u x  x / 2  0
*
*
*
c x  x / 2  c x  x se : u x  x / 2  0
*
x  x / 2
*
x  x
*
x  x / 2
• No caso de velocidade positiva (escoamento
para a direita):
c
t  t
x
 cx 
 c  c xx 
c
Vol xt  t  c xt Vol xt / t  Qxx / 2 C xx  Qxx / 2 C x     Axx / 2  x
    Ax x / 2  xx

x
x





Teste em problema unidimensional com
volume constante e caudal uniforme
Ci
Ci+1
Ci-1
c
t  t
x
 cx 
 c  c xx 
c
Vol xt  t  c xt Vol xt / t  Qxx / 2 C xx  Qxx / 2 C x     Axx / 2  x
    Ax x / 2  xx

x
x





Se o volume for constante e o caudal e a difusividade forem uniformes fica, em
upwind explícito:
c
t  t
i



 cit
cit1  cit
cit1  2cit  cit1
u

t
x
x 2
Rearranging the equation
t  t
i
c
 ut t  t  ut 2t  t  t  t

 2 ci 1  1 

c   2 ci 1
2  i
x x 
 x x 

 x 
• Using 3 vectors to store the coefficients that relate the
new (at t+dt) with values at time t in the point and in
the neighbouring points one can write:
cit t  diCit1  1  ei Cit  fiC txx
ei  di  fi 
Explicit, Upwind, Cr = 1, Dif=0
t  t
i
c
 ut t  t  ut 2t  t  t  t

 2 ci 1  1 

c   2 ci 1
2  i
x x 
 x x 

 x 
In this table each line representa a time instant (the first is the initial solution) and each
column displays the solution in one point. The right column is the addition of the solution
along one line (i.e is the total mass).
Time step i-3
0
1
2
3
4
i-2
0
0
0
0
0
i-1
0
0.00
0.00
0.00
0.00
Grid point number
i
i+1
0
1
0.00
0.00
0.00
0.00
0.00
0.00
0.00
0.00
i+2
0
1.00
0.00
0.00
0.00
i+3
0
0.00
1.00
0.00
0.00
Cr=(Espaço percorrido num intervalo de tempo)/(passo espacial)
Cr=1, implica uma célula por passo => a solução é exacta
Total
amount
0
0
0
0
0
1
1
1
0
0
Explícito, Upwind, Cr= 0.5, Dif=0
t  t
i
c
Time step i-3
0
1
2
3
4
5
6
7
8
9
10
 ut t  t  ut 2t  t  t  t

 2 ci 1  1 

c   2 ci 1
2  i
x x 
 x x 

 x 
i-2
0
0
0
0
0
0
0
0
0
0
0
i-1
0
0.00
0.00
0.00
0.00
0.00
0.00
0.00
0.00
0.00
0.00
Grid point number
i
i+1
0
1
0.00
0.50
0.00
0.25
0.00
0.13
0.00
0.06
0.00
0.03
0.00
0.02
0.00
0.01
0.00
0.00
0.00
0.00
0.00
0.00
i+2
0
0.50
0.50
0.38
0.25
0.16
0.09
0.05
0.03
0.02
0.01
i+3
0
0.00
0.25
0.38
0.38
0.31
0.23
0.16
0.11
0.07
0.04
0
0
0
0
0
0
0
0
0
0
0
Total
amount
1.00
1.00
1.00
0.88
0.69
0.50
0.34
0.23
0.14
0.09
0.05
A mancha espalha-se e por isso temos difusão numérica. Porque será?
É porque violámos a definição de concentração. Como se resolve? Com uma malha mais fina
Como aconteceu?
t0
0
0
1
0
0
0
0
0
0
t0+Δt
0
0
0.5
0.5
0
0
0
0
0
t0+2Δt
0
0
0.25
0.5
0.25
0
0
0
0
t0+3Δt
0
0
0.125
0.375
0.375
0.125
O modelo é estável: os erros que aparecem diminuem no tempo.
O modelo tem difusão numérica: a concentração vai baixando
apesar de a difusão física ser nula.
Explícito, Upwind, Cr=2
t  t
i
c
 ut t  t  ut 2t  t  t  t

 2 ci 1  1 

c   2 ci 1
2  i
x x 
 x x 

 x 
t0
0
0
1
0
0
0
0
0
0
t0+Δt
0
0
-1
2
0
0
0
0
0
t0+2Δt
0
0
+1
-4
4
0
0
0
0
t0+3Δt
0
0
-1
10
-16
8
Temos um modelo instável: os erros aparecem e crescem. Porquê? Num
modelo explícito Cr≤1. Os coeficientes têm que ser positivos.
As instabilidades são consequências da
violação de princípios físicos
• Quando as propriedades aumentam num
instante, nos instantes seguintes também só
podem aumentar.
• Quando Cr>1 o coeficiente de Ci fica negativo.
• Neste caso, durante um intervalo de tempo o
volume que sai de uma célula é maior do que
o que lá estava no início (Usando volumes finitos
é fácil ver que isso é a causa do problema).
• (ver Patankar, Fluid Flow)
Condição de estabilidade
t  t
i
c
 ut t  t  ut 2t  t  t  t

 2 ci 1  1 

c   2 ci 1
2  i
x x 
 x x 

 x 
Condição de estabilidade:
t 
 ut
2 2 0
1 
x
x 

2ª Aula Advecção - Difusão
Diferenças Centrais. Método
implícito. Método QUICK
Outra opção: Valores médios nas faces
=>Diferenças Centrais
c 
 c x  c x  x 


2


c 
 c x  c x  x 


2


*
x  x / 2
*
x  x / 2
*
*
c*xx / 2  c*xx / 2   cx  cxx  cx  cxx *   cxx  cxx *
x


2x
2x




2x


Diferenças Centrais Explícitas
c
t  t
x
cit  t
 cx  cx x 
 cx x  cx 
 c Ax  AuCx x  Cx  x     Ax x / 2 
    Ax x / 2 

 x 
 x 
t  t  ut t  t
 ut t  t 

 2 C i1  1  2 2 C i   
 2 C i1
x 
 x x 

 x x 
t
x

A equação continua a poder ser escrita na forma:
cit t  diCit1  1  ei Cit  fiC txx
ei  di  fi 
1D explicit central differences Courant=1
t  t
i
c
 ut t  t  2t  t  ut t  t

 2 ci 1  1 
c  
 2 ci 1
2  i
x 
 2x x 

 2x x 
Time step i-3
0
1
2
3
4
5
6
7
8
9
10
11
i-2
0
0
0
0
0
0
0
0
0
0
0
0
Grid point number
i-1
i
i+1
0
0
1
0.00
-0.50
1.00
0.25
-1.00
0.50
0.75
-1.13
-0.50
1.31
-0.50
-1.63
1.56
0.97
-2.13
1.08
2.81
-1.16
-0.33
3.93
1.66
-2.29
2.94
5.59
-3.76
-1.00
8.52
-3.26
-7.14
7.52
0.31
-12.54
0.38
i+2
0
0.50
1.00
1.13
0.50
-0.97
-2.81
-3.93
-2.94
1.00
7.14
12.54
Total
amount
i+3
0
0.00
0.25
0.75
1.31
1.56
1.08
-0.33
-2.29
-3.76
-3.26
0.31
0
0
0
0
0
0
0
0
0
0
0
0
Modelo Instável. Porquê? Há um dos coeficientes que é sempre negativo.
Propriedade transportiva violada. Como se resolve?
1
1
1
1
1
1
1
1
1
1
1
1
Porque é instável?
• Por advecção (ou por difusão) quando as
propriedades aumentam num ponto, nos
pontos vizinhos só podem aumentar também.
• Isso implica que os coeficientes que
multiplicam as concentrações nos pontos
vizinhos têm que ser positivos.
• Só adicionando difusão é que isso pode
acontecer….
Condição de estabilidade para
diferenças centrais explícitas
Porque é que adicionando difusão o método fica estável?
Porque é que excesso de difusão torna o modelo instável?
Interpretação das diferenças centrais
• Porque é que as diferenças centrais são instáveis
sem difusão?
– Resp: Violam a propriedade transportiva. Um ponto
fica a saber o que está abaixo através da advecção, o
que é fisicamente impossível.
• Porque é que a difusão pode estabilizar as
diferenças centrais?
– Resp: Porque a difusão transporta a informação para
montante. No caso de a difusão ser importante a
advecção transporta efectivamente para jusante
coisas que foram transportadas para montante pela
difusão.
Continuação
• Poderão as diferenças centrais explícitas ser usadas quando a advecção é
dominante?
– Resp: Não. Nesse caso difusão transporta para montante muito menos do que
a advecção transporta para jusante (Reynolds da malha)
• Se a difusão for dominante é preferível usar diferenças centrais ou
upwind?
– Se a difusão for dominante as diferenças centrais são vantajosas porque têm
precisão de 2ª ordem e por isso introduzem menos difusão numéricas
• E se o algoritmo fosse implícito? Seria o algoritmo mais estável?
– Resp: Sim. Nesse caso a solução seria função dos valores das variáveis no
passo de tempo seguinte. Se a advecção tende a criar concentrações
negativas, a difusão aumenta automaticamente para porque o gradiente de
concentração aumenta.
• E se o método fosse upwind?
– Resp: nesse caso as concentração não pode ficar negativa. Em upwind a
concentração fica negativa se retirarmos de uma célula mais do que lá existe
para sair. Mas como em implícito o que sai é função da nova concentração, se
ela ficasse negativa isso significaria que sairia uma quantidade negativa e por
isso a concentração cresceria…..
Outros métodos para a advecção
• Upwind: Passa numa face o que está a montante.
• Diferenças centrais: Passa numa face a média do que está
dos dois lados.
• E se ajustássemos um polinómio de 2ª ordem a 3
pontos? Obteríamos o método QUICK: (Quadratic
Upstream Interpolation for Convective Kinematics):

Q  0  C
Qi  0  Ci  1  6 8 Ci 1  3 8 Ci  18 Ci 2
i
2
i  12
 6 8 Ci  3 8 Ci 1  18 Ci 1


• Tem precisão de terceira ordem. Tem mais problemas de
estabilidade (em situações particulares, nomeadamente
junto às fronteiras.
• Afinal qual é o melhor método?
3ª Aula Advecção - Difusão
Cálculo Implícito
Método implícito
c
t  t
x
c

 c Ax  AuC x  x  C x 
t
x
t  t

C  Cx 
 cxt
 u x x
t
x
t  t
x
t  t
c c

   Ax  x / 2  x x  x 
 x 
 2cx  cx  x 
c
   x  x

x 2


t 
 ut t  t  t  ut

 2 c xx  1 
 2 2 cxt  t
x
x 
 x x 

 diCi1  1  ei C
t t
i
c 
c
   Ax  x / 2  x  x x 
 x 
t  t
t  t
 t 
  2 cxt tx  c tx
 x 
 diCit1 t  1 ei Cit t  fiCit1 t  Cit
t t
t  t
t t
 fiCi1  TIi
ei  di  fi 
In this case the Independent
Term is just the concentration
at time t.
Porque serão os métodos implícitos
incondicionalmente estáveis?
• UPWIND
– No método explícito o que sai de uma célula é o que lá
está em “t”. No método implícito o que sai é o que lá vai
estar em “t+dt”.
– No método explícito, quando se retira de uma célula mais
do que lá está para sair, a concentração fica negativa.
– No método implícito não pode ficar porque o que sai é
função do que lá vai estar e por isso, se a concentração
pudesse ficar negativa, sairia uma quantidade negativa e
por isso a concentração iria aumentar e não diminuir. Isto
mostra que é impossível ficar com concentrações negativas
em upwind.
– E em diferenças centrais?
Porque são as diferenças centrais
implícitas mais estáveis que as explícitas?
• No caso das diferenças centrais, o que entra numa célula é
o que está a montante e o que sai é calculado em função
do que está a jusante ( em explícito viola a propriedade
transportiva da advecção).
• Em explícito, sem difusão a solução é instável (viola a
propriedade transportiva). Em implícito, o que sai de uma
célula é o que vai estar a jusante e o que entra é o que vai
estar a montante. Se a concentração a montante de uma
célula for nula, nessa zona ela vai ter que ficar negativa. No
entanto, o valor negativo a montante vai entrar na célula de
jusante e vai fazê-lo baixar, o que implica que vai ser
removido menos material de montante e por isso que a
concentração vai ser menos negativa.
Diferença entre métodos explícitos e
implícitos
c
Método implícito
Têm erros da mesma ordem de
grandeza. Se um é por excesso o
outro é por defeito. O método ideal
é a média dos dois.
Método Explícito
t1
t1+Δt
t
Método Semi-implícito (Crank –
Nicholson)
Método explícito:
cit t  diCit1  1  ei Cit  fiC txx
Método implícito:
 diCit1 t  1 ei Cit t  fiCit1 t  Cit
Método Semi-implícito (Crank – Nicholson):
1
1
1
1
 1 
 1 
 di C it1 t  1  ei Cit  t  f i C it1 t  di C it1  1  ei Cit  f i C it1
2
2
2
2
 2 
 2 
Requer o dobro das contas, mas deve ser mais preciso.
Summary
• The rationale used to obtain the algebraic
equations was based on the direct application of
the conservation principle to a finite volume,
assuming that it is small enough to guarantee
that the average value is representative of each
point of the volume and that the properties over
each face can be assumed as uniforms and also
that the time step is small enough to assume that
its value does not change over the time step.
• Could we get the same result starting from the
differential equations?
3ª Aula Advecção - Difusão
Séries de Taylor para obtenção das
equações algébricas.
Formas da equação
ck u j ck



t
x j
x j
 c 

  ( Fk  Pk )
 x 
j 

Escrevendo na forma da divergência dos fluxos:
ck


t
x j

 u j ck  c

x j


  ( Fk  Pk )


Onde o 1º termos do 2º membro é o simétrico da divergência dos fluxos, i.e.
o que entra menos o que sai.
Ou na forma convencional:
dck ck
ck


uj

dt
t
x j x j
 c 

  ( Fk  Pk )
 x 
j 

Séries de Taylor
t
cit  t
t
t
t
2
2
3
3
t n   n c 
 c  t   c  t   c 
t
 2  
 3   ....
 n 
 ci  t   
n!  t i
 t i 2!  t i 3!  t i
• Estão na base do método das diferenças
finitas, que são da mesma família dos Volumes
Finitos.
• Os Elementos Finitos/Elementos de fronteira
são a segunda principal família de métodos
numéricos.
O que representa a série de Taylor?
cit  t
t
t
t
2
3
3
n
n






c

t

c

t

c

t

c


t
 2  
 3   ....
 n 
 ci  t   
n!  t i
 t i 2!  t i 3!  t i
t
2
c
Outras
derivadas
Δc
Δc 1ª Derivada: Δc/ Δt
Δt
t1
t1+Δt
t
Como usar para calcular as derivadas?
t
t
t
 c  t   c  t   c 
 2  
 3   ....
 c  t   
2!  t i 3!  t i
n!
 t i
t
t  t
i
t
i
c
2
2
3
3
n
 c 
cit  t  cit  t    (t 2 )
 t i
t
t  t
 c  ci  ci
 (t )
  
t
 t i
t
t
 c
 n 
 t i
n
t
Método Explícito: A derivada é calculada à esquerda “em t” e tem precisão de 1ª ordem,
ou seja, as derivadas que foram ignoradas estão multiplicadas por (t )
A derivada ser calculada à esquerda significa “à esquerda do intervalo de tempo”, i.e.
em “t” e por isso o método é explícito. Todas as derivadas (i.e. todos os termos da
equação) são calculados em “t”.
O erro ser proporcional a (t ) significa que “o erro do cálculo aumenta quando o
passo de tempo aumenta.
Mas poderia ter feito calculado a
derivada à direita do intervalo de tempo
t  t
t  t
i
c c
t
i
 c 
 t  
 t i
t

2!
2
t  t
 c
 2 
 t i
2
t

3!
3
t  t
 c
 3 
 t i
3
t
 ....
n!
t  t
cit  cit  t
t  t
 c 
 
 t i
 c 
 t  
 t i
 (t 2 )
cit  t  cit

 (t )
t
Método Implícito: A derivada é calculada à direita “em t+dt” e tem precisão
de 1ª ordem, ou seja, todas as derivadas que foram ignoradas estão
multiplicadas por (t )
Isto significa que o erro do cálculo aumenta quando o passo de tempo
aumenta. Os métodos implícitos e explícitos têm a mesma precisão.
n
t  t
 c
 n 
 t i
n
Para calcular a derivada no centro do intervalo teria
que calcular os valores nos extremos a partir daquele
t  t / 2
t  t
i
c
t  t / 2
i
c
 c 
 t / 2 
 t i
t  t / 2
i
c c
 c 
 t / 2 
 t i
2

t / 2   2 c 



2!
t  t / 2
 t 2 

i
2!
t  t / 2
t
i
2

t / 2   2 c 



t  t / 2
 t 2 

i
3

t / 2   3c 



3!
 t 3 

i
3

t / 2   3c 



3!
t  t / 2
t  t / 2
 t 3 

i
n

t / 2   n c 


 ....
n!
 t n 

i
n

t / 2   n c 


 ....
n!
t  t / 2
t  t / 2
 t n 

i
Subtraindo uma da outra:
t  t / 2
cit  t
 c 
 cit   t  
 t i
t  t / 2
 c 
 
 t i

  t / 2

cit  t  cit
2

  t / 2
t
3


Neste método a derivada é calculada no centro do intervalo de tempo e
tem precisão de 2ª ordem. Dá a solução exacta até uma evolução
parabólica. As derivadas ignoradas estão multiplicadas por t / 22
O que representa a série de Taylor?
cit  t
c
t
t
t
2
3
3
n
n






c

t

c

t

c

t

c


t
 2  
 3   ....
 n 
 ci  t   
n!  t i
 t i 2!  t i 3!  t i
t
Método Explícito
2
Método Implícito
Outras
derivadas
Δc
1ª Derivada: Δc/ Δt
Δt
Método Diferenças Centrais
t1
t1+Δt
t
Derivadas espaciais
t
i  x
c
t
t
x
 c  x   c  x   c 
 2  
 3   ....
 c  x  
2!  x i
3!  x i
n!
 x i
t
2
2
3
3
t
i
n
 c 
t
t
ci  x  ci  x   (x 2 )
 t i
Derivada à direita, Método downwind, se
t
velocidade positiva
cit x  cit
 c 
 (x)
  
x
 x i
t
Neste método a derivada espacial num ponto é calculada a partir da
informação no ponto e da informação à direita.
Veremos mais adiante que este cálculo cria problemas se esta derivada for
usada para calcular o termo advectivo quando a velocidade é positiva.
t
 c
 n 
 x i
n
Derivadas espaciais
t
i  x
c
t
t
x
 c  x   c  x   c 
 2  
 3   ....
 c  x  
2!  x i
3!  x i
n!
 x i
t
2
2
3
3
n
t
i
 c 
cit x  cit  x   (x 2 )
 t i
t
 c
 n 
 x i
n
t
t
t
 c  ci  ci  x
 (x)
  
x
 x i
t
Derivada à esquerda: “Método upwind” se
velocidade positiva e downwind se fosse negativa.
Neste método a derivada espacial num ponto é calculada a partir da
informação no ponto e da informação à esquerda. Este método respeita a
propriedade transportiva da velocidade se esta for positiva, mas não se for
negativa. Nesse caso a derivada deveria ser calculada “à direita”.
Subtraindo uma equação da outra
t
i  x
c
2
c
2
3
t
3
t
3
n
t
x
 c  x   c  x   c 




 c  x  

 ....
2 
3 


2!  x i
3!  x i
n!
 x i
2
2
3
t
i
 c 
3
 2x   (x )
 t i
t
t
i  x
c
c
t
i  x
 c  c  c
  
2x
 x i
t
t
i  x
t
i  x
 (x )
Diferenças Centrais
2
t
 c
 n 
 x i
n
t
i
t
t
i  x
t
x
 c  x   c  x   c 
 2  
 3   ....
 c  x  
2!  x i
3!  x i
n!
 x i
t
n
t
 c
 n 
 x i
n
2ª Derivada
*
*
c
x
 c  x   c  x   c 
 2  
 3   ....
 c  x  
2!  x i
3!  x i
n!
 x i
c
x
 c  x   c  x   c 
 2  
 3   ....
 c  x  
2!  x i
3!  x i
n!
 x i
*
*
i  x
2
3
3
n
*
i
*
*
i  x
2
2
*
2
3
*
i
Adicionando:
*
ci* x  ci* x
*
2


c
*
2
 2ci  x  2   (x 4 )
 t  i
  c  ci* x  2ci*  ci* x
2
 2  


(

x
)
2
x
 x i
2
3
*
*
 c
 n 
 x i
n
n
*
 c
 n 
 x i
n
4ª Aula Advecção - Difusão
Equações algébricas. Erro de
truncatura, condições iniciais e
condições de fronteira.
Sumário da aula anterior
•
•
•
•
Na última aula vimos como obter equações algébricas a partir das equações
diferenciais, usando séries de Taylor.
Vimos que poderíamos obter facilmente discretizações com precisão de primeira
ou de segunda ordem no tempo e/ou no espaço e vimos o que queria dizer o erro
de truncatura.
Combinando este conhecimento com o que obtivemos quando analisamos o
problema com o método dos volumes finitos concluímos que nem sempre o
menor erro de truncatura significa menor erro dos resultados.
Para se obterem bons resultados é necessário garantir o respeito pelos princípios
físicos, nomeadamente:
– Conceito de Concentração, que tem que ser mais ou menos uniforme no interior da célula,
– A transportividade da advecção,
– Que uma célula não é despejada numa iteração (Cr ≤ 1).
•
Os métodos implícitos respeitam os processos físicos de forma semelhante aos
explícitos e são mais estáveis. OS métodos semi-implícitos são mais estáveis e têm
maior precisão que os explícitos.
Equações Algébricas
• Obtêm-se substituindo as derivadas pelas
aproximações:
t t
cx  cxt
cxt  x  cxt x
cxt  x  2cxt  cxt x
2
2
 t   u
  x  



x
t
2x
x 2
 
 
• Explícito, diferenças centrais. Precisão de 2ª ordem no
espaço e 1ª no tempo.
t t
cx  cxt
cxt tx/ 2  cxt tx/ 2
cxt tx/ 2  2cxt  t / 2  cxt tx/ 2
2
2
2
 t   u
  x  



x
t
2x
x 2
 
 
• Semi-implícito (Crank-Nicholson) diferenças centrais
espaço. Precisão de 2ª ordem no tempo e no espaço.
O que se paga pela precisão de 2ª ordem no tempo?
Como se obtém o valor em (t+Δt/2) ?
Fazendo a média…..
t  t / 2
cit  t
 c 
 cit  t / 2  t / 2 
 t i
t  t / 2
 c 
cit  cit  t / 2  t / 2 
 t i
2

t / 2   2 c 



2!
 t 2 

i
2

t / 2   2 c 



2!
t  t / 2
3

t / 2   3c 



3!
 t 3 

i
3

t / 2   3c 



t  t / 2
 t 2 

i
3!
t  t / 2
t  t / 2
 t 3 

i
n

t / 2   n c 


 ....
n!
• Adicionando as equações!
t  t / 2
2cit  t / 2  cit  cit  t
t  t / 2
i
c
  2c 
 t / 2   2 
 t i
2
 t n 

i
n

t / 2   n c 


 ....
n!
 .....
cit  cit  t
2

 t / 2 
2
• Substituindo estes termos nas equações
obtém-se a equação a resolver.
t  t / 2
t  t / 2
 t n 

i
Explícito Upwind
cxt  t  cxt
cxt  cxt x
cxt  x  2cxt  cxt x
2
2
 t   u
  x  



x
t
x
x 2
 
 
• Precisão de 1ª ordem no tempo e no espaço para advecção.
Segunda ordem para difusão.
• Esta equação pode ser organizada na forma:
c
t  t
i
 ut t  t
 ut 2t  t  t  t

 2 ci 1  1 

c   2 ci 1
2  i
x x 
 x x 

 x 
cit t  di cit1  1  ei cit  fi cit1  (F  P)
Forma geral da Equação
kdi cit1t  1  kei cit t  kfi cit1t  1  k di cit1  1  1  k ei cit  1  k  fi cit1  (F  P)
K=1=> implícito. K=0 => Explicito, k=0.5=> Crank-Nicholson:
Explicito, upwind:
c
t  t
i
 ut t  t
 ut 2t  t  t  t

 2 ci 1  1 

c   2 ci 1
2  i
x x 
 x x 

 x 
Números de Courant e de Difusão
Cr 
ut
x
 t
N º Dif  2
x
Sobre a precisão do cálculo
• No cálculo implícito e no cálculo explícito as derivadas são
calculadas nos extremos do intervalo de tempo. Estes métodos
ignoram todas as derivadas a partir da primeira: têm precisão de
primeira ordem ou “até à primeira ordem”.
• Os termos da série de Taylor ignorados estão multiplicados por (t )
• Quando a derivada é calculada no centro do intervalo de tempo as
derivadas só são ignoradas a partir da segunda. São métodos com
precisão de 2ª ordem, ou “até à 2ª ordem”. Se a função for uma
recta ou uma parábola o cálculo da derivada é exacto.
2
• Os termos da série de Taylor ignorados estão multiplicados por t / 2
• Mas (t ) >1 então quanto maior é a ordem de precisão do
cálculo, maior é o coeficiente dos termos ignorados. Porque é
que a precisão do cálculo aumenta?
Porque aumenta a precisão com o
expoente de (t ) ?
Porque os termos ignorados são da forma:
t / 2
n!
n 1
t  t / 2
  nc 
 n 
 t  i
O cálculo da derivada faz aparecer em denominador o intervalo de
tempo elevado n e o coeficiente está elevado a (n-1) e por isso o
produto é proporcional a (c) /(t )ou seja à primeira derivada
multiplicada pelo inverso do factorial de n e por isso quanto maior é o
valor do expoente do intervalo de tempo, menor é o valor dos temos
desprezados.
Esta conclusão é consistente como facto de as derivadas perderem
importância à medida que a ordem aumenta.
Condições Iniciais e de Fronteira
• Iniciais podem ser importantes ou não
• Fronteira idem. Como se impõem?
Ci
Ci-1
Ci+1
 diCit1 t  1 ei Cit t  fiCit1 t  Cit
Condições de fronteira
• Difusão:
– Requer o cálculo dos fluxos nas células de fronteira e
por isso requer a concentração no exterior em ambas
as fronteiras. Se não for conhecida a melhor solução é
normalmente gradiente nulo.
• Advecção
– Quando o escoamento entra no domínio transporta as
propriedades do exterior. As propriedades têm que
ser conhecidas no exterior. Se não forem conhecidas,
a simulação só pode fazer sentido se as fontes e os
poços ou os fluxos através do fundo e/ou da superfície
livre dominarem a solução.
Transporte de calor
• No caso do calor os fluxos através do fundo
são normalmente pouco importantes.
• Pelo contrário os fluxos através da superfície
livre são essenciais (radiação, calor sensível e
calor latente.
5ª Aula Advecção - Difusão
Condições de fronteira na interface
com a atmosfera.
Condições de fronteira: Fronteiras
abertas
• Num canal as fronteiras de entrada e de saída do
escoamento (fronteiras abertas) requerem duas
condições de fronteira por via da difusão e da
advecção no caso de o escoamento estar a entrar.
• A difusão envolve uma segunda derivada e por
isso precisa de duas condições de fronteira (uma
na entrada e outra na saída).
• A advecção envolve uma primeira derivada e por
isso requer uma condição de fronteira: valor da
propriedade à entrada (ou fluxo advectivo à
entrada, pois estão relacionados pelo caudal).
Condições de Fronteira: Fronteiras
sólidas
• Nas fronteiras sólidas só pode haver fluxos
difusivos, que dependem da propriedade que se
está a estudar.
• No caso de sedimentos poderíamos ter erosão e
deposição (este último processo envolve a
velocidade de queda dos sedimentos).
• No caso do calor o fluxo difusivo através da
fronteira sólida é igual ao fluxo que se propaga
através do solo. Admitindo que esse fluxo é baixo,
então poderemos admitir que os fluxos através
das paredes sólidas são desprezáveis.
Condições de fronteira: fluxos através
da superfície livre
• Os fluxos através da superfície livre dependem
também da propriedade que estamos a
considerar.
• No caso de gases e de vapores dependem das
pressões parciais na atmosfera e na água. Na
grande maioria das propriedades são nulos, mas
no caso do calor são determinantes.
• Poderemos ter fluxos de calor latente, sensível e
fluxos de calor por radiação directa do sol, difusa
da atmosfera e ainda radiação da água para a
atmosfera.
Fluxo de calor latente
• Depende da temperatura da água e da
humidade relativa do ar. No modelo MOHID é
calculado como:
Fluxo e calor sensível
Fluxo de calor por radiação
• Ver:
• Brock, T. D. (1981) - Calculating solar radiation
for ecological studies. Ecological Modelling.