How to convert from euler´s method to ODE45?

1 回表示 (過去 30 日間)
Cesar García
Cesar García 2015 年 12 月 20 日
i want to use ODE45 based in the next code.
% Parameters
um=0.54;
Km=0.03;
k1=0.4;
k2=1.24;
Csx=0.095;
Csp=0.15;
Cnx=0.20;
Rcsx=0.022;
Rcnx=0.26;
%Initial Conditions
X(1)=1.16;
S(1)=20.05;
N(1)=2.00;
P(1)=1.01;
V(1)=3;
Sin=60;
Nin=0;
%Specific growth rate, we need to be care here, because u is not a function of time and it looks like a DAE.
u(1)=um*((N(1)/S(1))/((N(1)/S(1))+Km));
%simulation time and step for euler tsim=25;
t(1)=0;
dt=0.00001;
i=1;
while t(i)<tsim
if t(i)>4.9 && t(i)<5.001
F1=0.8;
F2=0;
else
F1=0;
F2=0;
end
V(i+1)=V(i)+(F1+F2)*dt;
X(i+1)=X(i)+(u(i)*X(i))*dt-((F1+F2)/V(i))*X(i)*dt;
S(i+1)=S(i)-(Csx*u(i)*X(i)+Rcsx*X(i))*dt-(Csp*(k1*u(i)*X(i)+k2*X(i)))*dt+(F1/V(i))*Sin*dt-((F1+F2)/V(i))*S(i)*dt;
P(i+1)=P(i)+((k1*u(i)*X(i))+(k2*X(i)))*dt-((F1+F2)/V(i))*P(i)*dt;
N(i+1)=N(i)-((Cnx*u(i)*X(i))+(Rcnx*X(i)))*dt+(F2/V(i))*Nin*dt-((F1+F2)/V(i))*N(i)*dt;
if N(i+1)<0
N(i+1)=0;
end
u(i+1)=um*((N(i+1)/S(i+1))/((N(i+1)/S(i+1))+Km));
t(i+1)=t(i)+dt;
i=i+1;
end
I will appreciate your help.
Best Regads

回答 (0 件)

カテゴリ

Help Center および File ExchangeOrdinary Differential Equations についてさらに検索

タグ

Community Treasure Hunt

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

Start Hunting!

Translated by