comparing doubles for equality

[John Smith]

This code for the comparison of fp types is taken from the C FAQ.
Any problems using it in a macro?

/* compare 2 doubles for equality */
#define DBL_ISEQUAL(a,b) (fabs((a)-(b))<=(DBL_EPSILON)*fabs((a)))

Do the same issues involved in comparing 2 fp types for equality
apply to comparing a float to zero? E.g. is if(x == 0.0)
considered harmful?

[Tim Prince]

John Smith wrote:

This code for the comparison of fp types is taken from the C FAQ.
Any problems using it in a macro?

/* compare 2 doubles for equality */
#define DBL_ISEQUAL(a,b) (fabs((a)-(b))<=(DBL_EPSILON)*fabs((a)))

Do the same issues involved in comparing 2 fp types for equality
apply to comparing a float to zero? E.g. is if(x == 0.0) considered
harmful?

Depending on how you use this, you could make your application highly
dependent on external issues, such as whether you compile for an extra
precision mode, or disable gradual underflow.

[Dik T. Winter]

In article <I8Clh.533642$1T2.401780@pd7urf2noJohn Smith <JS****@mail.netwrites:

This code for the comparison of fp types is taken from the C FAQ.
Any problems using it in a macro?

/* compare 2 doubles for equality */
#define DBL_ISEQUAL(a,b) (fabs((a)-(b))<=(DBL_EPSILON)*fabs((a)))

Do the same issues involved in comparing 2 fp types for equality
apply to comparing a float to zero? E.g. is if(x == 0.0)
considered harmful?

If you insert 0.0 for 'a' you will see that the two give precisely the
same result. On the other hand, I think it is possible to construct
two floating point numbers 'a' and 'b', where the result is asymmetric.
But let that not deter you from using the macro, when it is asymmetric
you are in the outskirts of floating-point arithmetic.
--
dik t. winter, cwi, kruislaan 413, 1098 sj amsterdam, nederland, [phone]
home: bovenover 215, 1025 jn amsterdam, nederland; http://www.cwi.nl/~dik/

[Thad Smith]

John Smith wrote:

This code for the comparison of fp types is taken from the C FAQ.
Any problems using it in a macro?

/* compare 2 doubles for equality */
#define DBL_ISEQUAL(a,b) (fabs((a)-(b))<=(DBL_EPSILON)*fabs((a)))

This construction is misleading and I would never use it, because the
implied function, determining whether two doubles are equal, is not an
accurate description of the returned value.

Do the same issues involved in comparing 2 fp types for equality
apply to comparing a float to zero? E.g. is if(x == 0.0) considered
harmful?

Which issues are those? The test "if (x == 0.0)", in contrast, does not
have the inaccurate description that the DBL_ISEQUAL macro does.

In both cases, you should employ good analysis for floating point
comparisons.

--
Thad

[CBFalconer]

Tim Prince wrote:

John Smith wrote:

>This code for the comparison of fp types is taken from the C FAQ.
Any problems using it in a macro?

/* compare 2 doubles for equality */
#define DBL_ISEQUAL(a,b) (fabs((a)-(b))<=(DBL_EPSILON)*fabs((a)))

Do the same issues involved in comparing 2 fp types for equality
apply to comparing a float to zero? E.g. is if(x == 0.0)
considered harmful?

Depending on how you use this, you could make your application
highly dependent on external issues, such as whether you compile
for an extra precision mode, or disable gradual underflow.

That's why gcc has the Wfloat-equal switch. As far as the macro is
concerned, the a argument better not have any side effects.

--
Merry Christmas, Happy Hanukah, Happy New Year
Joyeux Noel, Bonne Annee.
Chuck F (cbfalconer at maineline dot net)
<http://cbfalconer.home.att.net>

[CBFalconer]

Thad Smith wrote:

John Smith wrote:

>This code for the comparison of fp types is taken from the C FAQ.
Any problems using it in a macro?

/* compare 2 doubles for equality */
#define DBL_ISEQUAL(a,b) (fabs((a)-(b))<=(DBL_EPSILON)*fabs((a)))

This construction is misleading and I would never use it, because
the implied function, determining whether two doubles are equal,
is not an accurate description of the returned value.

How about:

#define DUNEQUAL(a, b) (fabs((a)-(b)) (DBL_EPSILON)*fabs((a)))

with a caveat against passing an a with side effects.

--
Merry Christmas, Happy Hanukah, Happy New Year
Joyeux Noel, Bonne Annee.
Chuck F (cbfalconer at maineline dot net)
<http://cbfalconer.home.att.net>

[Thad Smith]

CBFalconer wrote:

Thad Smith wrote:

>John Smith wrote:

>>This code for the comparison of fp types is taken from the C FAQ.
Any problems using it in a macro?

/* compare 2 doubles for equality */
#define DBL_ISEQUAL(a,b) (fabs((a)-(b))<=(DBL_EPSILON)*fabs((a)))

This construction is misleading and I would never use it, because
the implied function, determining whether two doubles are equal,
is not an accurate description of the returned value.

How about:

#define DUNEQUAL(a, b) (fabs((a)-(b)) (DBL_EPSILON)*fabs((a)))

with a caveat against passing an a with side effects.

That misses the point. The only true equality test is given by a==b.
If you need an epsilon, fine, but you should make it explicit.
DBL_EPSILON isn't necessarily the correct choice for a particular
application. In fact, the two tests above may not provide any advantage
over a strict equality.

--
Thad

[CBFalconer]

Thad Smith wrote:

CBFalconer wrote:

>Thad Smith wrote:

>>John Smith wrote:

This code for the comparison of fp types is taken from the C FAQ.
Any problems using it in a macro?

/* compare 2 doubles for equality */
#define DBL_ISEQUAL(a,b) (fabs((a)-(b))<=(DBL_EPSILON)*fabs((a)))
This construction is misleading and I would never use it, because
the implied function, determining whether two doubles are equal,
is not an accurate description of the returned value.

How about:

#define DUNEQUAL(a, b) (fabs((a)-(b)) (DBL_EPSILON)*fabs((a)))

with a caveat against passing an a with side effects.

That misses the point. The only true equality test is given by
a==b. If you need an epsilon, fine, but you should make it
explicit. DBL_EPSILON isn't necessarily the correct choice for a
particular application. In fact, the two tests above may not
provide any advantage over a strict equality.

No, I think you miss the point. Floating point is inherently an
approximation, and the above says that anything within the
resolution (i.e. the approximation) is to be considered equal. Any
time you see (a == b) for floating point operands, it is probably a
bug. For example:

for (a = 1.0; a != 0; a -= 0.1) dosomethingwith(a);

will probably not behave. While:

for (a = 1.0; !DUNEQUAL(a, 0.0); a -= 0.1) dosomethingwith(a);

will behave. Of course:

for (a = 1.0; a 0; a -= 0.1) dosomethingwith(a);

may behave. But it is a crapshoot whether it does an extra cycle.

--
Merry Christmas, Happy Hanukah, Happy New Year
Joyeux Noel, Bonne Annee.
Chuck F (cbfalconer at maineline dot net)
<http://cbfalconer.home.att.net>

[SM Ryan]

John Smith <JS****@mail.netwrote:
# This code for the comparison of fp types is taken from the C FAQ.
# Any problems using it in a macro?
#
# /* compare 2 doubles for equality */
# #define DBL_ISEQUAL(a,b) (fabs((a)-(b))<=(DBL_EPSILON)*fabs((a)))
#
# Do the same issues involved in comparing 2 fp types for equality
# apply to comparing a float to zero? E.g. is if(x == 0.0)
# considered harmful?

Floats and doubles use a small number of bits for their values,
so exact equality to zero or any other representable value is
meaningful. Two problems are

(1) most machines do not use radix 10 representation; some decimal
values, such as 0.1 can only be approximated. That means if you
write (x==0.1) you're not really testing x against one tenth but
some other value that is close to one tenth. So the test is not
what you think it is.

(2) many computations with reals are iterative approximations to
an unknown value. Because of the precision limits, they often
cannot reach the exact correct value, so you have to accept when
it is close enough.

--
SM Ryan http://www.rawbw.com/~wyrmwif/
So basically, you just trace.

[Dik T. Winter]

In article <45***************@yahoo.comcb********@maineline.net writes:
....

#define DUNEQUAL(a, b) (fabs((a)-(b)) (DBL_EPSILON)*fabs((a)))

....

For example:
for (a = 1.0; a != 0; a -= 0.1) dosomethingwith(a);
will probably not behave. While:
for (a = 1.0; !DUNEQUAL(a, 0.0); a -= 0.1) dosomethingwith(a);

you wish to omit the exclamation mark here.

will behave. Of course:
for (a = 1.0; a 0; a -= 0.1) dosomethingwith(a);
may behave. But it is a crapshoot whether it does an extra cycle.

And now compare the three with:
for (a = 1.0; DUNEQUAL(0.0, a); a -= 0.1) dosomethingwith(a);

A better test would be:

#define DUNEQUAL(a, b) (fabs((a) - (b)) \
(fabs(a) fabs(b) ? DBL_EPSILON * fabs(a) : DBL_EPSILON * fabs(b)))

in this way there is no asymmetry.
--
dik t. winter, cwi, kruislaan 413, 1098 sj amsterdam, nederland, [phone]
home: bovenover 215, 1025 jn amsterdam, nederland; http://www.cwi.nl/~dik/

[Thad Smith]

CBFalconer wrote:

Thad Smith wrote:

>CBFalconer wrote:

>>Thad Smith wrote:
John Smith wrote:

This code for the comparison of fp types is taken from the C FAQ.
Any problems using it in a macro?
>
/* compare 2 doubles for equality */
#define DBL_ISEQUAL(a,b) (fabs((a)-(b))<=(DBL_EPSILON)*fabs((a)))
This construction is misleading and I would never use it, because
the implied function, determining whether two doubles are equal,
is not an accurate description of the returned value.
How about:

#define DUNEQUAL(a, b) (fabs((a)-(b)) (DBL_EPSILON)*fabs((a)))

with a caveat against passing an a with side effects.

That misses the point. The only true equality test is given by
a==b. If you need an epsilon, fine, but you should make it
explicit. DBL_EPSILON isn't necessarily the correct choice for a
particular application. In fact, the two tests above may not
provide any advantage over a strict equality.

No, I think you miss the point. Floating point is inherently an
approximation, and the above says that anything within the
resolution (i.e. the approximation) is to be considered equal. Any
time you see (a == b) for floating point operands, it is probably a
bug. For example:

for (a = 1.0; a != 0; a -= 0.1) dosomethingwith(a);

While I agree that this is likely to fail ...

for (a = 1.0; !DUNEQUAL(a, 0.0); a -= 0.1) dosomethingwith(a);

.... there is no guarantee that this second form is correct, for the
definition supplied above. Why? Because each subtraction may invoke a
roundoff. After several roundings, you maybe off by more than one lsb
for a number with the magnitude of 1, let alone, say, one lsb of 1e-7.

Note also that DUNEQUAL, is not commutative, which is fails the
principle of least astonishment.

>
for (a = 1.0; a 0; a -= 0.1) dosomethingwith(a);

may behave. But it is a crapshoot whether it does an extra cycle.

And, to make it robust, you could do

for (a = 1.0; a 0.1/2; a -= 0.1) dosomethingwith(a);

or, to use the best approximation of of each nominal value of a:

int i;
for (i = 10; i 0; i--) { dosomethingwith (i*0.1); }

or, for the pragmatic types

for (a = 10; a 0; a -= 1) dosomethingwith (a/10);

The last form is not guaranteed, but performs correctly on any floating
point system I am familiar with, except maybe logarithmic systems, due
to the limited number of bits required to exactly represent the integer
values.

--
Thad

[CBFalconer]

"Dik T. Winter" wrote:

cb********@maineline.net writes:
...

>>>#define DUNEQUAL(a, b) (fabs((a)-(b)) (DBL_EPSILON)*fabs((a)))

...

>For example:
for (a = 1.0; a != 0; a -= 0.1) dosomethingwith(a);
will probably not behave. While:
for (a = 1.0; !DUNEQUAL(a, 0.0); a -= 0.1) dosomethingwith(a);

you wish to omit the exclamation mark here.

>will behave. Of course:
for (a = 1.0; a 0; a -= 0.1) dosomethingwith(a);
may behave. But it is a crapshoot whether it does an extra cycle.

And now compare the three with:
for (a = 1.0; DUNEQUAL(0.0, a); a -= 0.1) dosomethingwith(a);

A better test would be:

#define DUNEQUAL(a, b) (fabs((a) - (b)) \
(fabs(a) fabs(b) ? DBL_EPSILON * fabs(a) : DBL_EPSILON * fabs(b)))

in this way there is no asymmetry.

Not needed. The use of epsilon is only needed when a and b are
nearly equal, thus either can be used to derive it. The problem of
a zero remains, and ties in with the use of extremely small
operands in general. To handle all these cases we probably should
use a function rather than a macro.

--
Merry Christmas, Happy Hanukah, Happy New Year
Joyeux Noel, Bonne Annee.
Chuck F (cbfalconer at maineline dot net)
<http://cbfalconer.home.att.net>

[pete]

CBFalconer wrote:

>
Thad Smith wrote:

CBFalconer wrote:

Thad Smith wrote:
John Smith wrote:

This code for the comparison of fp types is taken from the C FAQ.
Any problems using it in a macro?

/* compare 2 doubles for equality */
#define DBL_ISEQUAL(a,b) (fabs((a)-(b))<=(DBL_EPSILON)*fabs((a)))
This construction is misleading and I would never use it, because
the implied function, determining whether two doubles are equal,
is not an accurate description of the returned value.

How about:

#define DUNEQUAL(a, b) (fabs((a)-(b)) (DBL_EPSILON)*fabs((a)))

with a caveat against passing an a with side effects.

That misses the point. The only true equality test is given by
a==b. If you need an epsilon, fine, but you should make it
explicit. DBL_EPSILON isn't necessarily the correct choice for a
particular application. In fact, the two tests above may not
provide any advantage over a strict equality.

No, I think you miss the point. Floating point is inherently an
approximation, and the above says that anything within the
resolution (i.e. the approximation) is to be considered equal. Any
time you see (a == b) for floating point operands, it is probably a
bug.

I disagree.
It really all depends
on what you're doing with the floating point values.
If you're using qsort to sort an array of doubles,
then there's really no place for DBL_EPSILON in your compar function.

int compar(const void *arg1, const void *arg2)
{
return *(double *)arg2 *(double *)arg1 ? -1
: *(double *)arg2 != *(double *)arg1;
}

--
pete