Problems with oder45, for or if
Show older comments
I have problem with my program. Seemingly everything is working well but when i run ode45 I get wrong results. My program doesn't see change of Ginf. I really don't know why. Ginf in two first iterations is very high and i should see this on graph. but it doesn't see
function F=model_magda(t,y)
F=zeros(2,1);
global Vg Vi Rinf Cinf Gg gamma lambda alfa C1 a1 Rm Ki Kg Tinf UFR Gd Vglucose Vinsulin Ginf
for t=1:61
if t-1<Tinf
Ginf(t)=Rinf*Cinf;
else
Ginf(t)=0;
end;
Vglucose(t)=Vg-UFR*(t-1);
Vinsulin(t)=Vi-UFR*(t-1);
F(1)=(Gg-gamma*y(1)*y(2)-lambda*y(1)+Ginf(t)-Kg*(y(1)-Gd)+y(1)*UFR)/Vglucose(t);
F(2)=(Rm/(1+exp((C1-y(1))/a1))-alfa*y(2)-Ki*y(2))/Vinsulin(t);
t=t+1;
end
-------------------------------------------------
global Vg Vi Rinf Cinf Gg gamma lambda alfa C1 a1 Rm Ki Kg Tinf UFR m Gd Vglucose Vinsulin yo meani meang time to Ginf
UFR=3/240; %0.004166667; %[L/min]
m=78; %[kg]
Gd=5.615292712;
Rinf=60/1000 %[L/min]
Vinf=0.5*m*100/(33*1000); %[L]
Cinf=33/180.16*1000*1000/100; %[mmol/L]
Tinf=Vinf/Rinf %[min]
Minf=Tinf*Rinf*Cinf;
Minfg=Minf*180.16/1000; %Ginfa=Rinf*Cinf
%----------srednia---------------
Go=6.63;
Io=8.89;
Gg=0.621129599; %[mmol/min]
gamma=0.002993997; %{L^2/mU*min]
lambda=0.101455909; %[L/min]
alfa=0.095577674; %[L/min]
C1=25.51490467; %[mmol/L]
a1=1.194742573; %[mmol/L]
Rm=79.44266365; %[mU/min]
Vg=8.764594059; %[L]
Vi=6.673505545; %[L]
time=[ 0 10 20 30 40 60];
meang=[ 6.63 18.67 14.19 11.59 10.9 9.08];
meani=[ 8.89 25.07 17.54 13.05 10.99 8.76];
Ki=63/1000; %[L/min]
Kg=223.2/1000; %[L/min]
[to yo]=ode45('model_magda',[0:60],[Go Io]);
[to yo(:,1) yo(:,2)]
figure(1)
plot(to,yo(:,1),'b')
hold on
plot(time,meang,'b--x')
grid on;
xlabel('Time[min]');
ylabel('glucose concetration in blood [mmol/l]');
legend('simulation','dane');
figure(2)
plot(to,yo(:,2),'r')
hold on
plot(time,meani,'r--x')
grid on;
xlabel('Time [min]');
ylabel('insulin concetration in blood [mU/l]');
legend('simulation','dane');
figure(3)
plot(yo(:,1),yo(:,2),'k')
hold on
plot(meang,meani,'k--x')
legend('simulation','dane');
Answers (0)
Categories
Find more on Ordinary Differential Equations in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!