code with while loop that start running but doesn't end

조회 수: 5 (최근 30일)
Moufid Tarek
Moufid Tarek 2023년 1월 20일
답변: Moufid Tarek 2023년 1월 20일
The code is designed to resolve (via a while loop) a least square adjustment problem, unknowns are a,b,c, observations y and constant x values such that : y(x) = a* exp(b*x) +c.
The code starts running but never ends, I don't know what to change..
Thank you in advance for your precious help.
% data: x and L are both 25x1 vector (scalars)
data = importdata('exercise9task2.txt');
x= data(:,1);
L=data(:,2);
%break-off conditions
epsilon = 10^-6;
delta = 10^-8;
max_x_hat = 1;
check2=1;
%Number of iterations
iteration = 0;
syms a1 b1 c1
bx= a1* exp(x.*b1 ) +c1;
while check2 > delta && max_x_hat > epsilon
%Observations as functions of the approximations for the unknowns
L_0 = subs(bx,{a1,b1,c1},{X_0(1),X_0(2),X_0(3)});
%Vector of reduced observations
l = L-L_0;
%Design matrix with the elements from the Jacobian matrix J
A = jacobian(bx) ;
A= subs(A,{a1,b1,c1},{X_0(1),X_0(2),X_0(3)});
%Normal matrix
N = A'*A ;
%Vector of right hand side of normal equations
n = A'*l;
%Inversion of normal matrix / Cofactor matrix of the unknowns
% Q_xx =inv(N) ;
%Solution of the normal equations
x_hat = inv(N) *n;
%Update
X_0 = x_hat + X_0;
%Check 1
max_x_hat = max(abs(x_hat));
%Vector of residuals
v = A*x_hat - l;
%Vector of adjusted observations
L_hat = L+v;
%Objective function
% vTPv = v'*v;
%Functional relationships witho ut the observations
%Check 2
Check2 = max(abs(L_hat - (subs(bx,{a1,b1,c1},{X_0(1),X_0(2),X_0(3)}))));
%Update number of iterations
iteration = iteration+1;
end
  댓글 수: 3
Moufid Tarek
Moufid Tarek 2023년 1월 20일
Thank you for your help, it was something I didn't even notice, but the problem remains
Moufid Tarek
Moufid Tarek 2023년 1월 20일
It seems like it is not able to compute the value of Check2, I tried running each line on its own, and the check2 line runs but doesn't end

댓글을 달려면 로그인하십시오.

답변 (2개)

Image Analyst
Image Analyst 2023년 1월 20일
While loops should always have a failsafe so you don't get an infinite loop. Here is an example of how to put a failsafe in your loop:
% Demonstration of how to avoid an infinite loop by setting up a failsafe.
% Set up a failsafe
maxIterations = 100; % Way more than you think it would ever need.
loopCounter = 0;
% Now loop until we obtain the required condition: a random number equals exactly 0.5.
% If that never happens, the failsafe will kick us out of the loop so we do not get an infinite loop.
r = nan; % Initialize so we can enter the loop the first time.
while (r ~= 0.5) && loopCounter < maxIterations
loopCounter = loopCounter + 1;
fprintf('Iteration #%d.\n', loopCounter)
r = rand; % This will never equal 0.500000000000000000000000000000 exactly
end
% Alert user if we exited normally, or if the failsafe kicked us out to avoid an infinite loop.
if loopCounter < maxIterations
% Then the loop found the condition and exited early, which means normally.
fprintf('Loop exited normally after %d iterations.\n', loopCounter);
else
% Then the loop never found the condition and exited when the number of iterations
% hit the maximum number of iterations allowed, which means abnormally.
fprintf('Loop exited abnormally after iterating the maximimum number of iterations (%d) without obtaining the exit criteria.\n', maxIterations);
end
fprintf('All done after %d iterations.\n', loopCounter)
  댓글 수: 1
Moufid Tarek
Moufid Tarek 2023년 1월 20일
Thank you, it is very good to know that, I will keep that in mind.
My code, on the otherhand, is not good enough to perfom more than one loop, if I stop the code from running and check the variable iteration, it is equal to 0

댓글을 달려면 로그인하십시오.


Moufid Tarek
Moufid Tarek 2023년 1월 20일
Thank you all for your support. I found the solution, everywhere I used the equation "bx" I had to convert the output to "double"
For instance:
L_0 = subs(bx,{a1,b1,c1},{X_0(1),X_0(2),X_0(3)});
% should be
L_0 = double(L_0); % before being used in other computations
% in order to make it work

카테고리

Help CenterFile Exchange에서 Loops and Conditional Statements에 대해 자세히 알아보기

태그

제품


릴리스

R2021b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by