V-dfc-SistemaAlgebrico

Prévia do material em texto

Angela Nieckele – PUC-Rio
SOLUÇÃO DE SISTEMA 
ALGÉBRICO
O conjunto de equações algébricos 
ap fp= anb fnb + b
dá origem a um sistema algébrico que pode ser 
representado pelo produto de uma matriz de 
coeficiente A e um vetor incógnita f igual a um 
vetor conhecido b
A f = b
1
Angela Nieckele – PUC-Rio
ME2
Maior parte do tempo de uma simulação por 
volumes finitos, elementos finitos, diferenças 
finitas ou outro método numérico é gasto na 
resolução do sistema de equações obtido com a 
discretização.
Necessidade de métodos robustos e rápidos
Para resolver esse sistema pode-se utilizar
MÉTODOS DIRETOS
MÉTODOS ITERATIVOS
Angela Nieckele – PUC-Rio
3
1. Métodos Diretos
Consiste em manipular a matriz de coeficientes A 
e obter de uma vez o campo da variável de 
interesse.
A f = b f = A-1 b
A solução exata (a menos de erros de 
truncamento do computador) é determinada após 
um número finito de operações
Angela Nieckele – PUC-Rio
4
Para situações 1-D, a matriz de coeficientes A com 
muita freqüência só possui três diagonais 
diferentes de zero, sendo possível desenvolver 
algoritmos eficientes, que requerem pouco espaço 
de memória, como por exemplo o Algoritmo de 
Thomas (TDMA).
-aW aP -aE f = b
Angela Nieckele – PUC-Rio
5
Para situações 2-D, a matriz de coeficientes A com 
muita freqüência só possui cinco diagonais 
diferentes de zero. No caso 3D, sete diagonais 
diferentes de zero
-aS …0…-aW aP -aE …0… -aN f = b
Angela Nieckele – PUC-Rio
6
Para situações multi-dimensionais, em geral os 
algoritmos de inversão da matriz A de coeficientes 
são muito caros e requerem muito espaço de 
memória e tempo de execução.
 Requer mais memória de armazenamento
 Grande esforço computacional
 Mais robusto e Mais rápido 
Angela Nieckele – PUC-Rio
7
2. Métodos Iterativos
Fornece uma seqüência de soluções aproximadas que 
convergem quando o número de passos tende a infinito
Os problemas típicos de Fenômenos de transporte são 
problemas não lineares e portanto um procedimento 
iterativo deve ser utilizado. Já que o sistema de equações 
algébricas deverá ser resolvido diversas vezes atualizando 
os coeficientes, não vale a pena o esforço para inverter a 
matriz de coeficientes A e obter diretamente o campo de 
temperatura, pois o mesmo deverá ser corrigido, já que os 
coeficientes estarão errados.
 Menor necessidade de memória de armazenamento
 Problemas de convergência
Angela Nieckele – PUC-Rio
8
MÉTODOS DIRETOS
SISTEMAS TRIANGULARES












=





























nnnn
nn
nn
nn
b
b
b
b
x
x
x
x
u
uu
uuu
uuuu






3
2
1
3
2
1
313
21222
1111211
000
00
0
Se niuii ,,2,1,0 = as incógnitas podem ser facilmente calculadas
ii
n
ik
kkii
i
nn
nnnn
nnnnnnnn
nn
n
n
u
xub
x
u
xub
xbxuxu
u
b
x

=




=

==
=
1
,
1,1
,11
11,111,1
 :i linha
; :1-n linha
; :n linha

RETROSUBSTITUIÇÃO
Angela Nieckele – PUC-Rio
9
Se a matriz for triangular inferior:












=


























 nnnnannnn
b
b
b
b
x
x
x
x
llll
ll
ll
l






3
2
1
3
2
1
21
3231
2221
11
00
00
000
A solução é calculada da seguinte forma:
ii
i
ik
kkii
i
l
xlb
x
u
xub
xbxuxu
l
b
x


=

=

==
=
1
,
22
11,22
2222,211,2
11
1
1
 :i linha
; :2 linha
; :1 linha

SUBSTITUIÇÃO A FRENTE
NÚMERO DE OPERAÇÕES: 
=
=
n
i
nnnnin
1
2
2
1
)1(
2
1
)1(
Angela Nieckele – PUC-Rio
Eliminação de Gauss
 Consiste em rearranjar a Matriz A em uma matriz 
triangular superior. A solução é obtida por substituição 
regressiva
 Exemplo do procedimento:
Primeiro passo
o linha 2 – linha 1 multiplicada por ( a21 / a11 )
o Procedendo de forma análoga para as 
outras linhas, elimina-se todos os 
coeficientes da primeira coluna com 
exceção da primeira linha.
10
















=
































5
4
3
2
1
5
4
3
2
1
5554535251
4544434241
3534333231
2524232221
1514131211
b
b
b
b
b
aaaaa
aaaaa
aaaaa
aaaaa
aaaaa
f
f
f
f
f
nibmbb
niamaa
ni
a
a
m
iii
jiijij
i
i
,,,;
,,,;
,,,;
)(
)(



32
32
32
11
2
11
2
11
1
1
==
==
==
Angela Nieckele – PUC-Rio
11























































































=


























































































































































































1
11
1
1
11
31
3
1
11
21
2
1
3
2
1
1
11
1
0
11
11
1
1
1
11
31
312
11
31
32
0
11
11
31
31
1
11
21
212
11
21
22
0
11
11
21
21
1131211
b
a
a
b
b
a
a
b
b
a
a
b
b
a
a
a
aa
a
a
a
a
a
a
aa
a
a
aa
a
a
a
a
a
a
aa
a
a
aa
a
a
a
aaaa
n
m
n
n
m
mn
m
m
nn
nn
n


  

  

  

f
f
f
f
o A seguir, manipula-se as linhas de
forma análoga para eliminar os
coeficientes da segunda coluna
com exceção da primeira e
segunda linha e assim por diante
nibmbb
niamaa
ni
a
a
m
iii
jiijij
i
i
,,,;
,,,;
,,,;
)()()(
)()()(
)(
)(



43
43
43
2
22
23
2
22
23
2
22
2
2
2
==
==
==
Angela Nieckele – PUC-Rio
12
oA matriz resultante será do seguinte tipo





















=





































5
4
3
2
1
0000
000
00
0
f
f
f
f
f
fpode ser obtido por substituição regressiva.
Os elementos são denominados de Pivots
 O lado direito do sistema de equações é modificado da mesma 
forma que os coeficientes das equações
 Melhor tratar o sistema na forma matricial, com o lado direito do 
sistema sendo a coluna n+1 da matriz, conforme mostrado a seguir
)(
,
)()(
,,,
1
11
2
22
1
11


n
nn
aaa 
Angela Nieckele – PUC-Rio
13


















nnnannnn
nn
nn
nn
b
b
b
b
aaaa
aaaa
aaaa
aaaa






3
2
1
21
3133231
2122221
1111211niba ki
k
ni ,,2,1,
)()(
1, ==
end
end
end
1,1For 
,1For 
1,1For 
)()()1(
)(
)(
k
kjik
k
ij
k
ij
k
kk
k
ik
ik
amaa
nkj
a
a
m
nki
nk
=
=
=
=
=

end
end
1For 
0
11For 
1
)(
)(
,
)(
*
,
,,
i
ii
i
ni
i
k
i
ik
a
suma
asumsum
nik
sum
ni

=
=
=
=
=

f
f
ALGORÍTMO
ELIMINAÇÃO RETROSUBSTITUIÇÃO
Angela Nieckele – PUC-Rio
14
 O maior custo computacional ocorre no processo de eliminação
 Supor que o tempo de cada operação seja de 1 microsegundo
O tempo em segundos de cada parte do algoritmo é mostrado abaixo
st 610=
NÚMERO DE OPERAÇÕES:


=

1
1
3
3
1)1)((
n
k
nknkn 
=
=
n
i
nnni
1
2
2
1
2
1 )1()1(
ELIMINAÇÃO RETROSUBSTITUIÇÃO
n Eliminação Retrosubstituição
10 0.0050 s 0.0008 s
100 5 s 0.075 s
1000 5000 s 7.5 s
número de multiplicações é proporcional a N3
erros de truncamento se acumulam para sistemas grandes
Angela Nieckele – PUC-Rio
15
Caso o coeficiente da diagonal principal (pivot) vier a se anular, cuidados 
especiais precisam ser tomados:
 - uma solução é trocar a ordem das equações
 Para evitar falha catastrófica (divisão por zero) ou resultados errados
é necessário fazer uma escolha criteriosa dos PIVOTS usados na eliminação
 PIVOTAMENTO PARCIAL PIVOTAMENTO COMPLETO
PIVOTAMENTO PARCIAL
No passo k do processo de eliminação
k
r
k
rk
nikaa
r
k
ik
k
rk
 e linhasTrocar 
,max
que talinteiromenor o como Escolher 
)()(

=

PIVOTAMENTO COMPLETO
No passo k do processo de eliminação
k
r
k
skrk
njikaa
sr
k
ij
k
rs
 e colunas e , e linhasTrocar 
,,max
que talinteiros menores os como e Escolher 
)()(

=

s
Angela Nieckele – PUC-Rio
16
 A Eliminação Gaussiana deve ser feita sempre com PIVOTAMENTO
para garantir estabilidade do método
 Na grande maioria dos casos, PIVOTAMENTO PARCIAL é suficiente
e deve ser usada no lugar de PIVOTAMENTO COMPLETO
 PIVOTAMENTO COMPLETO não é muito usado devido ao grande
tempo computacional gasto no processo de busca do pivot.
 PIVOTAMENTO não é necessário em dois casos particulares
MATRIZ DIAGONAL DOMINANTE
MATRIZ SIMÉTRICA E POSITIVA-DEFINIDA


=
=
n
ij
j
ijii niaa
1
.,,2,1, 
0xAxxAA == ,0 e)( Tjiij
T aa
Angela Nieckele – PUC-Rio
17
O método apresentado para resolver as matrizes com 3 diagonais 
adjacentes (Algoritmo de Thomas) é um caso particular do método acima, 
otimizado para não se efetuar contas desnecessárias com os coeficientes 
nulos. Para o caso de 5 diagonais adjacentes, também pode-se otimizar o 
método de eliminação de Gauss de forma análogo ao utilizado para a 
matriz tri-diadonal.
SOLUÇÃO DE SISTEMA ALGÉBRICO 
PENTA-DIAGONAL
ai i = bi i+1 + ci i-1 + ei i+2 + fi i-2 + di
(I)
sendo c1= 0 ; f1= 0 ; f2= 0 
bN =0 ; eN = 0 ; eN-1= 0
A matriz de coeficientes 
penta-diagonal é
f1
f2
.
.
..
.
.
.
.
.fN
Angela Nieckele – PUC-Rio
18
Para derivar este algoritmo deve-se proceder da mesma forma que procedemos 
para derivar o algoritmo para matrizes tri-diagonais.
- Vamos assumir que podemos determinar o valor de i em função de i+1 e
i+2 
i = Ri i+2 + Pi i+1 + Qi (II)
Rescrevemos a equação (II) para i-1 e i-2
i-1 = Ri-1 i+1 + Pi-1 i + Qi-1 (III)
i-2 = Ri-2 i + Pi-2 i-1 + Qi-2 (IV)
Vamos agora substituir (III) em (IV)
i-2 = Ri-2 i + Pi-2 (Ri-1 i+1 + Pi-1 i + Qi-1) + Qi-2 (V)
Agora vamos substituir a equação (V) e (III) em (I)
ai i = bi i+1 + ci (Ri-1 i+1 + Pi-1 i + Qi-1) + ei i+2 +
+ fi (Ri-2 i + Pi-2 (Ri-1 i+1 + Pi-1 i + Qi-1) + Qi-2 ) + di
Angela Nieckele – PUC-Rio
19
Rearrumando a equação acima obtemos explicitando i e comparando com a 
equação (II) 
i = Ri i+2 + Pi i+1 + Qi
temos que
)(
)(
2121
21




=
iiiiiii
iiiii
i
RPPfPca
PfcRb
P
)(
)(
2121
1221




=
iiiiiii
iiiiiii
i
RPPfPca
QPQfQcd
Q
)( 2121  
=
iiiiiii
i
i
RPPfPca
e
R
Angela Nieckele – PUC-Rio
20
Procedimento para matriz penta-diagonal
calcula-se ; ; 
; ; 
Usar as relações recursivas para obter Pi ; Qi e Ri para 
3  i  N
especificar N = QN e N-1 = PN-1 N + QN-1
Usando i = Ri i+2 + Pi i+1 + Qi , obter N-2 ; ..... ; 2 ; 1
)(
)(
2121
21




=
iiiiiii
iiiii
i
RPPfPca
PfcRb
P
)(
)(
2121
1221




=
iiiiiii
iiiiiii
i
RPPfPca
QPQfQcd
Q
)( 2121  
=
iiiiiii
i
i
RPPfPca
e
R
1
1
1
a
b
P =
1
1
1
a
d
Q =
1
1
1
a
e
R =
122
122
2
Pca
Rcb
P


=
122
122
2
Pca
Qcd
Q


=
122
2
2
Pca
e
R

=
Angela Nieckele – PUC-Rio
21
Decomposição LU
Resolver o sistema  [ A ] f = b
Todo matriz não singular pode ser decomposta como o produto de uma 
triangular inferior L e uma triangular superior U 
[ A ] = [L] [U] então [ L ] [ U ] f = b
Uma vez feita a decomposição, a solução do sistema fica reduzida a solução 
de dois sistemas triangulares:
















=
















==
 nn
n
n
n
nnnn u
uu
uuu
uuuu
lll
ll
l
000
00
0
1
01
001
0001
333
22322
1131211
121
3231
21










U;L;LUA
,
Angela Nieckele – PUC-Rio
22
Define-se um vetor auxiliar z 
Procedimento: 
1o. Resolver [ L ] z = b por substituição progressiva, 
determinando z
2o. Resolver [ U ] f = z por substituição regressiva, 
determinando
   bUL =

z
f





















=





































5
4
3
2
1
0
00
000
0000
z
z
z
z
z





















=





































5
4
3
2
1
0000
000
00
0
f
f
f
f
f
Angela Nieckele – PUC-Rio
23
=
















5554535251
4544434241
3534333231
2524232221
1514131211
aaaaa
aaaaa
aaaaa
aaaaa
aaaaa
















1
01
001
0001
00001
54535251
434241
3231
21
llll
lll
ll
l
















55
4544
353433
25242322
1514131211
0000
000
00
0
u
uu
uuu
uuuu
uuuuu
1111 au = 1212 au =
11
21
21
u
a
l =
11
31
31
u
a
l =
12212222
ulau =
22
123131
32
u
ula
l

=
, 
, ...........
, 
, ...........
, ...........
, ..............
13212323 ulau =, 
jipara
u
ula
l
jj
j
k
jkkiji
ji 

=


=
1
1 jiparaulau
i
k
jkkijiji = 

=
1
1
abaixo da diagonal acima da diagonal
Angela Nieckele – PUC-Rio
Observação:
1) Note que nem sempre podemos decompor uma matriz A em uma 
matriz triangular inferior e outra superior. Logo nem sempre o 
método de decomposição L U pode ser utilizado.
2) Assim como o método de Eliminação de Gauss também pode ter 
problemas, quando os coeficientes da diagonal principal são 
pequenos, e artifícios especiais devem ser introduzidos para evitar 
problemas
3) Os problemas típicos de fenômenos de transporte são problemas 
não lineares e portanto um procedimento iterativo deve ser 
utilizado. Já que o sistema de equações algébricas deverá ser 
resolvido diversas vezes atualizando os coeficientes, muitas vezes 
não vale a pena o esforço para inverter a matriz de coeficientes A
e obter diretamente o campo de temperatura, pois o mesmo 
deverá ser corrigido, já que os coeficientes estarão errados.
24
Angela Nieckele – PUC-Rio
Métodos Iterativos
Os métodos iterativos consistem em resolver o sistema algébrico 
utilizando duas etapas: 
(a) prever, 
(b) corrigir.
Os métodos iterativos podem ser:
 ponto-a-ponto: Os métodos iterativos ponto-a-ponto consistem em 
visitar um ponto nodal do domínio de cálculo de cada vez e estimar 
o valor da temperatura para o ponto nodal em função da 
temperatura estimada dos pontos vizinhos. Repete-se o 
procedimento para todos os pontos nodais. Como exemplo, pode-se 
citar o método de Gauss Seidel.
 em blocos: Os métodos iterativos em blocos consistem em resolver 
simultaneamente um conjunto de pontos nodais do domínio, em 
função da temperatura estimada dos pontos vizinhos. O 
procedimento é repetido até que todos os blocos do domínio 
tenham sido resolvidos. Como exemplo, pode-se citar o algoritmo 
TDMA linha por linha.
25
Angela Nieckele – PUC-Rio
Método de Jacobi
Para use
onde é o valor dos f vizinhos da última iteração. Todos os pontos
nodais são percorridos dessa maneira, atualizando os valores de f.
Este método apresenta uma convergência muito lenta, não sendo
muito utilizado.
 Exemplo:
T1 = 0,5 T2 + 1,5 T2 = 0,25 T1 + 0,5
Converge lentamente para a solução correta
26
baa nbnbii =  ff
i
nbnb
i
a
ba 
=

*f
f
*
nbf
T1 
0 
 1,5 1,750 1,9375 ...... 2,0 
T2 0 0,5 0,875 0,9375 ...... 1,0 
 
0
Angela Nieckele – PUC-Rio
Método de Gauss-Seidel
Este é um algoritmo iterativo ponto-a-ponto para resolver o sistema de
equações algébricas de situações uni- ou multi-dimensionais.
Para use
Onde é o último valor disponível na memóriado computador do
ponto vizinho. Todos os pontos nodais são percorridos dessa maneira,
atualizando os valores de f.
 Exemplo:
T1 = 0,5 T2 + 1,5 T2 = 0,25 T1 + 0,5
Converge para a solução correta
27
baa nbnbii =  ff
i
nbnb
i
a
ba 
=

*f
f
*
nbf
 
T1 0 1,5 1,938 1,990 ...... 2,0 
T2 0 0,875 0,984 0,998 ...... 1,0 
 
Angela Nieckele – PUC-Rio
 Exemplo:
T1 = 4 T2 - 2 T2 = 2 T1 - 3
a solução diverge. 
Note que as equações são as mesmas que o exemplo anterior, porém 
rearranjadas
28
T1 0 -2 -30 -254 -2046 ........... 
T2 0 -7 -63 -511 -4095 .......... 
 
Angela Nieckele – PUC-Rio
29
O Critério de Scarborough 
 
 O método de Gauss-Seidel não converge sempre. 
 
 As possibilidades de convergência do método de Gauss Seidel podem 
ser determinadas com referência ao critério de Scarborough, que é 
uma condição suficiente para a convergência do método de Gauss 
Seidel. É preciso que: 
 
1

P
nb
a
a
 para todas as equações e 
 < 1 para pelo menos uma equação 
nb são os vizinhos 
com valor 
desconhecido
Angela Nieckele – PUC-Rio
 Exemplo:
 As quatro regras básicas apresentadas, dão origem a equações de 
discretização que satisfazem o Critério de Scarborough, uma vez 
que os coeficientes são positivos e ap   anb.
 Um SP positivo pode fazer com que ap <  anb
 Um anb negativo pode levar a ap = ( anb ) <  anb
violando o critério
 A condição de desigualdade do critério em pelo menos um ponto 
nodal é garantida pelas condições de contorno.
30
Angela Nieckele – PUC-Rio
Algoritmo TDMA linha por linha
 Este é um algoritmo iterativo por blocos para resolver o sistema de 
equações algébricas de situações multi-dimensionais.
• Selecione uma linha da malha.
• Suponha que os valores de f ao longo 
da linhas vizinhas são conhecidos, 
a partir de seus últimos valores.
• Resolva os valores de f ao longo da linha selecionada pelo método 
direto TDMA. 
• Repita este procedimento para todas as linhas em uma direção, 
repita, se desejado, para as linhas nas outras direções.
31
  
*b
baaaaa *
E
E
*
W
WSSNNPp ffff=f
  
*b
baaaaa *
S
S
*
NNEEWWPp ffff=f
Angela Nieckele – PUC-Rio
• No método de Gauss Seidel, as informações viajam um ponto nodal 
por iteração. No método TDMA linha-linha, a informação do 
contorno é transmitida de uma vez para o interior do domínio; 
consequentemente a convergência é mais rápida.
• Se os coeficientes em uma direção são muito maiores que aqueles 
nas outras direções, a TDMA ao longo da direção dos coeficientes 
dominantes levará a uma convergência muito mais rápida.
• Se Dx >> Dy e dx >> dy → aN e aS >>> aW e aE logo os 
valores ao longo da coluna terão um peso muito maior na obtenção 
da solução
• A direção da varredura também é importante. (Para escoamentos, 
a varredura deve ser de montante para jusante).
32
  
*b
baaaaa *
E
E
*
W
WSSNNPp ffff=f
  
*b
baaaaa *
S
S
*
NNEEWWPp ffff=f
Angela Nieckele – PUC-Rio
Método ADI
Alternating - Direction - Implicit
(Peaceman - Rachford) - 1955.
Este método adota uma formulação explícita em uma direção 
utilizando meio passo de tempo (Dt/2) ( ) e formulação 
implícita na outra, utilizando meio passo de tempo (Dt/2). Alterna-se 
as direções de uma iteração para o seguinte.
33
tao
P
DD /=
  
b
aaSaaa
aaSaaa
n
P
nn
nnn
P
PSNC
o
SSNN
WWEEPP
o
WE
fDff
fffD





 
=




 

2
2
212121 ///
  
b
aaSaaa
aaSaaa
n
P
nn
nnn
P
PWEC
o
WWEE
SSNNPP
o
SN
212121
111
2
2
/// 






 
=




 
fDff
fffD
Angela Nieckele – PUC-Rio
Método Fortemente Implícito 
(SIP = Strongly Implicit Procedure, Stone 1968)
 Para ilustrar este método, considere um sistema de equações 
algébricas com 
[A] u = C
[A] matriz de coeficientes ; u vetor desconhecido
C vetor conhecido
 Método Fortemente Implícito é uma técnica de fatorização da matriz
 O objetivo é substituir [A] por uma matriz [A + P] tal que a matriz 
modificada possa ser decomposta em matrizes triangulares, uma 
superior [U] e uma inferior [L].
34
Angela Nieckele – PUC-Rio
Método Fortemente Implícito 
 [A]  [A + P] tal que [A + P] = [L] [U] 
[A] u = C
[A + P] un+1 = C + [P] un
[A + P] = [L] [U]
[L] [U] un+1 = C +[P]un
definindo Vn+1 = [U] un+1. Obtemos um algoritmo de 2 passos
1o passo [L] Vn+1 = C + [P] un
2o passo [U] un+1 = Vn+1
passo 1 consiste de uma substituição para frente e o passo 2 para trás.
repetir até convergir 35
Angela Nieckele – PUC-Rio
Sobre-Relaxação e Sub-Relaxação
Se 
Definimos: 
onde é o valor de da iteração anterior, e  é o fator de 
relaxação.
 < 1: sub-relaxação.
 > 1: sobre-relaxação.
 A sobre-relaxação pode ser utilizada para acelerar a velocidade de 
convergência. Em geral, é utilizada junto com o método de Gauss-
Seidel, o qual por ser um método iterativo ponto-a-ponto, apresenta 
taxas de convergência muito baixas. Este método é muitas vezes 
denominado de SOR (Sucessive Over Relaxation).
36
b a ap nbnbp =  ff
~
*
ppp ) - (1 fff =
~
*
pf pf
Angela Nieckele – PUC-Rio
 Em problemas não lineares, é desejável reduzir as variações da 
variável dependente, então, sub-relaxação é recomendável.
 O fator de relaxação pode ser introduzido diretamente nas 
equações de discretização,
 resultando em 
37
*
p
p
nbnb
p ) - (1 
a
b a
 f
f
f 

=

*
p
p
nbnbp
p
 
a
 ) - (1 b a 
a
f

ff

= 
Angela Nieckele – PUC-Rio
Observações:
 Não existe regra geral para determinação de .
 Depende de: 
 Natureza do problema número de pontos nodais
 espaçamento
 método iterativo utilizado
 etc.
 Pode-se variar a durante o processo iterativo.
 Pode-se também adotar a diferentes para cada ponto nodal (não é 
recomendado)
38
Angela Nieckele – PUC-Rio
Convergência de problemas não lineares
 O desenvolvimento de um método numérico que respeite 
as quatro regras básicas, leva a conjunto de 
discretização que, para coeficientes constantes, possui 
convergência garantida ao utilizarmos métodos de 
solução do tipo ponto a ponto ou linha por linha.
 Se os coeficientes variam lentamente, ainda podemos 
garantir convergência. Fatores de sub-relaxação 
apropriados para as variáveis dependentes diminuem as 
variações das variáveis e portanto dos coeficientes.
39
Angela Nieckele – PUC-Rio
 Quantidades auxiliares como massa específica  e 
viscosidade m podem ser sub-relaxadas através de:
 = novo +(1-) velho onde  é o fator de 
relaxação.
 Ao sub-relaxarmos as quantidades apropriadas, 
podemos enfraquecer temporariamente o acoplamento 
entre o conjunto de equações. Por exemplo, a sub-
relação de  separa as equações de 
temperatura/concentração. A sub-relaxação da 
viscosidade turbulenta mt reduz o acoplamento entre as 
quantidades turbulentas e as equações do escoamento.
40
Angela Nieckele – PUC-Rio
 Quando o acoplamento entre as equações é devido à 
presença do termo de fonte (por exemplo, em 
convecção natural, o termo de empuxo acopla a 
equação e quantidade de movimento com a equação de 
energia), este também pode ser sub-relaxado. através 
de:
SC =  SC novo + (1 - ) SC velho onde  é o 
fator de relaxação.
SP =  SP novo + (1 - ) SP velho 
 Condições de contorno também podem ser sub-
relaxadas. Por exemplo, aumenta-se gradativamente a 
velocidade de um disco girante que impõe um “swirl” no 
escoamento.
41
Angela Nieckele – PUC-Rio
 Não existe nenhuma garantia que, para todas as não 
linearidades e acoplamentos, sempre conseguiremos 
convergência.
 Cada problema não linear e acoplado precisa de 
cuidados especiais.Alguns artifícios que podem auxiliar 
na obtenção de solução convergida:
(i) iniciar os cálculos com fatores de sub-
relaxação bem pequenos, aumentá-los gradativamente 
à medida que a solução converge
. (ii) utilizar como “chute inicial” a solução 
convergida de um caso mais simples. Por exemplo: 
inicializar o escoamento turbulento com a solução do 
escoamento laminar.
(iii) Aumentar gradativamente o parâmetro que 
governa o escoamento. Por exemplo: iniciar com 
número de Reynolds até obter o valor desejado. 42
Angela Nieckele – PUC-Rio
RELAXAÇÃO POR INÉRCIA
Rescreve-se a equação de discretização como
i  inércia Se i > 0  sub-relaxação
Se i < 0  sobre-relaxação
 Mesmos comentários feitos para  são válidos 
FORMULAÇÃO TRANSIENTE
se comporta como a inércia i
 Pode-se resolver um problema de regime permanente, utilizando a 
formulação transiente. Cada passo de tempo corresponde a uma 
iteração. Menor passo de tempo, maior sub-relaxação.
43
 =
*
)( PnbnbPP ibaia fff
 =
o
P
o
PcnbnbPp
o
Pnb aSaSaa fDffD )(
o
Pa
Angela Nieckele – PUC-Rio
 Este é um algoritmo para acelerar a convergência. Consiste em 
resolver as equações de conservação para um conjunto de blocos 
de forma a fornecer o nível correto da temperatura mais 
rapidamente para o interior do domínio.
 Um domínio 3-D é aproximado por um domínio 2-D e um domínio 
2-D é aproximado por um domínio 1-D. A solução do domínio de 
ordem inferior é mais fácil de ser obtida, fornecendo uma primeira 
estimativa para o domínio de ordem superior. 
44
Algoritmo de Correção por Blocos 
Angela Nieckele – PUC-Rio
 Sem perda de generalidade, vamos considerar uma situação bi-
dimensional.
 Considere que
 Substitua a definição acima na equação de discretização e some 
todos os blocos ao longo da coluna i.
 A equação resultante será 45
i.ji,j-1i,j1i,ji,ji-1,ji,j1,jii,ji,ji,j d T f T e T c T b Ta = 
*
i,jii,j T T T =
 











 









 








 









 =



 


j
i,j
j
*
i,j-1ii,j
j
*
1i,jii,j
j
*
i-1,ji-1i,j
j
*
1,ji1ii,j
j
*
i,jii,j
d TTf
 T Te T Tc
 TTb ) T T(a
Angela Nieckele – PUC-Rio
46
 





 
=


j
*
i,ji,ji,j
*
i,j-1i,j
*
1i,ji,j
*
i-1,ji,j
*
1,jii,j
j
i,ji-1
j
i,j1i
j
i,ji,ji,ji
T - a d T f T e T c Tb
 c T b T - f - ea T
 BLC T BLM T BLP TBL i-11ii = 
 
 


=
==
=

j
*
i,ji,ji,j
*
i,j-1i,j
*
1i,ji,j
*
i-1,ji,j
*
1,jii,j
j
i,j
j
i,j
j
i,ji,ji,j
T - a d T f T e T c Tb BLC 
c BLM b BLP 
 - f - ea BL 
;
logo
onde
Angela Nieckele – PUC-Rio
 A equação resultante para os blocos é 1-D, podendo ser 
resolvida pelo método direto TDMA.
 Note que o termo BLC corresponde ao resíduo da 
equação original, logo, se a solução estiver convergida, 
o resíduo será nulo, a solução do sistema algébrico 
fornecerá valor nulo para , consequentemente, 
nenhuma nova correção será efetuada.
 O algoritmo de correção por blocos, só fornece uma 
correção média para o conjunto de blocos, portanto, 
após sua utilização o algoritmo TDMA linha por linha 
deverá ser utilizado para determinar a distribuição da 
temperatura ao longo dos blocos.
 Pode-se inverter a direção do algoritmo de correção por 
blocos.
47
*
ijjij
T T T =
Angela Nieckele – PUC-Rio
 Como o algoritmo de correção por blocos corrige da mesma 
maneira todos os volumes de controle dentro de um bloco, em 
algumas situações particulares, quando a variável dependente 
possui restrições físicas, este algoritmo poderá fornecer resultados 
indesejáveis, não devendo ser utilizado.
 Exemplo: fração em massa de uma espécie química,mℓ . De 
acordo com sua definição, temos que entre zero e um. Se durante 
o processo iterativo, está condição for violada, o processo irá 
divergir. O algoritmo de correção por bloco garante que a 
conservação para o bloco (conjunto de volumes de controle) seja 
satisfeita, mas um volume de controle pode estar mais errado do 
que outro. Logo, ao corrigir todos os volumes de controle de um 
bloco com a mesma constante, a condição acima poderá ser 
violada.
48
Angela Nieckele – PUC-Rio
49
MÉTODO DE MULTIGRID
Métodos iterativos ponto a ponto (Gauss Seidel) , ou 
bloco por bloco (TDMA linha por linha) removem 
rapidamente os erros locais (alta freqüência) da solução, 
porém os erros globais (baixa freqüência) são reduzidos a 
uma taxa inversamente proporcional ao tamanho da 
malha. Conseqüentemente, para um número elevado de 
nós, a taxa de redução do resíduo torna-se extremamente 
baixa.
Técnicas de Multigrid permitem que o erro global seja 
tratado utilizando uma seqüência de malhas mais 
grosseiras. O método é baseado no princípio que o erro 
global (baixa freqüência) existente em um malha fina pode 
ser representado numa malha grosseira, onde este torna-
se acessível como um erro local (de alta freqüência).
Angela Nieckele – PUC-Rio
50
Considere o conjunto de equações discretizadas
lineares
A f e + b = 0
onde f e é a solução exata. Se f é a solução 
aproximada, antes da solução convergir haverá um 
resíduo r
A f + b = r
Deseja-se corrigir y com f, tal que a solução exata seja 
dada por
f e = f + y
O conceito básico do Método de Multigrid
Angela Nieckele – PUC-Rio
51
Substituindo na equação de discretização
A (f  y ) + b = 0
A y  (A f + b)= 0
logo
A y  r = 0
A qual é a equação para o termo de correção em função 
do operado A do nível fino original e do resíduo r.
Assumindo que os erros locais (alta freqüência) foram 
suficientemente amortecidos pelo esquema de relaxação 
do nível fino, a correção y será suave e portanto mais 
eficientemente resolvida no próximo nível de malha 
grosseira.
Angela Nieckele – PUC-Rio
52
Resolver as correções no nível grosseiro, requer transferir 
o resíduo do nível fino (restrição), computar as correções, e 
então transferir de volta do nível grosseiro (prolongação).
Podemos escrever as equações para as correções do nível 
grosseiro yH como
A H yH + R r = 0
onde A H é o operador no nível grosseiro e R é o operador 
de restrição responsável por transferir o resíduo do nível fino 
para o nível grosseiro. Após a solução da equação acima, a 
atualização da solução no nível fino é dada por
fnew = f + P yH
sendo P o operador de prolongação usado para transferir 
as correções do nível grosseiro para o nível fino.
Restrição e Prolongação
Angela Nieckele – PUC-Rio
53
MÉTODO DE MULTIGRID
Pode-se 
ainda 
combinar o 
ciclos V e 
W (ciclo 
flexível)
Ciclos do Multigrid
Angela Nieckele – PUC-Rio
54
Multigrid tradicional
- resolve equação para h, 2h, 4h
- interpola e extrapola as soluções para troca 
de níveis
Multigrid Correção Aditiva (Multigrid Algébrico)
- utiliza solução em malha grosseira para corrigir 
solução da malha mais fina
Métodos Multigrid
h dependem a
banb ap
nb
nbp
=
= ff
Angela Nieckele – PUC-Rio
55
MULTIGRID DE CORREÇÃO ADITIVA 
(Hutchinson, Galpin, Raithby, 1988)
Relação entre malha fina e malha grosseira.
Angela Nieckele – PUC-Rio
56
k
*
ijij δ T T =
 im, ip , i a a im, ip , ia a
 jm,jp , j a a jm,jp , j a a
b δa δaδa δa δa
i
S
ijm
S
kl
i
N
ijp
N
kl
j
W
imj
W
kl
j
E
ipj
E
kl
kl-1k, 
s
1k, 1k, 
N
klk-l, 
W
kl1, lk
E
klkl
P
kl
====
====
=



ˆˆ
ˆˆ
ˆˆˆˆˆˆ

jm, jpjim, ip, ia a
j i
P
ij
P
k l
===  SUM4, - SUM3 - SUM2 - SUM1 ˆ
Angela Nieckele – PUC-Rio
57
)a - a a a a (b b̂
jp 1,jm j ip, im, i , a SUM4
1-jp jm,j ip, im, i , a SUM3
jpjm, j ip, 1,im i ,a SUM2
jpjm,j 1,-ip im,i ,a SUM1
*P
ij
*
1-ij
S
ij
*
1ij
N
ij
*
1
W
ij
*
1
E
ij
j i
ijkl
j i
S
ij
j i
N
ij
j
W
ijj i
E
ij
ijjiji
i
fffff =
===
===
===
===





Angela Nieckele – PUC-Rio
58
Ppnbnb φ b - a φ a =
=


 Res
Res R
R R
ij
iji
1-ii 