Anais do 15O Encontro de Iniciação Científica e Pós-Graduação do ITA – XV ENCITA / 2009
Instituto Tecnológico de Aeronáutica, São José dos Campos, SP, Brasil, Outubro, 19 a 22, 2009.
MODELO COMPUTACIONAL PARA ANÁLISE DA CONDUÇÃO DE
CALOR
João Batista Caldas Neto
Divisão de Engenharia Mecânica, Instituto Tecnológico de Aeronáutica, São José dos Campos, SP, Brasil, 12228-900
Bolsista PIBIC-CNPq
[email protected]
Ezio Castejon Garcia
Divisão de Engenharia Mecânica, Instituto Tecnológico de Aeronáutica, São José dos Campos, SP, Brasil, 12228-900
[email protected]
Gustavo Adolfo Ronceros Rivas
Divisão de Engenharia Mecânica, Instituto Tecnológico de Aeronáutica, São José dos Campos, SP, Brasil, 12228-900
[email protected]
Resumo.
Este projeto trata do desenvolvimento de um modelo computacional para análise da condução de calor em sólidos. O
modelo atende as diversas condições de contorno como temperaturas impostas, fluxos e convecção. Foram estudados primeiramente
casos unidimensionais, dependentes ou não do tempo. Numa segunda etapa o modelo foi expandido para casos bi e tridimensionais,
sendo utilizado o método de volumes finitos na discretização das equações. O método TDMA foi utilizado para solução dos sistemas
de equações gerados. A malha computacional utilizada foi uniforme. Comparações e validações do modelo foram feitas utilizando o
código computacional Fluent. Os resultados do modelo para exemplos 1D apresentaram erros da ordem de 2% para malhas
grosseiras, em relação à solução analítica. Para aplicação bidimensional, a discordância dos valores foi da ordem de 0,1% em
comparação com malhas idênticas simuladas em Fluent. Os casos tridimensionais apresentaram distorções geométricas mais
acentuadas, sendo necessário um refinamento, que não foi realizado.
Palavras chave: condução de calor, simulação numérica, modelo computacional, TDMA.
1. Introdução
O estudo e o desenvolvimento de rotinas computacionais utilizadas para análise da condução de calor é importante
para o entendimento e otimização de programas para casos particulares, uma vez que os códigos comerciais apresentam
códigos fechados. Além disso, uma análise bem sucedida de um problema sem solução analítica, através de métodos
numéricos, é mais fiel à realidade quanto maior a experiência do operador em identificar a necessidade de
considerações a respeito de simetria e refinamento da malha apenas em regiões específicas da geometria, por exemplo,
de modo a garantir uma precisão razoável sem comprometimento do tempo de cálculo.
Visando um entendimento mais profundo dessas habilidades no campo da transferência de calor e de sua
programação atrelada, o trabalho desenvolvido teve como caráter principal a iniciação ao projeto científico e a
introdução às técnicas mais simples de simulação em CFD através do método de volumes finitos. O projeto foi
desenvolvido com formulações próprias das equações discretizadas a partir de problemas simples unidimensionais
evoluindo até a resolução tridimensional de distribuição de temperatura em geometrias de paralelepípedos com
condições de contorno mistas (fluxo e temperatura imposta) de maneira variada e volumes de controle com dimensões
arbitrárias, resultando na construção de uma rotina computacional para a solução de problemas genéricos.
O estudo das equações de difusão do calor e a discretização de equações por método dos volumes finitos foram
suportados pela construção de rotinas em MATLAB com a análise de desvios entre soluções numéricas e analíticas,
através de gráficos. A partir dos problemas bidimensionais, um método mais simples e rápido de resolução de uma
classe particular de sistemas lineares foi utilizado. Com o desenvolvimento da rotina TDMA de resolução de sistemas
lineares algébricos com matriz tridiagonal de coeficientes, resolveu-se problemas com condições de contorno
constantes, inicialmente em geometrias planas e posteriormente em paralelepípedos com a extensão do método
iterativo. Em um dos exemplos bidimensional, a equação de Laplace foi resolvida pelo método de separação de
variáveis, obtendo a solução como uma serie de Fourier. Nos dois casos bidimensionais mostrados, as respectivas
simulações em Fluent estão devidamente apresentadas.
Um trabalho semelhante já foi desenvolvido para um caso particular 3D (Nascimento, Belo, de Lima, 2002), e foi
estudado a título de revisão bibliográfica.
Nomenclatura
k
T
Condutividade térmica do material
Temperatura
Anais do XV ENCITA 2009, ITA, Outubro, 19-22, 2009
,
q
ρ
c
t
x
h
L
A
Ac
n
Fonte ou Sorvedouro de calor
Densidade do material
Calor específico do material
Tempo
Coordenada de espaço
Constante de convecção
Dimensão da geometria na direção da condução unidimensional do calor
Área da secção de condução unidimensional do calor
Área de resfriamento convectivo
Número de nós da malha
2. Equações Governantes
 
T
.(k T )  q   c
t
(1)
3. Análise da condução unidimensional do calor
Embora a simplicidade da condução unidimensional do calor não exija um tratamento matemático muito
complicado, muitos problemas em engenharia podem ser razoavelmente aproximados para estes casos. Como um
resumo do trabalho desenvolvido nesta primeira etapa, apenas três análises importantes serão citadas: Caso (I), a
propagação do calor em regime permanente em uma barra delgada homogênea de material isotrópico isolada com seção
transversal constante sem geração de calor; Caso (II), condução do calor em regime permanente entre as faces de uma
superfície infinita de material isotrópico a temperaturas constantes com espessura constante e geração uniforme de
calor; e Caso (III), problemas de transferência unidimensional em situações com termo fonte de calor variando com a
temperatura de forma linear e com uma extremidade isolada (aleta dissipadora).
Caso (I): A equação governante do problema (Eq. 1) se reduz a:
d 2T
0
dx 2
(2)
Integrando em um infinitésimo de volume qualquer:
dT
dx

final
dT
dx
0
(3)
inicial
Para discretização da malha, adotou-se a convenção mostrada na Fig. 1, como citado em (Veersteg & Malalasekera,
1995).
Figura 1. Convenção para malha de um volume de controle unidimensional.
Dessa forma, por diferenças finitas, para pontos intermediários da malha:
TE  TP TP  TW
T T

 0 , logo TP  E W
 xPE
 xWP
2
Para volumes de controle nas extremidades da malha:
(4)
Anais do XV ENCITA 2009, ITA, Outubro, 19-22, 2009
,
TP 
TE  2TA
3
e
TW  2TB
3
TP 
(5)
Onde o índice A indica a extremidade esquerda do domínio do problema e B, a extremidade direita.
Com a Eq. (4) e a Eq. (5), tem-se um sistema genérico para a distribuição de temperatura em um problema de
condução do calor unidimensional nas condições citadas.
3T1  T2  2TA
 3 1 0 0
T  2T  T  0
 1 2 1 0
2
3
 1

...
 0 1 2 1


Tn 1  2Tn  Tn 1  0  0 0 1 2
 ... ... ... ...
...


 0 0 0 0
TN 1  3TN  2TB
... 0   T1   2TA 
... 0   T2   0 
... 0   ...   ... 

   
... 0   Tn   ... 
... 1  ...   0 

  
1 3  TN   2TB 
(6)
Cabe salientar o fato da dependência exclusiva das temperaturas nas fronteiras para determinação da distribuição de
temperatura. Os resultados da resolução do sistema por um método numérico obviamente coincide com a solução
analítica do problema que é simples e dada pela Eq. (7).
T ( x) 
TB  TA
x  TA
L
(7)
Caso (II): A equação governante do problema (Eq. 1) se reduz a:
 dT 
 qV  0
 kA

 dx  e  w
(8)
Trata-se do problema anterior com termo fonte. Por volumes finitos, a equação discretizada é dada pela Eq. (9).
 TE  TP
 kA  x

  TP  TW
   kA  x
 

  qA x  0

(9)
Assim, para pontos intermediários da malha:
2
q  x 
1
TP   TW  TE 
2 
k




(10)
E para pontos nas extremidades da malha:
2
q  x 
1
TP   TE  2TA 
3 
k
2

q  x 
1
 e TP   TW  2TB 

3 
k





(11)
Portanto, com a Eq. (10) e a Eq. (11), constroe-se o sistema da Eq. (12), que representa um modo mais geral de
obtenção de uma solução também para o problema anterior, Caso (I).
 3 1 0 0
 1 2 1 0

 0 1 2 1

 0 0 1 2
 ... ... ... ...

 0 0 0 0
... 0   T1   2TA 
... 0   T2   0 
2
... 0   ...   ...  q  x 

   
... 0   Tn   ... 
k
... 1  ...   0 

  
1 3  TN   2TB 
1
1

1

1
1

1
(12)
Anais do XV ENCITA 2009, ITA, Outubro, 19-22, 2009
,
Observa-se, portanto, que existe a dependência no termo fonte da geometria da malha. Então, pôde-se fazer a
comparação entre a solução numérica e a solução analítica mostrada na Eq. (13).
  T  TA  q  L  x  
T x   B

 x  TA
L
2k 

(13)
A título de exemplo, resolveu-se o problema para as condições: L=0,02m; q=1000KW/m3; k=0,5W/m.K; n=5;
TA=100oC; TB=200oC.
Figura 2. Soluções numérica e analítica para o Caso II, segundo as condições descritas.
A observação principal a respeito do desvio em porcentagem da solução numérica é feita observando um erro
mínimo numa posição intermediária, ocorrendo os maiores erros nas extremidades.
Caso (III): A equação governante do problema (Eq. 1) se reduz a:
kA
d  dT 

  q1T  q2  0
dx  dx 
(14)
Ou seja, o termo fonte é uma função linear da temperatura de acordo com a Lei de Newton para o resfriamento,
sendo:
q1  hAc e q2  hAcT
(15)
Integrando em um volume de controle:
 dT 
kA  A
  q1TP  q2  A x  0 

 dx  e  w
dT
dx
  q1TP  q2 
ew
x
kA
0
(16)
Pelo método das diferenças finitas, lineariza-se o termo com derivada:
q2 x
1
1
 2 q1 x 
  x  kA  TP   x TE   x TW  kA


(17)
Dados numéricos do problema analisado: L=1m; A=0,01m2; n=10; q1=-25W/m3.K; q2=500W/m3; k=100W/m.K;
TA=100oC.
Particularizando o problema, para as extremidades, utilizando as condições de contorno de Dirichlet para a
extremidade esquerda e de Newman para a extremidade direita (extremidade isolada termicamente), tem-se:
2

q  x 
3 1

kA

2

 q  x  2
q  x 
 TP  TE  2TA  2
e 1  1


kA
kA


2

q  x 
 TP  TW  2

kA

(18)
Anais do XV ENCITA 2009, ITA, Outubro, 19-22, 2009
,
Com a Eq. (17) e a Eq. (18), montou-se o sistema:
 q2  x 2


 2TA 
 kA

 I 1 0 .. 0 0   T1  
2


    q2  x 

T

1
J

1
...
0
0
2

 

kA
 0 1 J ... 0 0   ...  

...

   

T
...
...
...
...

1
0

 n 

...
 0 0 0 1 J 1  ...  


   
...

 0 0 0 0 1 K   TN  
2

q

x


2


kA


2

q  x 
Sendo: I   3  1

kA

2


q  x 
 , J  2 1


kA



 q  x 2
 e K  1  1


kA


(19)




(20)
A solução analítica para este problema é dada pela Eq.21.
q2
 Ta
 q
 q
q1
T ( x) 
cos  1  x  L    2
 kA
 q
 q 

 1
cos  1 L 
kA


(21)
Neste caso, fez-se a análise de aproximação da solução numérica com a solução analítica, para uma malha com 5
nós e refinando a malha para 10 nós, conforme mostrado nas Fig. 3 e 4.
Figura 3. Soluções numérica e analítica e erro associado a cada ponto para o caso da aleta com malha de 5 nós.
Anais do XV ENCITA 2009, ITA, Outubro, 19-22, 2009
,
Figura 4. Solução numérica e analítica e erro associado a cada ponto para o caso da aleta com malha de 10 nós.
Duas observações foram imediatas com a plotagem do erro do método numérico. Mesmo tendo diminuído com o
refinamento da malha, o maior erro ocorreu em um volume de controle de uma extremidade, como já havia ocorrido em
exemplos estudados anteriormente. A região mais próxima da condição de contorno de Dirichlet apresentou maiores
erros que perto da condição de contorno de Newman. Certamente, a primeira condição de contorno aproxima uma
derivada de segunda ordem.
4. Análise da condução bidimensional do calor
A análise de condução bidimensional do calor é particularmente útil em engenharia para geometrias com condução
de calor desprezível ao longo de sua espessura. Apenas dois casos estudados foram mostrados. Caso (I): Problema de
Dirichlet no retângulo com temperaturas constantes nas fronteiras, e Caso (II): Problema envolvendo condições de
contorno de fluxo e temperatura imposta.
Caso (I): A Eq. (22) é a simplificação da Eq. (1).
      
k T  k T  S  0
x  x  y  y 
(22)
A convenção adotada para as análises feitas nos casos de condução bidimensional do calor é mostrada na Fig.5.
Figura 5. Convenção para volume de controle bidimensional.
Integrando em um volume de controle e por diferenças finitas, segundo a convenção da Fig.5, tem-se a Eq. (23).
Anais do XV ENCITA 2009, ITA, Outubro, 19-22, 2009
,
 T   T    T   T  
   kA
   qV  0 
 kA
   kA
    kA
 x e  x  w   y n  y  s 
 2k y 2k x

k y
k y
k x
k x

 q1xy  TP 
TW 
TE 
TN 
TS  q2 xy

y
x
x
y
y
 x

(23)
O cálculo das temperaturas foi feito percorrendo linhas da malha, de maneira iterativa. Para cada coluna da malha
resolveu-se um sistema pelo método TDMA. Através da Eq. (23), tem-se as equações para os pontos internos. Os
volumes de controle nas fronteiras foram analisados separadamente, observando um padrão de alteração da equação
principal Eq. (23) de acordo com as condições de Dirichlet. Para ilustração, as equações discretizadas para a primeira
coluna é mostrada no sistema da Eq. (24) e para a segunda coluna, na Eq. (25), adotando a seguinte notação dos
coeficientes.
 any ,1

 aS
 ...

 0

 any ,2

  aS
 ...

 0
 aN
aP
...
...
 aN
...
aP1
...
...
...
...
 aS
...
...
...
 aS
0   Tny ,1 
 Tny ,2 




...  Tny 1,1

  a  Tny 1,2   q 
E
2 1
 ... 
aN   ... 
 



a1,1   T1,1 
 T1,2 
0  Tny ,2 
 Tny ,3 




...  Tny 1,2 
Tny 1,3 
 aE 
a
 ...  W
aN  ... 




a1,2 
 T1,2 
 T1,3 
 Tny ,1 


 Tny 1,1   q 
2 2
 ... 


 T1,1 
(24)
(25)
...
0 
 any ,1 aN


aS aP1
...
... 

Onde: 
é a matriz de coeficientes das temperaturas calculadas em linha; e q2 1 é o vetor
...
...
... aN 


 0
... aS a1,1 

com os termos fonte para cada equação.
A resolução do problema foi simulada em MATLAB e comparada com resultado de simulações em Fluent de
acordo com os dados: dimensões da geometria retangular: 0,3 x 0,4m; dimensões do volume de controle: 0,01 x 0,01m;
condutividade térmica: k=1000 W/m.K. As temperaturas dos contornos estão mostradas na formulação matemática da
Eq. (26), que é a simplificação da Eq. (1).
 2T  0
T ( x, 0)  400
T ( x , 0, 4)  100
T (0, y )  300
T (0,3 , y )  200
0  x  0,3 0  y  0, 4
0  x  0,3
0  x  0,3
0  y  0, 4
0  y  0, 4
(26)
Anais do XV ENCITA 2009, ITA, Outubro, 19-22, 2009
,
Figura 6. (a) Gráfico plotado pelo programa construído em MATLAB. (b) Simulação usando o programa Fluent.
Para a solução analítica do problema, fez-se a decomposição em 4 problemas com condições homogêneas.
Utilizando o método de separação de variáveis, resolveu-se separadamente os problemas de Sturm-Liouville regular,
encontrando as autofunções e os respectivos autovalores. Com as soluções fundamentais e o princípio de superposição
generalizado, encontrou-se a solução formal, após a última etapa de determinação do coeficiente da série de Fourier por
uma extensão periódica ímpar. O campo escalar de temperatura encontrado analiticamente é mostrado na Eq. (27).
T ( x, y ) 
1
 (2n  1) (0.4  y)   (2n  1) x 
sinh 
 sin 

 n 1
0.3
0.3
 (2n  1) 0.4 

 

(2n  1)sinh 

0.3


480


(27)
Caso (II): Foi analisado o mesmo problema anterior, apenas substituindo três condições de temperatura imposta
por condições de fluxo, de acordo com a Eq. (28).
 2T  0
T ( x, 0.4)  100
0  x  0.3 0  y  0.4
0  x  0.3
Ty ( x, 0)  0
0  x  0.3
Tx (0, y )  500000
0  y  0.4
Tx (0.3, y )  200
0  y  0.4
(28)
Figura 7. (a) Gráfico plotado pelo programa construído em MATLAB. (b) Gráfico plotado pelo Fluent.
5. Análise da condução tridimensional do calor
Da mesma forma como foi feito no caso bidimensional, o cálculo da matriz tridimensional de temperaturas é feito
em planos separadamente, mas de maneira iterativa, de modo que planos adjacentes interferem progressivamente na
distribuição de temperatura no sólido.
Figura 7. Convenção para volume de controle tridimensional.
Anais do XV ENCITA 2009, ITA, Outubro, 19-22, 2009
,
As equações discretizadas são semelhantes ao caso bidimensional, pois foram construídas pelo mesmo método. A
equação para um ponto interno é mostrada na Eq. (29), sendo os coeficientes definidos na Tab. 2. Ou seja, o termo
independente do sistema gerado, análogo aos mostrados nas Eqs. (24) e (25), tem a contribuição das temperaturas da
terceira dimensão.
aPTP  aW TW  aETE  aS TS  aN TN  aBTB  aT TT  q2
Tabela 1. Coeficientes da Eq. (29).
aE
aW
aS
k y z
k y z
k xz
x
x
y
aN
k xz
y
aB
k xy
z
aT
k xy
z
(29)
aP
aE +aW +aS +aN +aB +aT-q1ΔV
Para os pontos em fronteira, os coeficientes são alterados (somar os termos mostrados nas Tabs. 3 e 4) de maneira a
admitir as condições de contorno do problema nas faces. Os volumes de controle em uma aresta ou em uma quina do
paralelepípedo devem ter os respectivos coeficientes alterados simultaneamente pelas faces as quais faz fronteira.
Tabela 2. Alteração dos coeficientes da Eq. (29) caso haja temperatura imposta como condição de contorno.
Face
aP
q2
S
2k xz
2k xz
S
y
y
N
2k xz
2k xz
N
y
y
W
2k y z
2k y z
W
x
x
E
2k y z
2k y z
E
x
x
B
2k xy
2k xy
B
z
z
T
2k xy
2k xy
T
z
z
Tabela 3. Alteração dos coeficientes da Eq. (29) caso haja fluxo de calor como condição de contorno.
Face
q2
S
FluxoS.x.z
N
FluxoN.x.z
W
FluxoW .y.z
E
FluxoE.y.z
B
T
FluxoB.x.y
FluxoT .x.y
Importante explicitar a ordem de iteração da rotina (Anexo II): os planos paralelos xy da malha são percorridos no
sentido B -> T (valores crescentes de z), e nos planos, as linhas são percorridas no sentido W -> E (valores crescentes
de x) e nestas linhas, as equações dos pontos são resolvidas no sentido S -> N (valores crescentes de y). O
conhecimento da ordem de iteração é importante para o operador do código atribuir as condições de contorno de modo
a obter uma convergência mais rápida.
6. Análise da condução de calor em um caso real
Como forma de exemplificar a análise da condução tridimensional de calor, o modelo foi aplicado utilizando parte
do projeto de cooperação-científica TKMCL-ITA intitulado “Aquecimento de matriz de forjamento por radiação na
banda infravermelha provinda de lâmpadas”. O projeto simula o comportamento térmico durante o pré-aquecimento de
matrizes de forjamento. Foram realizados ensaios com diferentes potências das lâmpadas medindo-se a temperatura da
Anais do XV ENCITA 2009, ITA, Outubro, 19-22, 2009
,
superfície da matriz em função do tempo, conforme mostrado no gráfico da Fig. 9. O objetivo do referido trabalho é
estimar a potência mais adequada das lâmpadas para aquecimento das matrizes em menor tempo.
Figura 8. Ensaio para medida de temperaturas com termopares de junção na superfície do bloco de aço.
Temperatura na Superfície da Matriz (oC)
360
340
320
300
280
260
240
220
200
180
160
140
120
100
Superf-Matriz Potência 60%
Superf-Matriz Potência 70%
Superf-Matriz Potência 80%
Superf-Matriz Potência 90%
80
60
40
20
0
0
2
4
6
8
10
12
14
16
18
20
22
24
26
28
30
Tempo (min)
Figura 9. Temperaturas na superfície da matriz para diferentes níveis de potência em função do tempo medidas
com os termopares de junção mostrados na Fig. 8.
Para a simulação usando o programa Fluent as propriedades termofísicas estão mostradas na Tab. 4 e as dimensões
da geometria são x=500mm; y=44mm; z =160mm. O valor de convecção natural foi h = 15 W/m2°C e a temperatura
ambiente 27°C. As simulações estão mostradas nas Figs. 8 e 9.
Tabela 4. Variação das Propriedades Termofísicas do material com a temperatura
Propriedade
Temperatura (°C)
Valor adotado para a
propriedade
20
400
Densidade (kg/m3)
7800
7700
7750
Condutividade Térmica – k (W/m°C)
25
29
27
Calor Específico (J/kg.K)
442
444
443
Figura 8. (a) Campo de temperatura (°C) na matriz para 1 minuto com potência de 60%. (b) Campo de temperatura
(°C) na matriz para 11 minutos com potência de 60%.
Anais do XV ENCITA 2009, ITA, Outubro, 19-22, 2009
,
Figura 9. (a) Campo de temperatura (°C) na matriz para 1 minuto com potência de 90%. (b) Campo de temperatura
(°C) na matriz para 5 minutos com potência de 90%.
Como a rotina construída não abrange o regime transiente, calculou-se a distribuição da matriz de temperaturas
para o bloco de aço considerando o valor de temperatura na superfície da matriz em 1min como sendo constante, caso
semelhante ao mostrado nas Figs. 8(a) e 9 (a). Apenas como ilustração da aplicação do modelo criado, tomou-se como
condições de contorno a temperatura de 70oC (dada pelo gráfico da Fig.9) na superfície superior do bloco e a
temperatura ambiente (medida em 27oC) nas outras faces restantes do paralelepípedo. Esses dados foram rodados na
rotina criada, tendo como saída o gráfico da Fig. 10. Observação: a maior temperatura de contorno foi adotada na face
inferior, já que o método iterativo se inicia por convergência nesta face.
Figura 10. Difusão de temperatura no bloco de aço (k=27 W/m°C) calculada pela rotina.
7. Conclusões
O trabalho realizado com o modelo desenvolvido por discretização apresentou bons resultados, principalmente para
casos uni e bidimensionais, por se utilizar de uma linguagem robusta de programação, a qual permitiu o procedimento
recorrente de operações matriciais. Para as aplicações tridimensionais, pôde-se perceber um distanciamento um pouco
pronunciado do resultado esperado, sendo necessário o refinamento da malha utilizada e levando a um tempo
computacional expressivo, já que a malha gerada foi uniforme em todos os casos. Os resultados do modelo para
exemplos 1D apresentaram erros da ordem de 2% para malhas grosseiras, em relação à solução analítica. Para aplicação
bidimensional, a discordância dos valores foi da ordem de 0,1% em comparação com malhas idênticas simuladas em
Fluent. Os casos tridimensionais apresentaram distorções geométricas mais acentuadas, sendo necessário um
refinamento, que não foi realizado. Ainda, o método de expansão da resolução dos problemas de 2D para 3D necessitou
muitas alterações de coeficientes, tornando extenso o código. Por fim, todo o projeto foi embasado na construção de
malhas retangulares uniformes, assim não houve comparações entre iterações feitas entre malhas diferentes.
8. Sugestões de prosseguimento do trabalho
Seguindo as conclusões apresentadas, otimizações do código podem ser feitas de maneira a tornar menos lento o
procedimento computacional. Atalhos podem ser criados na identificação de casos simples e relações de simetria, já que
generalidade do programa implicou a análise de comparações desnecessárias e criação de variáveis ociosas.
Anais do XV ENCITA 2009, ITA, Outubro, 19-22, 2009
,
Implementações podem ser feitas ainda fazendo-se as atualizações dos coeficientes dos nós de contorno de maneira
matricial, diminuindo, assim, o numero de linhas do código.
9. Agradecimentos
Agradecimento, principalmente, ao CNPq e à CAPES pela concessão das bolsas e também ao Professor Sandro.
10.
Referências
Nascimento, J. J. da S.; Belo, F. A.; de Lima, A. G. B., (2002). “Solução numérica da equação de difusão
transiente aplicada a sólidos paralelepípedos”. II Congresso Nacional de Engenharia Mecânica.
Veersteg, H. K., & Malalasekera, W. (1995). "An introduction to computational fluid dynamics", The finite volume
method. New York: John Wiley & Sons Inc. p.96
Incropera, F. P., & DeWitt, D. P. (2002). "Fundamentals of heat and mass transfer". Fifth Edition. New York: John
Wiley & Sons Inc.
Tikhonov, A. N., & Samarskii, A. A.(1963). "Equations of mathematical physics". New York: Dover Publications
Inc.
Miller, Kenneth S. (1953)."Partial differential equations in engineering problems". New York: Prentice-Hall Inc.
Trim, Donald W.(1989). "Applied partial differential equations". Boston: PWS-Kent Publishing Company.
Kakaç, S. and Yener, Y., Heat Conduction, 2 Ed., Hemisphere Publishing Corporation, Washington, 1985.
Incropera, F.P.; DeWitt, D.P. Transferência de Calor e de Massa. Apêndice A Propriedades Termofísicas da Matéria. 5ª
Ed., LTC, Rio de Janeiro, 2003, pag. 643.
Cooperação-Científica TKMCL-ITA, Doc: 01/2008, Projeto: “Aquecimento de Matriz de Forjamento por Radiação
na Banda Infravermelha provinda de Lâmpadas”.