Bugs in the FAQ....

[fungus]

I reported this bug years ago but it's still not fixed.

On this page: http://c-faq.com/fp/fpequal.html

it says to use the following code:

if (fabs(a - b) <= epsilon * fabs(a))

What happens when a is zero and b is (eg.) epsilon/2?

In this case a and b should compare as equal, but they don't...

:-)
--
<\___/>
/ O O \
\_____/ FTB. Remove my socks for email address.

[jacob navia]

fungus wrote:

I reported this bug years ago but it's still not fixed.

On this page: http://c-faq.com/fp/fpequal.html

it says to use the following code:

if (fabs(a - b) <= epsilon * fabs(a))

What happens when a is zero and b is (eg.) epsilon/2?

In this case a and b should compare as equal, but they don't...

:-)
--
<\___/>
/ O O \
\_____/ FTB. Remove my socks for email address.

If you read the next line it will say:

.... where epsilon is a value chosen to set the degree of ``closeness''
(and where you know that a will not be zero).

So a can't be zero
--
jacob navia
jacob at jacob point remcomp point fr
logiciels/informatique
http://www.cs.virginia.edu/~lcc-win32

[fungus]

On Jul 14, 4:50*pm, jacob navia <[e-mail]:

>
If you read the next line it will say:

... where epsilon is a value chosen to set the degree of ``closeness''
(and where you know that a will not be zero).

So a can't be zero

Ah, Ok. It didn't used to say that....

So is that the "fix" then? :-)
--
<\___/>
/ O O \
\_____/ FTB. Remove my socks for email address.

[jacob navia]

fungus wrote:

On Jul 14, 4:50 pm, jacob navia <[e-mail]:

>If you read the next line it will say:

... where epsilon is a value chosen to set the degree of ``closeness''
(and where you know that a will not be zero).

So a can't be zero

Ah, Ok. It didn't used to say that....

So is that the "fix" then? :-)
--
<\___/>
/ O O \
\_____/ FTB. Remove my socks for email address.

Well, it says that in the link you gave us.

Obviously it would be better if the code made no such an assumption like

if( a == 0 || fabs(a - b) <= epsilon * fabs(a))

but all the faqs have bugs and this one is not an exception. If you
reread the text it looks like the sentence within the parenthesis was
hastily added after someone discovered that...
well you know :-)

--
jacob navia
jacob at jacob point remcomp point fr
logiciels/informatique
http://www.cs.virginia.edu/~lcc-win32

[Willem]

fungus wrote:
) I reported this bug years ago but it's still not fixed.
)
) On this page: http://c-faq.com/fp/fpequal.html
)
) it says to use the following code:
)
) if (fabs(a - b) <= epsilon * fabs(a))
)
) What happens when a is zero and b is (eg.) epsilon/2?
)
) In this case a and b should compare as equal, but they don't...

Why should a and b compare as equal in that case ?
Consider the case where a = 0.001 and b = epsilon/2, for example.

However, a much more serious problem with the above code is that it
is not commutative. compare(a,b) is not always equal to compare(b,a).

To be commutative, shouldn't the code be something like:

if (fabs(a-b) <= epsilon * max(fabs(a), fabs(b)))

(Which, incidentally, solves the a=0 problem also)
SaSW, Willem
--
Disclaimer: I am in no way responsible for any of the statements
made in the above text. For all I know I might be
drugged or something..
No I'm not paranoid. You all think I'm paranoid, don't you !
#EOT

[fungus]

On Jul 14, 5:04*pm, jacob navia <[e-mail]:

>
Obviously it would be better if the code made no such an assumption

I changed it to:

if (((a == 0) && (b<epsilon)) || (fabs(a - b) <= epsilon * fabs(a)))
There's still a _tiny_ bug because:

(epsilon * fabs(a)) != (epsilon * fabs(b))

But it's a lot better....

well you know :-)

Yeah, I know...

--
<\___/>
/ O O \
\_____/ FTB. Remove my socks for email address.

[fungus]

On Jul 14, 5:13*pm, Willem <[e-mail]:

>
Why should a and b compare as equal in that case ?

Because "epsilon" is pretty close to zero...?
--
<\___/>
/ O O \
\_____/ FTB. Remove my socks for email address.

[Keith Thompson]

jacob navia <ja***@nospam.comwrites:

fungus wrote:

>On Jul 14, 4:50 pm, jacob navia <[e-mail]:

>>If you read the next line it will say:

... where epsilon is a value chosen to set the degree of ``closeness''
(and where you know that a will not be zero).

So a can't be zero

Ah, Ok. It didn't used to say that....
So is that the "fix" then? :-)

Well, it says that in the link you gave us.

Obviously it would be better if the code made no such an assumption like

if( a == 0 || fabs(a - b) <= epsilon * fabs(a))

but all the faqs have bugs and this one is not an exception. If you
reread the text it looks like the sentence within the parenthesis was
hastily added after someone discovered that...

[...]

Yes, but that's not the right fix. It assumes that a and b are
"equal" whenever a==0.

--
Keith Thompson (The_Other_Keith) ks***@mib.org <http://www.ghoti.net/~kst>
Nokia
"We must do something. This is something. Therefore, we must do this."
-- Antony Jay and Jonathan Lynn, "Yes Minister"

[Eric Sosman]

fungus wrote:

On Jul 14, 5:04 pm, jacob navia <[e-mail]:

>Obviously it would be better if the code made no such an assumption

I changed it to:

[The "it" is part of one of the FAQ's suggestions on how to
compare floating-point numbers for approximate equality]

if (((a == 0) && (b<epsilon)) || (fabs(a - b) <= epsilon * fabs(a)))
There's still a _tiny_ bug because:

... because it says 0 and -1E30 are approximately equal?

--
Er*********@sun.com

[user923005]

On Jul 14, 7:47 am, fungus <[e-mail]:

I reported this bug years ago but it's still not fixed.

On this page:http://c-faq.com/fp/fpequal.html

it says to use the following code:

if (fabs(a - b) <= epsilon * fabs(a))

What happens when a is zero and b is (eg.) epsilon/2?

In this case a and b should compare as equal, but they don't...

:-)

It depends on whether you are trying to calculate relative or absolute
error.

This is one possible relative error comparison for floating point
types:

#include <float.h>
#include <math.h>

int double_compare (double d1, double d2)
{
if (d1 d2)
if ((d1 - d2) < fabs (d1 * DBL_EPSILON))
return 0;
else
return 1;
if (d1 < d2)
if ((d2 - d1) < fabs (d2 * DBL_EPSILON))
return 0;
else
return -1;
return 0;
}

int float_compare (float d1, float d2)
{
if (d1 d2)
if ((d1 - d2) < fabsf (d1 * FLT_EPSILON))
return 0;
else
return 1;
if (d1 < d2)
if ((d2 - d1) < fabsf (d2 * FLT_EPSILON))
return 0;
else
return -1;
return 0;
}

On the other hand, maybe you do not want relative error.
An example might be if you have a floating point function that is
given to a root search algorithm.
When the root has a y value of DBL_EPSILON you have found the zero and
you do not want to keep searching or refining. So the above method
would be incorrect for that application.
I don't think that there is a simple answer to the problem.