Debugging ode45. Need Help Please.

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

댓글 수: 2

Torsten
Torsten 2015년 1월 29일
Where do you define y0(1) abd y0(3) ?
Best wishes
Torsten.
Torsten
Torsten 2015년 1월 29일
In both routines, you will have to define y0 as global.
Best wishes
Torsten.

댓글을 달려면 로그인하십시오.

답변 (0개)

카테고리

도움말 센터File Exchange에서 General Applications에 대해 자세히 알아보기

태그

질문:

2015년 1월 29일

댓글:

2015년 1월 29일

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by