setdiff not working for a particular value, bug?

조회 수: 6(최근 30일)
Set 1:
R = 1.2;
F = 1.78;
D = 1.29;
M = 0.2
M = 0.2000
Set 2:
R = 1.2;
F = 2.5;
D = 1.59;
M = 0.5
M = 0.5000
And here's the function:
V= R/(F*(1/(M)-1)+R);
M = 0.5000
X = setdiff( 1:0.005:2,D);
err = (abs((((((D-(1-V).*X)./V)).*(X./((((D-(1-V).*X)./V)).*(1./M-1)+X))+(1-(X./((((D-(1-V).*X)./V)).*(1./M-1)+X))).*X)-D)./D)).*100;
[~,imin]=min(err);
New=X(imin);
For set 1, it gives the right value for New = 1.205 but for set 2 it gives New = D = 1.59.
Why is that?

채택된 답변

Walter Roberson
Walter Roberson 2021년 12월 9일
not a bug.
D = 1.59;
vec = 1:0.005:2;
[~, idx] = min(abs(vec - D));
vec(idx)
ans = 1.5900
vec(idx) - D
ans = -2.2204e-16
You can see that the closest value to D in vec is about 2.22e-16 less than D.
Your colon operator is operating like this:
start = 1;
stop = 2;
incr = 0.005;
if start > stop
vec = [];
else
current = start;
vec(end+1) = current;
new = current + incr;
if new > stop; break; end
current = new;
end
Now if (1.59 - 1) were exactly divisibile by the binary double precision representation of 0.005 and assuming that start (1.0) is exactly representable in binary double precision, then under those circumstances 1.59 exactly would eventually be reached.
However... the binary double precision representation of 0.005 is 0.005000000000000000104083408558608425664715468883514404296875 exactly. Binary double precision is not able to represent 0.001 exactly . And those differences add up.
The reason that binary double precision is not able to represent 0.001 exactly is the same as the reason that a finite decimal floating point representation cannot exactly represent 1/3 or 1/7 . Suppose that you had a system that stored 10 decimal digits, then 1/3 would be 0.3333333333 . Add another of the same and you would get 0.6666666666 . Add the third value on and you would get 0.9999999999 -- which is not 1.0000000000 ! Algebraic rational 1/3 can only be exactly represented in a finite number of digits if the digit base is divisible by 3, which is not the case for decimal (base 10). Algebraic rational 1/10 can only be exactly represented in a finite number of digits i the digit base is divisible by 10, which is not the case for binary (base 2.)
You need to learn this rule:
Never compare floating point numbers for equality if the floating point numbers might have been calculated different ways. Even just very slight differences in computation like .3 - .2 - .1 versus .3 - .1 - .2 make a difference in floating point calculations.
  댓글 수: 1
Walter Roberson
Walter Roberson 2021년 12월 9일
By the way, because of floating point round off if you take 1/49*49 then the result is not exactly 1.0, so do not trust visual examination.

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

추가 답변(1개)

Pelajar UM
Pelajar UM 2021년 12월 9일
OK. Seems like it's a tolerance issue. It can be solved by using smaller steps:
X = setdiff( 1:0.000001:2,D);
  댓글 수: 1
Stephen23
Stephen23 2021년 12월 9일
"It can be solved by using smaller steps:"
No, that does not solve the problem. If you want to write robust code then you need to avoid testing for exact equivalence of binary floating point numbers (i.e. avoid SETDIFF, EQ, ISMEMBER, etc.)
Instead compare the absolute difference against a tolerance (selected to suit your data):
abs(A-B)<tol

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

Community Treasure Hunt

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

Start Hunting!

Translated by