How to convert from euler´s method to ODE45?

1 次查看（过去 30 天）

显示更早的评论

Cesar García 2015-12-20

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/261110-how-to-convert-from-euler-s-method-to-ode45

在 MATLAB Online 中打开

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;