An iteration Code causing error
1 次查看(过去 30 天)
显示 更早的评论
Thanks guys for the response so far
function []=model3()
%parameters
h=0.0006944444; V=10;
tfinal=30;
Qe = 10; Qo = Qe; Y=0.4;
Xo=0; Xe=Xo; u1=0.1; S=500; Snoe1=100; Ks=10; Kno=0.5; b=0.02;
%initial conditions
X1(i)=200; Sse(i)=5; Snoe(i)=2;
t(i)=0;
%define the ODE function hundle
f1=@(t,X1) ((Qo*Xo-Qe*Xe-V*X1/10)/V)+u1*(Snoe/(Ks+Snoe))*(Snoe/(Kno+Snoe))*X1-b*X1;
f2=@(t,Snoe) (Qo*Snoe1-Qe*Snoe)/V-(1-Y)/(2.86*Y)*u1*(Sse/(Ks+Sse))*(Snoe/(Kno+Snoe))*X1;
f3=@(t,Sse) (Qo*S-Qe*Sse)/V-u1/Y*(Sse/(Ks+Sse))*(Snoe/(Kno+Snoe))*X1;
for i=1:ceil(tfinal/h)
%update t
t(i+1)=t(i)+h;
%update x
kX11=f1(t(i), X1(i));
kX12=f1(t(i)+0.5*h, X1(i)+0.5*kX11*h);
kX13=f1(t(i)+0.5*h, X1(i)+0.5*kX12*h);
kX14=f1(t(i)+ h, X1(i)+ kX13*h);
X1(i+1)=X1(i)+(h/6)*(kX11+2*kX12+2*kX13+kX14);
%update Snoe
kn1=f3(t(i), Snoe(i));
kn2=f3(t(i)+0.5*h, Snoe(i)+0.5*kn1*h);
kn3=f3(t(i)+0.5*h, Snoe(i)+0.5*kn2*h);
kn4=f3(t(i)+ h, Snoe(i)+ kn3*h);
Snoe(i+1)=Snoe(i)+(h/6)*(kn1+2*kn2+2*kn3+kn4);
%update Sse
ks1=f3(t(i), Sse(i));
ks2=f3(t(i)+0.5*h, Sse(i)+0.5*ks1*h);
ks3=f3(t(i)+0.5*h, Sse(i)+0.5*ks2*h);
ks4=f3(t(i)+ h, Sse(i)+ ks3*h);
Sse(i+1)=Sse(i)+(h/6)*(ks1+2*ks2+2*ks3+ks4);
end
end
1 个评论
Jan
2018-12-3
You forgot to post the error message. It is much easier to fix a problem than to guess, what the problem is.
回答(2 个)
Jan
2018-12-3
The code stops in:
X1(i)=200; Sse(i)=5; Snoe(i)=2;
t(i)=0;
While i is not defined yet, it is the imaginary unit. I guess you meant 1 instead of i.
X1(1)=200; Sse(1)=5; Snoe(1)=2;
t(1)=0;
0 个评论
YT
2018-12-3
编辑:YT
2018-12-3
You say you're getting an error, but you didnt specify it in the description of your question. I ran the script and I assume you meant the following error
Unable to perform assignment because the left and right sides have a different number of elements.
which first occured at
Snoe(i+1)=Snoe(i)+(h/6)*(kn1+2*kn2+2*kn3+kn4);
In this case Snoe will be a 1-by-X array, but you are trying to add an vector to it Snoe(i)+(h/6)*(kn1+2*kn2+2*kn3+kn4).
I think the easiest way to fix it, is by doing
Snoe(i+1,:) %also for X1/Sse
Extra:
I would initialize these matrices before the loop so
X1 = zeros([ceil(tfinal/h)+1 V]);
Snoe = zeros([ceil(tfinal/h)+1 V]);
Sse = zeros([ceil(tfinal/h)+1 V]);
for i=1:ceil(tfinal/h)
...
end
1 个评论
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!