Loop problem

% I have got this simple program runnning a loop to satisfy the condition 'nvvni = 0' stop the program
PT= 50;
c = [0.20,0.10,0.10;0.20,0.20,0.20;0 0 0];
P= [190 72.2 51.6; 20.44 15.57 4.956;0 0 0];
K = P/PT;
A = c(1,1)*(K(1,1)-1)+ c(1,2)*(K(1,2)-1)+c(1,3)*(K(1,3)-1)+c(2,1)*(K(2,1)-1)+c(2,2)*(K(2,2)-1)+c(2,3)*(K(2,3)-1);
B = (c(1,1)*(K(1,1)-1)/K(1,1))+(c(1,2)*(K(1,2)-1)/K(1,2))+(c(1,3)*(K(1,3)-1)/K(1,3))+(c(2,1)*(K(2,1)-1)/K(2,1))+(c(2,2)*(K(2,2)-1)/K(2,2))+(c(2,3)*(K(2,3)-1)/K(2,3));
%disp(B)
for i=1:10
n=0;
nvni = (A )/(A-B);
n=n+1;
nvvni = (c(1,1)*(K(1,1)))/(nvni*(K(1,1)-1)+1)+(c(1,2)*(K(1,2)))/(nvni*(K(1,2)-1)+1)+(c(1,3)*(K(1,3)))/(nvni*(K(1,3)-1)+1)+(c(2,1)*(K(2,1)))/(nvni*(K(2,1)-1)+1)+(c(2,2)*(K(2,2)))/(nvni*(K(2,2)-1)+1)+(c(2,3)*(K(2,3)))/(nvni*(K(2,3)-1)+1);
if nvvni == 0;
disp(nvni);
break
else
nvvvi = (c(1,1)*(K(1,1)^2)/(nvni*(K(1,1) + 1)^2))+(c(1,2)*(K(1,2)^2)/(nvni*(K(1,2) + 1)^2))+(c(1,3)*(K(1,3)^2)/(nvni*(K(1,3) + 1)^2))+(c(2,1)*(K(2,1)^2)/(nvni*(K(2,1) + 1)^2))+(c(2,2)*(K(2,2)^2)/(nvni*(K(2,2) + 1)^2))+(c(2,3)*(K(2,3)^2)/(nvni*(K(2,3) + 1)^2));
nvvni = nvni - (nvvni/nvvvi);
nvni =nvvni;
end
end

1 Comment

Jan
Jan on 19 Feb 2012
Please format your code properly. I've done this for you this time.
What is your question?

Sign in to comment.

Answers (1)

Jan
Jan on 19 Feb 2012
It is very unlikely that the result of such a complicated formula is exactly 0.0 due to rounding errors. Better check for a certain limit:
if abs(nvvni) < 1e-14
The size of the limit depends on the physical meaning of the formula and cannot be defined uniquely. An analysis of the sensitivity is recommended in addition: How sensitive is nvvni to tiny variations of the inputs, e.g. when using c+eps.

1 Comment

Walter Roberson
Walter Roberson on 19 Feb 2012
http://matlab.wikia.com/wiki/FAQ#Why_is_0.3_-_0.2_-_0.1_.28or_similar.29_not_equal_to_zero.3F

Sign in to comment.

Categories

Find more on Loops and Conditional Statements in Help Center and File Exchange

Tags

Asked:

on 19 Feb 2012

Edited:

on 16 Oct 2013

Community Treasure Hunt

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

Start Hunting!