An iteration Code causing error

2 views (last 30 days)
WILLIAM BAYA MWARO
WILLIAM BAYA MWARO on 3 Dec 2018
Edited: Jan on 3 Dec 2018
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 Comment
Jan
Jan on 3 Dec 2018
You forgot to post the error message. It is much easier to fix a problem than to guess, what the problem is.

Sign in to comment.

Answers (2)

Jan
Jan on 3 Dec 2018
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;

YT
YT on 3 Dec 2018
Edited: YT on 3 Dec 2018
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 Comment
Jan
Jan on 3 Dec 2018
Edited: Jan on 3 Dec 2018
After I fixed the 'i' to '1' to define the inital values, the code ran fine. The pre-allocation is an important advice.

Sign in to comment.

Categories

Find more on Programming 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!