Debugging ode45. Need Help Please.
Ältere Kommentare anzeigen
Here's my code. I'm not sure why it's not running:
function CSTR tspan=[0 10.]; y0=[1.0; 1.0; 100; 80]; % initial concentrations of A, B,Temp Reactor, Temp Cooling Jacket
% A=1.0 M
% B=1.0 M
[t,y]=ode45(@CSTRfunction,tspan,y0);
plot(t,y, 'LineWidth',3);
title('CSTR')
xlabel('Time (mins)'),ylabel('Variables A,B,T, TCJ')
legend('[A]','[B]','T'),grid
end
function dydt=CSTRfunction(time,y)
V= 100;
VR= 100;
E1= -9758.3;
E2= -9758.3;
E3= 8560;
k10= 1.287*10^12;
k20= 1.287*10^12;
k30= 9.043*10^9;
HAB= 4.2;
HBC= -11;
HAD= -41.85;
rho= .9342;
Cp=3.01;
Kw= 4032;
A= .215;
mk= 5;
CpK= 2;
Qk=-1113.5;
k1=k10*exp(E1/(y(3)+273.15));
k2=k20*exp(E2/(y(3)+273.15));
k3=k30*exp(E3/(y(3)+273.15));
dCAdt=V/VR*(y0(1)-y(1))-k1*y(3)*y(1)-k3*y(3)*(y(1))^2;
dCBdt=-V/VR*y(2)+k1*y(3)*y(1)-k2*y(3)*y(2);
dTdt=V/VR*(y0(3)-y(3))-1/rho/Cp*(k1*y(3)*y(1)*HAB+k2*y(3)*y(2)*HBC+k3*y(3)*(y(1))^2*HAD)+Kw*A/rho/Cp/VR*(y(4)-y(3));
dTCJdt=1/mk/CpK*(Qk+Kw*A*(y(3)-y(4)));
% dydt
dydt=[dCAdt,dCBdt,dTRdt,dTRKdt]
end
Antworten (0)
Kategorien
Mehr zu General Applications finden Sie in Hilfe-Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!