Clear Filters
Clear Filters

For Loop Not Working??

114 views (last 30 days)
Kevin
Kevin on 9 Jul 2014
Edited: dpb on 20 Aug 2023
Hi,
I am designing a programme which will operate using a for loop starting at a value of 4 and ending at 10 as shown in the code below:
for TSR=4:10;
I am unsure as to whether I should put the following line of code in the end of my program:
TSR=TSR+i;
I have also secified that i is a counting index of 1. If I dont put it in, the loop seems to keep interating for a TSR value of 4. But if I do put it in the for loop surpasses the final value of 10 and increases continuously. Could anybody please help with this problem?
Thanks, Kevin

Answers (7)

dpb
dpb on 9 Jul 2014
Your code snippet is incomplete; you have no i defined as shown what more as "a counting index of 1" .
The for index variable in your case is TSR, there's no i involved and, "Yes, Virginia, you do not artificially index the for loop variable inside a for construct."
doc for % and all will be explained...read all the doc including the _Tips_ section at the bottom.
Then try the following exercise at the command line --
for i=4:10
disp(i)
i=i+1;
end
and explain the results.
NB: while often used as "traditional" loop indices, in Matlab i and j are predefined as the imaginary sqrt(-1). You can override them as above and it is often done even in Matlab, but be aware of what you are doing when do so particularly if there's any chance whatsoever that there are real (so to speak :) ) complex values lurking...

dpb
dpb on 9 Jul 2014
...I put an end if statement at the end which seems to have solved it...
Solved what???
for TSR=4:10 % TSR from 4 to 10
while abs(Check)>=tolerance
...
if TSR==10
else
TSR=TSR+i;
end
end
end
The solution is as documented to NOT try to update the loop variable inside the loop at all.
Remove the code snippet
if TSR==10
else
TSR=TSR+i;
end
entirely; it isn't needed and doesn't do anything, anyway. As the demo I suggested you run shows, the value of the loop index (TSR here in your case; i in mine) is overwritten anyway by the next loop iteration. Since after you updated TSR it isn't ever referenced until the beginning of the next loop iteration, that local copy that is changed is never referenced. To see this even more clearly, modify the sample I gave before to--
for i=4:10,disp(i),i=i+2;disp(i),end
and observe the results.
As the above conclusively demonstrates, whatever the actual problem you had before, it was NOT owing to the cause you've surmised; there's insufficient data provided to tell what, precisely, you did incorrectly before as you didn't post an example code that reproduced the problem to see...
I'll make some additional comments on the code in making it more "Mablaby"...
Rewrite
r_local=[(R./N)*1; (R./N)*2; ... (R./N)*8; (R./N)*9];
as
r_local=bsxfun(R./N),[1:9]); r_local=r_local(:);
Use inbuilt constant pi instead of crude approximation for pi...multiple places.
C=((8*pi.*r_local)./(B.*Cl_design)).*(1-cos(Phi));
for TSR=4:10
while abs(Check)>=tolerance
axial_induction_old=axial_induction;
TSR_local=TSR*(r_local./R);
Phi=(2/3)*atan(1./TSR_local);
axial_induction=1./(((4.*(sin(Phi).^2))./(sigma.*Cl_design.*cos(Phi)))+1);
angular_induction=(1-(3*axial_induction))./((4.*axial_induction)-1);
relative_wind=(1-axial_induction)./((1+axial_induction).*TSR);
F=2/pi*acos(exp(-((B/2.*(1-r_over_R))./(r_over_R.*sin(relative_wind)))));
C_T=(sigma.*((1-axial_induction).^2).*((Cl.*cos(relative_wind))+ ...
(Cd.*sin(relative_wind))))./((sin(relative_wind)).^2);
if C_T<0.96
axial_induction=1./(1+(4.*F.*(sin(relative_wind).^2))./ ...
(sigma.*Cl.*cos(relative_wind)));
else
axial_induction=1./(((4.*F.*cos(relative_wind))./(sigma.*Cl))-1);
end
D=(8/(TSR.*N)).*(F.*(sin(Phi).^2).*(cos(Phi)- ...
((TSR_local).*(sin(Phi)))).*(sin(Phi)+ ...
((TSR_local).*(cos(Phi)))).*(1-(Cd./Cl).*atan(Phi)).*(TSR_local.^2));
Cp=sum(D);
Diff=axial_induction-axial_induction_old;
Check=max(Diff(:));
store_Phi(:,i)=Phi;
store_TSR_local(:,i)=TSR_local;
store_axial_induction(:,i)=axial_induction;
store_angular_induction(:,i)=angular_induction;
store_relative_wind(:,i)=relative_wind;
store_Check(:,i)=Check;
store_Diff(:,i)=Diff;
store_Cp(:,i)=Cp;
store_TSR(:,i)=TSR;
end
end

Ben11
Ben11 on 9 Jul 2014
Edited: Ben11 on 9 Jul 2014
That is very strange; the default increment for for-loops is 1 so you should not have to add the TSR = TSR+i. By the way is it really +i or +1 ? Maybe you have something in your code which makes TSR smaller than 10 at every iteration so the loop keeps going. Can you post your whole code so we can better help?
  2 Comments
dpb
dpb on 9 Jul 2014
Edited: dpb on 9 Jul 2014
... Maybe you have something in your code which makes TSR smaller than 10 at every iteration so the loop keeps going.
That's not how Matlab loops work...the demo I posted for Kevin illustrates this just as well with a slight modification...try
for i=4:10,disp(i),i=i-2;disp(i),end
at the command line.
Again I commend
doc for
to your attention--read the "Tips" section carefully...
Clearly if his code had an infinite loop he had "something" in his code that caused it but whatever it was, it was not the for loop index itself. As noted earlier, it's not possible to tell what his error was that caused the symptom; he didn't post any code that demonstrated the claimed behavior.
Ben11
Ben11 on 9 Jul 2014
Yep you're so right I need more coffee I guess thanks for correcting my non sense!

Sign in to comment.


Kevin
Kevin on 9 Jul 2014
Thanks for the help. I put an end if statement at the end which seems to have solved it. This is my entire code:
clc
clear
filename = 'book2.xlsx';
Input = xlsread(filename)
N = xlsread(filename,'C7:C15')
r = xlsread(filename,'D7:D15')
C = xlsread(filename,'E7:E15')
S = xlsread(filename,'F7:F15')
Non = xlsread(filename,'G7:G15')
% Inputs
R=0.4; % Radius of Rotor
B=3; % Number of blades
V=2; % Fluid velocity
Rho=998; % Fluid Density
N=9; % Number of Blade Elements
Cp_estimate=0.5; % Estimate power coefficient
Alpha_design=4; % Design alpha
Cl_design=1.04; % Design lift coefficient
% Variables
TSR=4;
Cp=0;
i=1; % Counter
alpha_new=0;
axial_induction=0;
tolerance=0.1; % Tolerance Value
Check=1;
axial_induction_old=0; % Initial value for old axial induction factor
r_local=[(R./N)*1; (R./N)*2; (R./N)*3; (R./N)*4; (R./N)*5; (R./N)*6; (R./N)*7; (R./N)*8; (R./N)*9];
r_over_R=r_local./R;
TSR_local=r_over_R.*TSR;
Phi=(2/3)*atan(1./TSR_local);
C=((8.*3.14.*r_local)./(B.*Cl_design)).*(1-cos(Phi));
sigma=(B*C)./(3.14.*r_local.*2);
Cl=[1.3; 1.1; 1; 0.9; 0.86; 0.83; 0.8; 0.75; 0.5]; % Lift Coefficients
Cd=[0.027; 0.024; 0.02; 0.019; 0.018; 0.016; 0.013; 0.012; 0.01]; % Drag Coefficients
F=(2/3.14).*acos(exp(-(((B/2).*(1-(r_over_R)))./((r_over_R).*sin(Phi)))));
for TSR=4:10; % TSR from 4 to 10
while abs(Check)>=tolerance
axial_induction_old=axial_induction;
TSR_local=TSR*(r_local./R); % Local Tip Speed Ratio
Phi=(2/3)*atan(1./TSR_local); % Angle of Relative Fluid
axial_induction=1./(((4.*(sin(Phi).^2))./(sigma.*Cl_design.*cos(Phi)))+1);
angular_induction=(1-(3*axial_induction))./((4.*axial_induction)-1);
relative_wind=(1-axial_induction)./((1+axial_induction).*TSR);
F=(2/3.14).*acos(exp(-(((B/2).*(1-(r_over_R)))./((r_over_R).*sin(relative_wind))))); % Tip Loss Factor
C_T=(sigma.*((1-axial_induction).^2).*((Cl.*cos(relative_wind))+(Cd.*sin(relative_wind))))./((sin(relative_wind)).^2);
if C_T<0.96
axial_induction=1./(1+(4.*F.*(sin(relative_wind).^2))./(sigma.*Cl.*cos(relative_wind)));
else
axial_induction=1./(((4.*F.*cos(relative_wind))./(sigma.*Cl))-1);
end
D=(8./(TSR.*N)).*(F.*(sin(Phi).^2).*(cos(Phi)-((TSR_local).*(sin(Phi)))).*(sin(Phi)+((TSR_local).*(cos(Phi)))).*(1-(Cd./Cl).*atan(Phi)).*(TSR_local.^2));
Cp=sum(D);
Diff=axial_induction-axial_induction_old;
Check=max(Diff(:));
store_Phi(:,i)=Phi;
store_TSR_local(:,i)=TSR_local;
store_axial_induction(:,i)=axial_induction;
store_angular_induction(:,i)=angular_induction;
store_relative_wind(:,i)=relative_wind;
store_Check(:,i)=Check;
store_Diff(:,i)=Diff;
store_Cp(:,i)=Cp;
store_TSR(:,i)=TSR;
if TSR==10
else TSR=TSR+i;
end
end
end
figure(1)
plot(store_Cp,store_TSR)
hold all
title('Cp vs Tip Speed Ratio')
xlabel('TSR')
ylabel('Cp')

dpb
dpb on 10 Jul 2014
I seriously recast your code below...I may have ended up w/ a mismatched parens or two as I tried to remove many of the superfluous pairs of them and the unneeded dot operators to make it simpler to try to parse the actual code.
The biggest change, however, was to preset a vector of TSR values over which to iterate and then replaced the loop index to select from that vector. Then, the results are stored in a (preallocated) set of arrays by column from the first to last instead of off in the yonder somewhere as before.
Then, it appears you reversed the role of x- and y- axes in your plot command so I swapped them around.
filename = 'book2.xlsx';
data = xlsread(filename); % only read the file once...
N = data(7:15,3); % assumes numeric array is at origin;
r = data(7:15,4); % adjust row/col offset to the actual
C = data(7:15,5); % locations based on what text is also in
S = data(7:15,6); % the spreadsheet
Non = data(7:15,7);
%
% Inputs
R=0.4; % Radius of Rotor
B=3; % Number of blades
V=2; % Fluid velocity
Rho=998; % Fluid Density
N=9; % Number of Blade Elements
Cp_estimate=0.5; % Estimate power coefficient
Alpha_design=4; % Design alpha
Cl_design=1.04; % Design lift coefficient
%
% Variables
TSR=4; % Initial tip speed ratio
Cp=0; % Initial power coefficient
i=1; % Counter
alpha_new=0; % Initial value for alpha new
axial_induction=0; % Initial axial induction factor
tolerance=0.01; % Tolerance Value
Check=1; % Initial check value
axial_induction_old=0; % Initial value for old axial induction factor
r_local=bsxfun(R./N),[1:9]); r_local=r_local(:);
r_over_R=r_local./R;
TSR_local=r_over_R.*TSR;
Phi=(2/3)*atan(1./TSR_local);
C=(8*pi*r_local./(B*Cl_design)).*(1-cos(Phi));
sigma=B*C./(2*pi*r_local);
Cl=[1.3; 1.1; 1; 0.9; 0.86; 0.83; 0.8; 0.75; 0.5]; % Lift Coefficients
Cd=[0.027; 0.024; 0.02; 0.019; 0.018; 0.016; 0.013; 0.012; 0.01];
F=2/pi*acos(exp(-(B/2*(1-r_over_R)./(r_over_R.*sin(Phi)))));
axial_induction=1/(4*sin(Phi).^2./(sigma.*Cl_design.*cos(Phi))+1);
angular_induction=(1-3*axial_induction)./(4*axial_induction-1);
TSRvec=4:10; % set the desired TSR values
[nr nc]=[length(r_local) length(TSRvec)]; % result sizes for preallocation
store_Phi=zeros(nr, nc);
store_TSR_local=zeros(nr, nc);
store_axial_induction=zeros(nr, nc);
store_angular_induction=zeros(nr, nc);
store_relative_wind=zeros(nr, nc);
store_Check=zeros(nr, nc);
store_Diff=zeros(nr, nc);
store_Cp=zeros(nr, nc);
store_TSR=zeros(nr, nc);
for i=1:length(TSRvec) % and iterate over them
TSR=TSRvec(i); % set local temporary from vector
while abs(Check)>=tolerance
axial_induction_old=axial_induction;
TSR_local=TSR*(r_local./R);
Phi=(2/3)*atan(1./TSR_local);
relative_wind=(1-axial_induction)./((1+axial_induction)*TSR(i));
F=2/pi*acos(exp(-(B/2*(1-r_over_R)./(r_over_R.*sin(relative_wind)))));
C_T=(sigma.*((1-axial_induction).^2).*(Cl.*cos(relative_wind)+ ...
Cd.*sin(relative_wind)))./((sin(relative_wind)).^2);
if C_T(TSR,:)<0.96
axial_induction=1./(1+(4.*F.*(sin(relative_wind).^2))./ ...
(sigma.*Cl.*cos(relative_wind)));
else
axial_induction=1./(((4.*F.*cos(relative_wind))./(sigma.*Cl))-1);
end
D=(8/(TSR*N))*(F.*sin(Phi).^2.*(cos(Phi)-(TSR_local.*sin(Phi))).* ...
(sin(Phi)+(TSR_local.*cos(Phi))).*(1-Cd./Cl.*atan(Phi)).*TSR_local.^2);
Cp=sum(D);
Diff=axial_induction-axial_induction_old;
Check=max(Diff(:));
store_Phi(:,i)=Phi;
store_TSR_local(:,i)=TSR_local;
store_axial_induction(:,i)=axial_induction;
store_angular_induction(:,i)=angular_induction;
store_relative_wind(:,i)=relative_wind;
store_Check(:,i)=Check;
store_Diff(:,i)=Diff;
store_Cp(:,i)=Cp;
store_TSR(:,i)=TSR;
end
end
%
figure
plot(store_TSR,store_Cp)
title('Cp vs Tip Speed Ratio')
xlabel('TSR')
ylabel('Cp')
There is one serious problem in your code that undoubtedly is not doing what you intend--since I don't know what the intent is for certain I can only point it out for your attention--
if C_T(TSR,:)<0.96
axial_induction=1./(1+(4.*F.*(sin(relative_wind).^2))./ ...
(sigma.*Cl.*cos(relative_wind)));
else
axial_induction=1./(((4.*F.*cos(relative_wind))./(sigma.*Cl))-1);
end
In Matlab, the test above will be TRUE iff all values in the vector C_T are true so the above branch will only trigger when that occurs. If it is need for the case to be taken for any one value, then you would need
if any(C_T(TSR,:))<0.96
As noted, you'll need to fix this logic up as you intend.
  1 Comment
Amin Mohammed
Amin Mohammed on 22 May 2018
Well explained, I'm facing similar issue, but still not able to work it out. I'll try to take the advantage of this answer. In answer provided by @dpb, in line 56:
relative_wind=(1-axial_induction)./((1+axial_induction)*TSR(i));
Should read:
relative_wind=(1-axial_induction)./((1+axial_induction)*TSR);

Sign in to comment.


Kevin
Kevin on 9 Jul 2014
Hi dpb,
Thank you so much for the help. I'm still relatively new to the world of Matlab and your advice really makes my life a bit easier.
Best regards, Kevin
  1 Comment
dpb
dpb on 10 Jul 2014
This part of it would be easier yet if you would take and apply some of the advice given... :)

Sign in to comment.


Maximilian
Maximilian on 19 Aug 2023
For the for loop below it works when I do not include the (e^(-n/r) part and I replace
it with (e*(-n/r). How can I make it work with the raised to the power part.
e = exp(1);
N= 1001;
Ic= zeros(1,N);
Vs = 1;
for n = 1:N
Vc(n) = 1 - (e^(-n/r));
Ic(n) = (Vs - Vc(n))/R;
end
  7 Comments
Walter Roberson
Walter Roberson on 19 Aug 2023
I tend to forget about eps(0) so what I actually used was log(eps(realmin)) which is identical to log(eps(0))
dpb
dpb on 19 Aug 2023
Chuckle...

Sign in to comment.

Categories

Find more on Loops and Conditional Statements 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!