MATLAB Answers

Why if condition does not store my variables?

1 view (last 30 days)
Meva
Meva on 3 Dec 2016
Commented: Meva on 4 Dec 2016
Hello,
Please forgive me my post is long. But I debugged and I came up with no solution. I want to store my u, p and gap arrays at time t=2, t=3, t=4, and t=5, then plot all. To do this, I put if commands associated with the times above as can be seen from the code below. At first, Matlab gave me the plots. However,now matlab says undefined variable u2 which is the first stored variable.
clear all;close all;clc;
Tsim=4;
ii=101;
dt=10^(-5);
tstart = 1;
dx=0.01;
Nt=(Tsim/dt)+1;
% Nt = 1;
ii_vect = 1:ii;
x = (ii_vect-1)*dx;
%----------------------------------------------------------------
b = 0.5; %|
kappa = 1; %|
gbar = .1; %|
M = 1; %|
g = gbar/M; %|
Gamma = 0.2; %|
F = kappa.*x.*(1-x); %|
h = kappa*b*(1-b); % h0 = F(b); %|
theta = kappa*(1-2*b); % htheta0 = F'(b); %|
%----------------------------------------------------------------
pi=4.*atan(1);
hD=0;
thetaD=0;
gap0 = 0.5;
cb=gap0;
ub=zeros(1,ii); u=zeros(1,ii); p=zeros(1,ii); q=zeros(1,ii); gap=zeros(1,ii);
gap=gap0+(-F+h+theta.*(x-b));
ub=(cb-hD.*(x-b)-thetaD/2.*(x-b).^2)./gap;
dc=-10^(-5);
%______________ TIME ITERATION____________%
for nt=1:Nt
flag=0;
mmm=1;
m=1;
c=cb;
t=tstart+(nt-1)*dt;
while (flag==0)
u=(c-hD.*(x-b)-thetaD/2.*(x-b).^2)./gap;
p(1)=0.5*(1-u(1)^2);
q(1)=0;
for i=2:ii
q(i)=q(i-1)-dx*(u(i-1)-ub(i-1))/dt;
p(i)=0.5*(1-u(i)^2)+q(i);
end
st(m)=p(ii)-0;
m=m+1;
if (m==2)
c=c+dc;
end
if (m==3)
c=(c*st(1)-(c-dc)*st(2))/(st(1)-st(2));
end
if (m==4)
mmm=mmm+1;
if (mmm==2)
m=1;
else
flag=1;
end
end
end
sumint1=zeros(1,ii); sumint2=zeros(1,ii);
sumint1=0.5*(p(1)); sumint2=0.5*(p(1))*(-b);
for k=2:ii-1
xx=(k-1)*dx;
sumint1=sumint1+(p(k));
sumint2=sumint2+(p(k))*(xx-b);
end
hDDOT=(sumint1*dx - M*g)/M;
hD=hD+dt*hDDOT;
h=h+dt*hD;
thetaDDOT=sumint2*dx/(Gamma*M);
thetaD=thetaD+dt*thetaDDOT;
theta=theta+dt*thetaD;
hL=h-theta*b;
hR=h+theta*(1-b);
gap=gap0+(-F+h+theta.*(x-b));
%____THESE IF CONDITIONS ARE NOT RECOGNISED !___%
if t==2 %%(nt==(2-tstart)/dt +1)
u2=u;
p2=p;
gap2=gap;
end
if t==3 %%(nt==(3-tstart)/dt +1)
u3=u;
p3=p;
gap3=gap;
end
if t==4 %%(nt==(4-tstart)/dt +1)
u4=u;
p4=p;
gap4=gap;
end
if t==5 %%(nt==(5-tstart)/dt +1)
u5=u;
p5=p;
gap5=gap;
end
end
ub(:)=u(:);
cb=c;
%_______PLOTTING_______%
figure(1)
subplot(2,2,1)
plot(x,u2,'k-')
hold on
plot(x,u3,'k:')
hold on
plot(x,u4,'k--')
hold on
plot(x,u5,'k-.')
hold on
legend('t=2', 't=3', 't=4','t=5')
subplot(2,2,2)
plot(x,p2,'k-')
hold on
plot(x,p3,'k:')
hold on
plot(x,p4,'k--')
hold on
legend('t=2', 't=3', 't=4','t=5')
subplot(2,2,[3,4])
plot(x,gap2,'k-')
hold on
plot(x,gap3,'k:')
hold on
plot(x,gap4,'k--')
hold on
plot(x,gap5,'k-.')
legend('t=2', 't=3', 't=4','t=5')

  0 Comments

Sign in to comment.

Accepted Answer

the cyclist
the cyclist on 3 Dec 2016
I have not gone through your code in detail, but I have an educated guess.
When you do your tests such as
if t==2
t is a floating-point variable which it looks to me might not be exactly equal to 2 (or whatever), because you are calculating it via a for loop that will not lead to exact integers. [It might be equal to 2.000000000001, for example.]
I expect you need to either calculate t in a way that ensure it is exactly equal to the integer you expect, or make your comparison in a way the includes a tolerance.

  7 Comments

Show 4 older comments
the cyclist
the cyclist on 3 Dec 2016
Yes, 1.e-8 represents 10^(-8).
My best guess regarding your point about the gaps is that if you set the tolerance too large, then there could be many points which will trigger the condition, some of which are too far away from the values of t you want. Only the last one gets saved, not the one that is closest to the actual value of t.
Meva
Meva on 3 Dec 2016
Very many thanks @cyclist ! So, would you please mind giving me a recommendation for time iteration then?
t
is the floating point now. So how can I set t iteration and how can I declare t so that it would not be a floating point ? Or is there any alternative way to do this:
Tsim=4;
dt=10^(-5);
tstart = 1;
Nt=(Tsim/dt)+1;
for nt=1:Nt
t=tstart+(nt-1)*dt
...
end
Thanks so much. It will solve my problem and anyone having trouble with floating points.
Meva
Meva on 4 Dec 2016
Hello again @cyclist,
Having look at this on the web, the best solution I have come across is
Nt=400001;
instead of the Nt line in above. Just that. I was wondering that if I have any better way?

Sign in to comment.

More Answers (0)

Sign in to answer this question.


Translated by