EXERCICIOS RESOLVIDOS DE CONTROLE E SIMULAÇÃO DE PROJETOS

Prévia do material em texto

ALUNA: VANESSA WANDERSEE CUNHA OSTROSKI
MÉTODO DE EULER
EXEMPLO 1 (Euler)
clc;
clearall;
%PARÂMETROS
E = 0.9; %emissividade
cte_rad = 5.6e-8; %constante de Boltzmann W/m2 K4
D = 0.7e-3; %diâmetro da esfera m
Cp = 400; %capacidade calorífica J/Kg K
ro = 8500.0; %densidade kg/m3
h1 = 400.0; %coeficiente convectivo W/m2 K
Tp = 400.0 + 273.15; %temperatura da parede K
T_ar = 200.0 + 273.15; %temperatura do ar K
To = 25.0 + 273.15; %temperatura inicial K
Nt = 500; %pontos no tempo
tf = 6.0; %tempo final s
dt = tf / Nt;
lambida1 = (6*h1)/(D*ro*Cp);
lambida2 = (6*E*cte_rad)/(D*ro*Cp);
%condição inicial
T1(1) = To ;
T11(1) = To - 273.15;
%loop do tempo
for n = 1:Nt;
T1(n+1) = T1(n) + dt * (lambida1 * (T_ar - T1(n)) + (lambida2 * (Tp^4 - T1(n)^4)));
T11 (n+1) = T1(n+1) - 273.15;
end
t = (0.0:1.0:Nt) * dt;
plot (t,T11,'g')
holdon
%Resposta 1: T = 218,18 °C no regime estacionário.
%Resposta 2: Não é possível determinar o tempo necessário, pois como têm-se um regime estacionário a temperatura inicial mantêm-se constante.
EXERCÍCIO 5.8 (Euler)
clc;
clear all;
%PARÂMETROS E CONSTANTES
E = 0.8; %emissividade
cte_rad = 5.6e-8; %constante de Boltzmann W/m2 K4
Cp = 900; %capacidade calorífica a pressão constante J/Kg K
p = 2800.0; %densidade kg/m3
q = 1.25e4; %fluxo térmico de calor gerado W/m2
h = 10.0; %coeficiente convectivo W/m2 K
A = 0.040; %área m2
L = 0.007; %espessura da base do ferro m
Tp = 25.0 + 273.15; %temperatura da parede (vizinhança) K
T_ar = 25.0 + 273.15; %temperatura do ar K
To = 25.0 + 273.15; %temperatura inicial K
Nt = 1000; %número de pontos no tempo
tf = 6.0; %tempo final de simulação (s)
dt = tf / Nt;
lambida1 = (A*h)/(ro*Cp*(A*L))
lambida2 = (A*E*cte_rad)/(ro*Cp*(A*L))
lambida3 = (q*A)/(ro*Cp*(A*L))
%condição inicial
T(1) = To;
T11(1) = To - 273.15;
%loop do tempo
for n = 1:Nt;
T(n+1) = T(n) + dt * (-((lambida1 * ( T(n) - T_ar)) + (lambida2 * (T(n)^4 - Tp^4)))) + lambida3;
T11 (n+1) = T(n+1) - 273.15;
end
t = (0.0:1.0:Nt) * dt;
plot (t,T11,'g')
holdon
%Resposta: t = 4,68 s, a placa atingiu a T = 135,2991 °C.
 
EXERCÍCIO 5.15 (Euler)
clc;
clear all;
%PARÂMETROS
Cp = 2200.0; %capacidade calorífica J/Kg K
ro = 1200.0; %densidade kg/m3
hi = 10000.0; %coeficiente convectivo interno W/m2 K
he = 2000.0; %coeficiente convectivo externo W/m2 K
V = 2.25; %volume m3
Ds = 0.02; %diâmetro da serpentina m
R = Ds /2;
Ls = 11; %comprimento da serpentina m
To = 300.0; %temperatura do ar K
Ts = 500.0; %temperatura da serpentina K
Nt = 50000; %pontos no tempo
tf = 3600; %tempo final de simulação (valor estimado por mim)
dt = tf / Nt;
U = 1 / ((1/hi)+(1/he));
lambida1 = (U * pi * Ls * Ds)/(V*ro*Cp);
tempo = 24883 * dt / 60
%condição inicial
T1(1) = To;
%loop do tempo
for n = 1:Nt;
T1(n+1) = T1(n) + dt * (lambida1 * ( Ts - T1(n)));
end
t = (0.0:1.0:Nt) * dt; %para transformar ponto discreto em intervalo de tempo
plot (t,T1,'g')
grid on %linhas de grade
xlabel('Tempo (s)') %escrever o título do eixo
ylabel ('Temperatura (K)') %escrever o título do eixo
title ('Temperatura x Tempo') %escrever o título da figura
%Resposta1: t = 30 min, a placa aqueceu de 300 a 400 K.
%Resposta2: t = 60 min, têm-se um comprimento de serpentina menor.
MÉTODO DAS DIFERENÇAS FINITAS
CONDUÇÃO BARRA (MDF)
clc;
clearall;
%parametros
dx= 0.01;
dt = 0.01;
Nt = 100; % números de pontos no tempo
Nx = 20; % número de pontos na posição
T0 = 20;
T1 = 100;
T2 = 50;
Fo = 0.5; %constante cinética 1/s
% vetor para o plot da posição ao longo da barra
x = ( 0.0:1.0:Nx -1 );
% vetor para o plot do tempo
t = ( 0.0:1.0:Nt -1);
% condição inicial primeiro passo no tempo
fori=1: Nx-1
T(1,i)=T0;
end
% plotinitialvalues
figure(1)
plot(T,'b*-') % linha azul com pontos
xlabel ('direção x')
ylabel('temperatura')
title ('Placa plana condição inicial')
gridoff
holdoff
% loop do tempo
for n = 1:Nt;
T(n,1)= T1;
T(n, Nx) = T2;
T(Nt +1, 1) = T1; 
T(Nt +1, Nx) = T2 
for i= 2:Nx -1 
T(n+1,i)= Fo*(T(n,i+1) + T (n,i-1)) + T(n,i)*(1-2*Fo); % E.M.M
end
end
%plotperfildinâmico de temperatura
figure (2)
plot (x,T)
xlabel('deslocamento')
ylabel('deslocamento')
title ('temperatura')
gridon
holdoff
% plot de superficie
n = (0.0:1.0:Nt);
figure (3)
surf(x,n,T)
colormaphsv
xlabel('posição');
ylabel('tempo');
zlabel('temperatura');
title('superficie');
Outra maneira:
clc;
clear all;
 
T1 = 100.0;
T2 = 50.0;
T3 = 30.0;
L = 1.0; %metros
t_final = 100.0;
Fo = 0.5; %valor fornecido pelo problema (se alterarmos para 1.0 ou 0.55, o gráfico diverge. Isso ocorre porque na equação geral temos (1-2Fo), logo o valor final terá que ser positivo, para isso (1-2Fo>0 ou 1-2Fo=0).
Nt = 100; %
Nx = 20; %
dt = t_final/Nt;
dx = L / Nx;
 
%condição inicial
for i = 1: Nx+1;
 T(1,i) = T3;
end
for n = 1: Nt;
 T(n,1) = T1;
 T(n,Nx+1) = T2;
 T(Nt+1,1) = T1;
 T(Nt+1,Nx+1) = T2;
 
 for i = 2: Nx;
 T(n+1,i) = Fo*(T(n,i+1)+T(n,i-1))+(1-2*Fo)*T(n,i);
 end
end
 
x = (0.0:Nx)*dx;
plot(x,T)
%EXERCÍCIO 5.105 
clc;
clear all;
%PARÂMETROS
ro = 7800.0; %densidade kg/m3
Cp = 700.0; %calor específico J/Kg K
k = 30; %condutividade térmica W/m K
To = 1400.0; %temperatura inicial °C
T_ar = 50.0; %temperatura do ar °C
h = 5000.0; %coeficiente convectivo W/m2 K
L = 0.1; %espessura da placa m
Nt = 2000;
Nx = 10;
tf = 180.0;
dt = tf / Nt;
dx = L / Nx;
lambida1 = (h*dx)/k;
lambida2 = (k / (ro*Cp))*(dt/dx^2);
%condição inicial
for i=1:Nx+1
T(1,i) = To;
end
%loop do tempo
for n = 1:Nt;
T(n+1,1) = lambida2*(T(n,2) + T(n,1)) + (1-2*lambida2)*T(n,1) ;
aa =lambida1*(T_ar - T(n,Nx+1)) + T(n,Nx+1);
T(n+1,Nx+1) = lambida2*(aa + T(n,Nx)) + (1-2*lambida2)*T(n,Nx+1);
for i=2:Nx
T(n+1,i) = lambida2*(T(n,i+1) + T(n,i-1)) + (1-2*lambida2)*T(n,i) ;
end
end
tempo = 1900*dt
l = tempo *15e-3
t = (0.0:1.0:Nt)*dt;
plot (t,T)
grid on
%Resposta1: t = 171 s, para a placa ser resfriada até 200 °C.
%Resposta2: L = 2,56 m, para a placa se mover a 15 mm/s.
%TROCADOR DE CALOR CASCO E TUBO CONTRACORRENTE
%Turbulento casco e tubo,Ti(t,x) e Te(t,x): convecção e troca com a parede
%agua quente passa nos cascos
clc;
clear all;
L = 2.802; %comprimento sugerido do trocador
tf = 50.0;
n_passes = 2; % número de passes
n_tubos = 148; 
D = 0.0127; %diametro sugerido do tubo interno
Dext = 0.38;%diâmetro sugerido do tubo externo
R = D/2;
Rext = Dext/2;%raio sugerido do tubo interno
U = 2800; % coeficiente global de troca térmica
p_q = 979.13 ;% densidade quente de motor a 60ºC
p_f = 1015.23 ;% densidade frio kg/m3 á 25ºC
cp_q = 4188.9; %calor específico fluido quente 90ºC
cp_f = 4031.9; %calor específico água á 25ºC
Nx = 10; %número de pontos no espaço
Nt = 270000;
Tqe = 90.0; %temperatura do fluido quente na entrada do trocador
Tfe = 32.0; %temperatura da fluido frio na entrada do trocador
vm_q = 0.0217; %vazão volumétrica Q m3/s 
vm_f = 0.02778/(n_tubos); %vazão volumétrica Q m3/s 
v_f = (4*vm_f)/(pi*D^2) % velocidade do frio, m/s
v_q = (vm_q)/(pi*(Rext^2)-(n_passes*n_tubos*R^2)) % velocidade da quente, m/s
TEMPO_f = L /(v_f); 
TEMPO_q = L /v_q;
areaext=pi*((Rext^2)-(n_passes*n_tubos*R^2));
areaint=pi*R^2;
 
%calculo das constantes
dx = L/Nx;
dt = tf / Nt;
cte2 = (4*U*dt)/ (D *p_q *cp_q); %constante de quente 
cte1 = (2*R*U*n_passes*n_tubos*dt)/ (((Rext^2)-(n_passes*n_tubos*R^2)) *p_f *cp_f); %constante de frio
cte3 = (v_q *dt)/ (dx); %constante de óleo 
cte4 = -(v_f *dt)/ (dx); %constante de água
 
% condição inicial (primeiro passo no tempo)
for i = 1:Nx+1;
Tq(1,i) = Tqe;
Tf(1,i) = Tfe;
end
 
%loop no tempo
for n = 1:Nt; 
 Tq(n,1) = Tqe;
 Tq(Nt+1,1) = Tqe;
 
 Tf(n,Nx+1) = Tfe;% cc1 
 Tf(Nt+1,Nx+1) = Tfe;Tf(n+1,1) = -cte4*(Tf(n,2)-Tf(n,1))+(cte2*(Tq(n,1)-Tf(n,1)))+Tf(n,1);% M.M. T água 
 
 for i = 2:Nx;
 Tq(n+1,i) = -cte3*(Tq(n,i)-Tq(n,i-1))+(cte1*(Tf(n,i)-Tq(n,i)))+Tq(n,i);% M.M. T óleo 
 Tf(n+1,i) = -cte4*(Tf(n,i+1)-Tf(n,i))+(cte2*(Tq(n,i)-Tf(n,i)))+Tf(n,i);% M.M. T água 
 
 end
 Tq(n+1,Nx+1) =-cte3*(Tq(n,Nx+1)-Tq(n,Nx))+(cte1*(Tf(n,Nx+1)-Tq(n,Nx+1)))+Tq(n,Nx+1);% M.M. T óleo 
end
x = (0.0:1.0:Nx)*dx;
 
for i = 1:Nx+1;
 Tq1(i) = Tq(10,i);
 Tf1(i) = Tf(10,i);
 Tq2(i) = Tq(Nt+1,i);
 Tf2(i) = Tf(Nt+1,i);
end
%plot perfil dinâmico da concentração Ca1
figure(1)
plot (x,Tq1,'r*-')
hold on
plot (x,Tf1,'b*-')
hold on
plot (x,Tq2,'m')
hold on
plot (x,Tf2,'g')
xlabel('direção axial, [m]')
ylabel('temperatura Tq e Tf, [ºC]')
title('Perfil espacial da temperatura de fluido quente e frio - CONTRACORRENTE')
grid on
legend('Tquente, t=10','Tfrio,t=10','Tquente, t=Nt+1','Tfrio,t=Nt+1')
%figure(2)
%plot (x,Tq)
%figure(3)
%plot (x,Ta)