While Loop not converging??

Question

Kevin 2014-7-22

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/142888-while-loop-not-converging

评论： Geoff Hayes 2014-7-22

Hi,

I am using a while loop in order to perform a tolerance test. When tha calculated check value becomes greater than an assigned tolerance value, the loop should stop. However my check value appears to contain imaginary numbers which stops it from converging. Does anybody know how I could fix this? Or can anybody see any simple mistakes I have made as I am still relatively mnew to the world of Matlab? I have attached my code for clarity:

% Inputs
R=0.4;                                          % Radius of Rotor
B=3;                                            % Number of blades
U=1;                                            % 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;                                    % Design lift coefficient
% Variables
i=1;                                            % Counter
alpha_new=0;                                    % Initial value for alpha new
tolerance=0.00001;                                 % Tolerance Value
axial_induction=[0;0;0;0;0;0;0;0;0];
Check=1;                                        % Initial check value
axial_induction_old=0;                          % Initial value for old axial induction factor
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
r_local=R/9*(1:9)';
r_over_R=r_local / R;
ii=0;                                               % Counting Interval for integral values
for TSR=0:0.1:10                                     % TSR from 1 to 10
    disp(TSR)                                                                    
    Check=1;
ii=ii+1;
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) ./ (pi.*r_local.*2);
axial_induction= 1 ./ (((4.*(sin(Phi).^2)) ./ (sigma.*Cl_design.*cos(Phi)))+1);
angular_induction= (1-(3*axial_induction)) ./ ((4.*axial_induction)-1);
angular_velocity=TSR.*U./R;
     while abs(Check)>=tolerance
relative_wind = atan(U.*(1-axial_induction))./(angular_velocity*R).*(1+angular_induction);
F=(2/pi) .* 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);     
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
for i=1:length(C_T)
   if  C_T(i) > 0.96;
      axial_induction(i) = 1 / (((4*F(i)*cos(relative_wind(i))) / (sigma(i)*Cl(i)))-1);
     else
         axial_induction(i) = 1 / (1+(4*F(i)*(sin(relative_wind(i))^2)) / (sigma(i)*Cl(i)*cos(relative_wind(i)))); 
     end;
  end;
D=(8./(TSR.*9)).*(F.*(sin(Phi).^2).*(cos(Phi)-((TSR_local).*(sin(Phi)))).*(sin(Phi)+((TSR_local).*(cos(Phi)))).*(1-(Cd./Cl).*cot(Phi)).*(TSR_local.^2));
Cp=sum(D);
Diff=axial_induction-axial_induction_old;
Check=max(Diff(:));
Check
     end
store_sigma(:,ii)=sigma;
store_Phi(:,ii)=Phi;
store_TSR_local(:,ii)=TSR_local;
store_axial_induction(:,ii)=axial_induction;
store_angular_induction(:,ii)=angular_induction;
store_axial_induction_old(:,ii)=axial_induction_old;
store_relative_wind(:,ii)=relative_wind;
store_Check(:,ii)=Check;
store_Diff(:,ii)=Diff;
store_Cp(:,ii)=Cp;
store_TSR(:,ii)=TSR;
store_F(:,ii)=F;
end

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Geoff Hayes 2014-7-22

在 MATLAB Online 中打开

Kevin - you state that When tha calculated check value becomes greater than an assigned tolerance value, the loop should stop whereas your code is doing the opposite

while abs(Check)>=tolerance

continuing so long as the check value is greater than the assigned tolerance.

Note that when TSR is zero (on the first iteration of the outer for loop) there is a division by zero error which leads to a couple of vectors being set to all Inf and (later) all NaN values. You may want to skip this first case.

Your code is getting imaginary numbers at the calculation of F

 F=(2/pi) .* acos(exp(-(((B/2) .* (1-(r_over_R))) ./ ...
      ((r_over_R) .* sin(relative_wind)))));

Has this equation been written correctly?

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Michael Haderlein 2014-7-22

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/142888-while-loop-not-converging#answer_145946

If you have imaginary numbers, maybe the problem comes from the acos in the calculation of F? acos is rad based (not degree). Could that be the error? If you want to use degree, use the acosd function.

If that's not the problem, you should check stepwise where the imaginary numbers appear first (use breakpoints and F10).

Best,

Michael