Análise numérica envolvendo problema de transmissão de calor por condução

Prévia do material em texto

UNIVERSIDADE FEDERAL DE PERNAMBUCO - UFPE
CENTRO DE TECNOLOGIA E GEOCIÊNCIAS - CTG
DEPARTAMENTO DE ENGENHARIA QUÍMICA - DEQ
CURSO DE ENGENHARIA QUÍMICA
ANÁLISE NUMÉRICA ENVOLVENDO PROBLEMA DE TRANSMISSÃO DE CALOR POR CONDUÇÃO 
CARLOS BRIAN OLIVEIRA DE CARVALHO
ELAIANE KARINE DA SILVA BARBOSA
FABIANA MARIA GOES OLIVEIRA SANTOS
MIRELA RIBEIRO EMBIRASSÚ DE ARRUDA
OSCAR DA CUNHA FERREIRA NETO
RECIFE - 2017
INTRODUÇÃO
Ao longo do curso de transferência de calor foram apresentados modelos matemáticos e suas equações governantes que descrevem a transferência de calor por condução assim como suas soluções analíticas. É claro observar que as soluções analíticas geralmente são possíveis somente para casos relativamente simples. 
Contudo, essas soluções têm uma função importante na análise de transferência de calor, pois fornecem insight em problemas complexos de engenharia que podem ser simplificados usando certas suposições.
Entretanto, muitos problemas práticos envolvem geometrias e condições de contornos complexas, propriedades termofísicas variáveis e não podem ser resolvidos analiticamente, mas por métodos numéricos ou de computador que incluem métodos de diferença finita, elemento finito, elemento de contorno, entre outros. 
Além de fornecer um método de solução para esses problemas mais complexos, a análise numérica é mais eficiente em termos de tempo total necessário para encontrar uma solução. 
Outra vantagem é que as alterações nos parâmetros dos problemas podem ser mais facilmente realizadas, permitindo que um engenheiro determine o comportamento de um sistema térmico ou otimize um sistema térmico muito mais facilmente.
Com base nisto, será realizada uma análise numérica envolvendo transferência de calor em aleta de seção variável aplicados aos problemas 3.8 e 3.9 do livro Princípios de Transferência de Calor (Kreith et. al, 7ª edição), utilizando o método iterativo de Gauss-Seidel através de código desenvolvido no Matlab.
Método iterativo de Gauss-Seidel
O método consiste em uma pequena alteração do método de Jacobi, pois os valores atualizados das variáveis são utilizados imediatamente na mesma iteração durante o cálculo das demais variáveis. Ou seja, chuta-se um valor inicial para cada variável e calcula-se uma por uma cada interação. O critério de convergência é baseado em um erro de diferença entre as temperaturas calculadas em cada interação. 
 
METODOLOGIA
2.1 Apresentação dos problemas
3.8 – Desenvolva uma equação das diferenças para o volume de controle para condução estável unidimensional em uma aleta com área transversal variável A(x) e perímetro P(x). O coeficiente de transferência de calor do estabilizador vertical para o ambiente é uma constante 0 e a ponta do estabilizador é adiabático.
3.9 – Usando os resultados do problema 3.8, encontre o fluxo de calor na base da aleta para as seguintes condições:
2.2 Considerações do problema
Regime permanente;
Condução unidimensional na direção x;
Condutividade do material constante ao longo da aleta;
Coeficiente de convecção (h) é constante;
Sem geração interna de calor.
2.3 Dedução das equações nodais
Considerando um volume de controle mostrado na Figura 1, foi realizado um balanço de energia. 
Figura 1 – Volume de controle do problema
 =+ 
Multiplicando por , obtém-se:
Isolando :
 (Equação 1)
Para a resolução do problema, é necessário definir as condições de contorno. Sabe-se que na base da aleta a temperatura é especificada, logo T1=T0, onde T1 é a temperatura do nó 1 e T0 é a temperatura da parede. Para o último nó foi preciso realizar um balanço de energia a fim de obter a segunda condição de contorno. A Figura 2 descreve o processo realizado.
Figura 2 – Balanço de energia no último nó
Como a extremidade da aleta é adiabática não há condução na face direita do volume de controle. 
 = 
	Isolando TN na equação, tem-se:
	Conhecendo as condições de contorno do primeiro e do último nós, aplica-se a equação 1 para determinação das temperaturas dos nós internos, as quais serão apresentadas a seguir:
Onde e 
Em seguida, os termos C foram isolados como apresentados abaixo:
	A partir das equações obtidas foi possível implementar um código numérico em Matlab para resolução do problema proposto.
	
CÓDIGO EM MATLAB
% Equação base: 
% T(i) = [C(i) - T(i-1)*Area(i) - T(i+1)*Area(i+1)]/denominador(i)
% Sendo C(i) e denominador(i)constantes de auxílio 
% C(i) = -dx^2*P(i)*h0*Tinf/k
% denominador(i) = (-Area(i) - Area(i+1) - dx^2*P(i)*h0/k)
% Método de Gauss-Seidel
% Matriz estritamente diagonal dominante para convergência |aii| > soma|aij|
 
%DADOS--------------------------------------------
k = 34; %W/m K
L = 0.05; %m
dx = 0.005; %m
x = 0:dx:L; %posição em cada nó
Area = 0.047*(1 - 1/3*sin(x/L)); %m^2 Area em cada posição do nó
P = 0.3*Area.^0.5; %m Perímetro
h0 = 110; %W/m^2 K
T0 = 93; %[C] - temperatura fixa no nó 1
Tinf = 27; %[C] - temperatura ambiente
n = L/dx + 1; %nós
%-----------------------------------------------
%criando o termo independente C(i)
for i = 2:n-1 último e primeiro nó já estão predefinidos
 C(i) = -dx^2*P(i)*h0*Tinf/k; %termo independente do numerador
end
 
A = zeros(n,n); %cria matriz de zeros e define tamanho 11x11
 
%Cria os termos do denominador e monta a matriz de temperaturas
for i = 2:n-1
 denominador(i) = (-Area(i) - Area(i+1) - dx^2*P(i)*h0/k);
 for j = 2:n-1
 if i==j
 A(i,i) = denominador(i);
 A(i,i-1) = Area(i);
 A(i,i+1) = Area(i+1);
 end
 end
end
A(1,1) = 1; %coeficiente 1 já está definido T1 = T0;
A(n,n) = 1; %já está definido T11 = T10;
 
b = zeros(n,1); %cria zeros e define 1 coluna com 11 linhas
b(2:n-1) = C(2:n-1); % Coloca os termos independentes em b
b(1) = T0; 
b(n) = C(n-1); %vetor de coeficientes
 
Ab = [A b]; %Junta as duas matrizes --> equações e coeficientes
As = Ab; % Se, caso a matriz não fosse dominante, alterava ordem aqui
 
T(2:n) = 80; % Condições iniciais, chute
T(1) = T0;
 
erro = zeros(n,1); %define matriz de zeros iniciais
 
for iter = 1:105 % definição de iteração i=1:iteração escolhida
 for k = 2:n-1
 Tantigo = T(k);
 numerador(k) = As(k,end) - Area(k) * T(k-1) - Area(k+1)*T(k+1);
 T(k) = numerador(k)/denominador(k);
 erro(k)= abs(T(k) - Tantigo);
 iteracao_T(iter,k) = T(k); %matriz iterativa das temperaturas
 iteracao_erro(iter,k) = erro(k); %matriz iterativa dos erros
 end
 disp(['Iteracao ',num2str(iter),'; Erro = ', num2str(max(erro))]);
end
%-------------------------------------------------
T(n) = (h0*P(n)*(1/2*dx^2)*Tinf - T(n-1)*(h0*P(n)*(1/2*dx^2)-k*Area(n)))/(k*Area(n)); %T11 = T10 -->definição para ponta adiabática
transpose(T) % mostra a temperatura no workspace
%Resultados---------------------------------------
fluxo_base = -k*(T(1)-T(2))/(x(1)-x(2)) %fluxo de calor no nó 1, x = 0
Thor = (T); % T horizontal para mostrar no imagesc
figure(1)
imagesc(T)
colorbar
figure(2)
plot(T)
FLUXOGRAMA PARA SOLUÇÃO ITERATIVA
RESULTADOS E DISCUSSÃO
A resolução do sistema linear de equações visto anteriormente foi feita baseada na construção de uma matriz tridiagonal 11x11 dos coeficientes, como mostrada na figura 3.
 
Figura 3 – Matriz tridiagonal elaborada pelo código
A partir dos cálculos iterativos feitos com auxílio do código desenvolvido em Matlab foi possível obter a distribuição de temperatura ao longo da aleta de seção variável. 
Observou-se que a temperatura nos nós diminui ao longo do comprimento da aleta como já era esperado, pois a mesma perde calor para o ambiente. Como a temperatura da base (T1) é maior quea do meio, ao longo do processo de condução ocorre o resfriamento da aleta, como pode ser visto na escala de cores da Figura 4.
Figura 4 – Escala de temperatura ao longo do comprimento da aleta
A Tabela 1 apresenta os valores das temperaturas nodais utilizados posteriormente para plotar o gráfico mostrado na Figura 5.
	Nó
	Temperatura (°C)
	1
	93.0000
	2
	91.8582
	3
	90.6838
	4
	89.4749
	5
	88.2302
	6
	86.9488
	7
	85.6300
	8
	84.2743
	9
	82.8826
	10
	81.4569
	11
	81.4446
Tabela 1 – Temperaturas nodais
Figura 5 – Representação gráfica das temperaturas em função dos nós
O gráfico acima também mostra a variação das temperaturas nodais por todo comprimento da aleta. Duas observações foram feitas a partir de sua leitura: entre os nós 1 e 10 vê-se um declínio linear de temperatura devido ao fato da perda de calor por convecção que se dá através da redução da seção reta da aleta. 
Uma segunda observação feita foi que entre os nós 10 e 11 a variação de temperatura é muito pequena – da ordem de , sendo assim praticamente estável. Esta baixa variação é causada pela extremidade da aleta ser adiabática.
O critério de parada do programa foi feito a partir do número de iterações a fim de que o erro fosse o mínimo possível – na ordem de . A figura 6 apresenta os erros obtidos em função do número de iterações. A partir da quarta iteração observa-se uma queda drástica no erro. 
Figura 6 – Representação gráfica do erro
Por fim foi calculado o fluxo de calor na base da aleta como pedido no problema. O fluxo foi calculado através da seguinte expressão:
A partir do programa foi obtido o valor de .
CONCLUSÃO
Conclui-se, portanto, que a ferramenta utilizada na resolução do problema facilita a obtenção de um erro mínimo. A partir do programa, pode-se refinar o tamanho da malha a fim de se ter mais dados de temperaturas. 
A análise numérica de transferência de calor auxilia a resolver problemas onde não seria possível a resolução de maneira analítica.