EXEMPLOS RESOLVIDOS DO LIVRO AMOS GILAT

Prévia do material em texto

UNIVERSIDADE DO ESTADO DO AMAPÁ
ENGENHARIA QUIMICA
 MODELAGEM E SIMULAÇÃO DE PROCESSOS	
					
	
MODELAGEM E SIMULAÇÃO DE PROCESSOS
					Acadêmico: Muller dos Santos Silva
Macapá-Ap
2017
EXEMPLOS RESOLVIDOS DO LIVRO AMOS GILAT
Exemplo 8-1
Script
function [x, y] = edoEULER(EDO,a,b,h,yINI)
N = (b - a)/h;
for i = l:N
x(i + 1) = x(i) + h;
y(i + 1) = y(i) + feval(EDO(x(i) ,y(i)))*h;
end
clear
a = 0 ; b = 2.5 ; h = 0.1 ; yINI = 3 ;
xp = a:0.5:b;
yp=70/9*exp(-0.3*xp)-43/9*exp(-1.2*xp);
plot(x,y,'--r' ,xp,yp)
xlabel ('x') ; ylabel ('y')
function dydx = Cap8Exmp1EDO(x,y)
dydx = -1.2*y + 7*exp(-0.3*x);
end
Comentários:
Desculpa-me por não desenvolver mais a questão pelo não entendimento do script.
No método de Euler Explicito nesta questão mostra e neste método depende do valor de h, porque quanto menor o valor de h, menor será o erro, comparando outros métodos mostrado no livro sobre aproximações o método de Euler explicito mostra-se como de segunda ordem.
Exemplo 8-2
Script
%Solução de uma EDO de primeira usando o método de implícito de Euler
clear all
a = 0; b = 0.5; h = 0.002;
N = (b - a)/h;
n(l) = 2000; t(l) = a;
for i=l:N
t(i + 1) = t(i) + h;
x = n(i);
% Tem início o método de Newton.
for j = 1:20
 num =x + 0.800*x (3/2) *h - 10.0*n (1) * (1 - exp (-3*t (i + 1))) *h - n (i) ;
denom = 1+0.800*1.5*x (1/2) *h;
xnew = x - num/ denom;
if abs((xnew - x)/x) < 0.0001
break
else
end
end
x = xnew;
end
end
x = xnew;
if j == 20
fprintf('Não foi possivel calcular a solução númerica em t=%gs', t(i))
break
end
% Termina o metodo de Newton.
n(i + 1) = xnew;
end
plot(t,n)
axis([0 0.5 0 2000]), xlabel('t (s) ') , ylabel('n')
Comentários:
Neste exemplo usa-se o método de Euler implícito pois depende dos dois lados da equação f(t,n) não for o suficiente, podemos resolver a equação isolando os termos, mas como vemos no problema que não é possível e temos que recorrer a métodos vistos nas aulas passadas para determinar raízes da equação, como o método de Newton que está no script e o método da bisseção. 
Exemplo 8-8
Script
function [t, x, y]=Sis2EDOsRK4(EDO1,EDO2,a,b,h,xl,yl)
t(l) = a; x(l) = xl; y(l) = yl;
n=(b - a)/h;
for i = l:n
 t(i+l) = t(i) + h;
 tm = t(i) + h/2;
 Kxl = EDOl(t(i) ,x(i) ,y(i));
 Kyl = EDO2(t(i) ,x(i),y(i));
 Kx2 = EDOl(tm,x(i)+ Kxl*h/2,y(i)+ Kyl*h/2);
 Ky2 = EDO2(tm,x(i)+ Kxl*h/2,y(i)+ Kyl*h/2);
 Kx3 = EDOl(tm,x(i)+ Kx 2*h/2,y(i)+ Ky2*h/2);
 Ky3 = EDO2(tm,x(i)+ Kx2*h/2,y(i)+ Ky2*h/2);
 Kx4 = EDOl(t(i + 1),x(i)+ Kx3*h,y(i)+ Ky3*h);
 Ky4 = EDO2(t(i + 1),x(i)+ Kx3*h,y(i)+ Ky3*h);
 x(i+l) = x(i) + (Kxl + 2*Kx2 + 2*Kx3 + Kx4)*h/6;
 y(i+l) = y(i) + (Kyl + 2*Ky2 + 2*Ky3 + Ky4)*h/6;
end
function dxdt = PenduloDtetaDt(t,x,y)
dxdt = y;
function dydt = PenduloDwDt(t,x,y)
c = 0.16; m = 0.5; g = 9.81; L = 1.2;
dydt = -(c/m)*y - (g/L)*sin(x);
[t, x, y] = Sis2EDOsRK4('PenduloDtetaDt,PenduloDwDt',0,18,0.1,pi/2,0);
plot(t,x)
xlabel('Tempo(s)')
ylabel ('ânngulo teta (rad)')
Comentários:
Neste exemplo usa-se o método de Ruge-Kutta de segunda, terceira e quarta ordem, em que superiores a quarta ordem se mostram menos eficazes, sem modificações.
Diferentes dos métodos anteriores consistem em aproximações locais em função de h, calcula-se os Ks primeiramente, a diferença desse método a sua o acréscimo dessas equações para aproximações e tornar os valores os mais precisos possíveis.
	
Exemplo 8-10
Script
function TempChapa(T1,Vol,Area,Tamb)
tInt = [ 0 180 ];
[Time Temp] = edo45 (@TaxaTemp,tInt,Tl,[ ],Vol,Area,Tamb);
plot(Tempo,Temp)
xlabel ('Tempo(s)')
ylabel('Temperatura(K)')
function dTdt=TaxaTemp(t, T, Vol, As, Tamb)
Rho= 300; Cv = 900; h = 30; Epsi=0.8; Sigma= 5.67e-8;
ARVC =As/(Rho*Vol*Cv) ;
SigEps =Epsi*Sigma;
dTdt = -ARVC*(SigEps*(T^4 - Tamb^4) + h*(T - Tamb));
Comentários:
Neste exemplo é um pouco diferente dos outros por conter comandos ou argumentos opcionais e também usando a função residente ode45 em problemas não-rígidos, é a melhor função a ser usada como primeira tentativa. Método de passo simples baseado em métodos de Runge-Kutta explícitos de quarta e quinta ordem.
E como a questão trata de transferência de calor na chapa haverá uma variação de volume em certo tempo.
Exemplo 9-1
Script
clear all
a=0; b=0.l; TINI=473; wINil=-1000; h=0.001;
[x, Tl, w]=Sis2EDOsRK2('edoCap9ExmpldTdx','edoCap9Exmpldwdx',a,b,h,TINI,wINI1);
n = length(x);
fprintf('A temperatura em x=0.1 é %5.3f, para o valor inicial dt/dx= %4.lf\n' ,Tl(n) ,wINil)
wINI2 = -3500;
[x, T2, w] =Sis2EDOsRK2('edoCap9ExmpldTdx','edoCap9Exmpldwdx',a,b,h,TINI,wINI2);
fprintf ('A temperatura em x=0.1 é %5.3f, para o valor inicial dt/dx= %4. lf\n' ,T2 (n) ,wINI2)
wINI3 = wINil + (293 - Tl (n)) * (wINI2 - wINil) I (T2 (n) - Tl (n));
[x, T3, w] =Sis2EDOsRK2('edoCap9ExmpldTdx','edoCap9Exmpldwdx',a,b,h,TINI,wINI3);
fprintf ('A temperatura em x=0.1 é %5.3f, para o valor inicial dt/dx = %4. lf\n' ,T3 (n) ,wINI3)
plot (x, Tl, '-k' , x, T2, '-k' , x, T3, '-r')
xlabel ('Distância (m)') ; ylabel ('Temperatura(K)')
function dTdx = edoCap9ExmpldTdx(x,T,w)
dTdx = w;
function dwdx = edoCap9Exmpldwdx'(x,T,w)
he = 40; P = 0.016; eps = 0.4; k = 240; Ac = 1.6E-5; Seg = 5.67E-8;
Ts = 293;
kAc = k*Ac;
Al = hc*P/kAc; A2 = eps*Seg*P/kAc;
dwdx = Al*(T - Ts) + A2* (T^4 - Ts^4) ;
Comentários:
Neste exemplo usa-se o método do tiro que é o método do tiro é um método numérico utilizado para a resolução de PVCs. Nesse método, um PVC é convertido em dois PVIs pela transformação de uma EDO de segunda ordem em doas EDOs de primeira ordem, cada uma com a sua condição inicial. Uma EDO de segunda ordem é convertida em duas EDOs de primeira ordem pela mudança de variável. Igualmente a questão anterior resolve-se pelo método de Runge-kutta usando a função Sis2EDOsRK2 e como dita nas questões anteriores um dos métodos para chegar um valor possivelmente aproximado da condição fornecida.
Exemplo 9-2
Script
clear all
a = 0; b = 0.1; TINI = 473; h = 0.001; Yb = 293;
tol = 0.01; imax = 15;
wH = -1000;
[x, T, w] = Sis2EDOsRK2('edoCap9ExmpldTdx','edoCap9Exmpldwdx',a,b,h,TINI,wH);
n = length(x);
wL = -3500;
[x,T,w]=Sys20DEsRK2('edoCap9ExmpldTdx','odeCap9Exmpldwdx',a,b,h,TINI,wL);
for i = 1: imax + 1
 wi = (wH + wL)/2;
 [x,T,w]=Sis2EDOsRK2('edoCap9ExmpldTdx','edoCap9Exmpldwdx',a,b,h,TINI,wi);
 E = T (n) - Yb;
 if abs(E) < tol
break
end
if E > 0
wH = wi;
else
 wL = wi
end
end
if i > imax
fprintf('A solução não foi obtida em %i iterações.\n' ,imax}
else
plot(x,T}
xlabel('Distância(m} '); ylabel('Temperatura(K}'}
fprintf('A temperatura calculada x = 0.1 é %5.3f K.\n' ,T(n))
fprintf('A solução foi obtida em %2.0f iterações.\n' ,i}
end
Comentários:
Igualmente a questão 9-1 usando o método do tiro e agora recorrendo ao método de Newton e o método da bisseção e pode-se se observar bem a o uso dos métodos verificando os gráficos das temperaturas, onde as aproximações com as 9 iterações.