compare two float values

U

user923005

This version has better formatting, which will allow us to see the
dirt.
That makes the output much easier to grok.

#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;
}

#ifdef UNIT_TEST

#include <stdio.h>

int main ()
{
double u = 1.0;
double v = u + 1e-14;
double big = 1.0;
double a = 1.0, b = 0.99999999999999999999999999;
printf ("compare(%.*g, %.*g) = %d.\n", DBL_DIG+5, a, DBL_DIG+5, b,
double_compare (a, b));
a *= -1;
b *= -1;
printf ("compare(%.*g, %.*g) = %d.\n", DBL_DIG+5, a, DBL_DIG+5, b,
double_compare (a, b));

big *= 1 << 16;
big *= 1 << 16;
big *= 1 << 16;
big *= 1 << 16;


printf ("compare(%.*g, %.*g) = %d.\n", DBL_DIG+5, u, DBL_DIG+5, v,
double_compare (u, v));


u *= big;
v *= big;

printf ("compare(%.*g, %.*g) = %d.\n", DBL_DIG+5, u, DBL_DIG+5, v,
double_compare (u, v));

u = DBL_EPSILON;
v = -DBL_MAX;
printf ("compare(%.*g, %.*g) = %d.\n", DBL_DIG+5, u, DBL_DIG+5, v,
double_compare (u, v));

return 0;
}
/*
compare(1, 1) = 0.
compare(-1, -1) = 0.
compare(1, 1.00000000000001) = -1.
compare(18446744073709552000, 18446744073709736000) = -1.
compare(2.2204460492503131e-016, -1.7976931348623157e+308) = 1.
*/
#endif
 
K

Keith Thompson

Mark McIntyre said:
By some different means. Two mathematically identical floats may
compare inequal by obvious methods.

Depends on what you mean by "mathematically identical".
Doubtful, but he also probably expected the below to print "true".

#include <stdio.h>
#include <math.h>
int main(void)
{
printf("%f %f %s\n",
exp(log(4.444)), 4.444,
exp(log(4.444))==4.444?"true":"false");
return 0;
}

$ ./a.out
4.444000 4.444000 false

Here you have two floating-point values that result from calculations
that, if they were carried out with complete mathematical exactness,
would be equal. But the actual values are not equal, because
floating-point calculations aren't completely exact. In this case,
you can see that they're inexact by using a different format. For
example:

#include <stdio.h>
#include <math.h>
int main(void)
{
printf("%.32f\n%.32f\n%s\n",
exp(log(4.444)),
4.444,
exp(log(4.444))==4.444 ? "true" : "false");
return 0;
}

4.44399999999999906208358879666775
4.44399999999999995026200849679299
false

If you want to know whether two floating-point values are "equal", the
first thing you have to ask yourself is just what that means. The
usual recommendation is to test whether the values are close enough --
but how close is close enough depends on the nature of the
calculations that generated the values. If the tolerance is too
tight, numbers that are logically equal might appear to be unequal; if
it's too loose, numbers that are logically distinct might appear to be
equal.

And sometimes exact equality (the built-in "==" operator) is just what
you want. For example, with a comparison function that checks whether
two numbers are within some epsilon of each other, you can have X == Y
and Y == Z, but X != Z. Whether that's a problem depends on the
algorithm.

Floating-point arithmetic is very probably harder than you think it is
unless your name is David Goldberg, and perhaps even if it is.
Speaking of which, see David Goldberg's paper "What Every Computer
Scientist Should Know About Floating-Point Arithmetic",
<http://docs.sun.com/source/806-3568/ncg_goldberg.html>.
 
J

Joachim Schmitz

user923005 said:
This version has better formatting, which will allow us to see the
dirt.
That makes the output much easier to grok.

#include <float.h>
#include <math.h>
int float_compare (float d1, float d2)
{
if (d1 > d2)
if ((d1 - d2) < fabsf (d1 * FLT_EPSILON))
<snip>

Where is fabsf()? Doesn't seem to be avaible on the platform I'm using,
which claims to be C89.
What do use instead?

Bye, Jojo
 
C

CBFalconer

Joachim said:
user923005 said:
user923005 said:
[snip]

What if d1 == DBL_EPSILON and d2 == -LARGE_NUMBER? I think you
are trying to do too much in one routine. I also think I have
the foul-up condition fouled. :)

This version has better formatting, which will allow us to see
the dirt. That makes the output much easier to grok.

#include <float.h>
#include <math.h>
int float_compare (float d1, float d2) {
if (d1 > d2)
if ((d1 - d2) < fabsf (d1 * FLT_EPSILON))
<snip>

Where is fabsf()? Doesn't seem to be avaible on the platform I'm
using, which claims to be C89. What do use instead?
From N869, for C99. For C90 use doubles and fabs.

7.12.7.2 The fabs functions
Synopsis
[#1]
#include <math.h>
double fabs(double x);
float fabsf(float x);
long double fabsl(long double x);
 
K

karthikbalaguru

Please quote relevant parts of the post to which you are replying. This
is good practise on Usenet.

Comparing floats for equality is fraught with danger as many values are
not perfectly representible in binary floating point. You might find
that you need to allow for a small variation on either side of the
comparison value, just like the margin of error in scientific
experiments.

You might want to read:

<http://www.cygnus-software.com/papers/comparingfloats/Comparing flo...>
<http://docs.sun.com/source/806-3568/ncg_goldberg.html>

I like the links. Nice Collection :):)

Karthik Balaguru
 
U

user923005

And surprise, surprise, it is: 14.5, <http://c-faq.com/fp/fpequal.html>.

Unfortunately, it's hopelessly broken:
C:\tmp>reldif
RelDif of 1 and -1 is 2

C:\tmp>type reldif.c
#define Abs(x) ((x) < 0 ? -(x) : (x))
#define Max(a, b) ((a) > (b) ? (a) : (b))

double RelDif(double a, double b)
{
double c = Abs(a);
double d = Abs(b);

d = Max(c, d);

return d == 0.0 ? 0.0 : Abs(a - b) / d;
}

int main(void)
{
printf("RelDif of %.20g and %.20g is %.20g\n", 1.0, -1.0,
RelDif(1.0, -1.0));
return 0;
}
 
U

user923005

Unfortunately, it's hopelessly broken:
C:\tmp>reldif
RelDif of 1 and -1 is 2

C:\tmp>type reldif.c
#define Abs(x)    ((x) < 0 ? -(x) : (x))
#define Max(a, b) ((a) > (b) ? (a) : (b))

double          RelDif(double a, double b)
{
    double          c = Abs(a);
    double          d = Abs(b);

    d = Max(c, d);

    return d == 0.0 ? 0.0 : Abs(a - b) / d;

}

int             main(void)
{
    printf("RelDif of %.20g and %.20g is %.20g\n", 1.0, -1.0,
RelDif(1.0, -1.0));
    return 0;



}

Nevermind. It works the way it was designed to work.
 
M

Mark McIntyre

Keith said:
Depends on what you mean by "mathematically identical".

Er, I mean mathematically identical, as in the same.
Here you have two floating-point values that result from calculations
that, if they were carried out with complete mathematical exactness,
would be equal. But the actual values are not equal, because
floating-point calculations aren't completely exact. In this case,
you can see that they're inexact by using a different format.


printf("%.32f\n%.32f\n%s\n",

Note that other than as a demo of the issue, I consider this meaningless
since one is printing greater than the precision of the type.
Floating-point arithmetic is very probably harder than you think

Did you confuse me with the OP?
unless your name is David Goldberg,

Its not, but I did spend several years using numerical analysis as part
of my post-doc research so I hope I'm fairly familiar with FP....
 
U

user923005

What's wrong with that? What answer did you expect, and why?

Nothing is wrong with it. I just expected it to work like a standard
compare() function. Once I actually understood what was intended, it
is clear that it works as designed.
 
J

Jack Klein

float t, check;
...
if ( t > check)
printf("%f is greater than %f\n" t, check);
else
printf("%f is less than %f\n" t, check);

Don't compare for equality, won't work for float (or double), rather compare

That's a little strong. It won't work for floating point values in
some circumstances. It most certainly will in others.
the difference against some delta that suits your needs for accuracy

--
Jack Klein
Home: http://JK-Technology.Com
FAQs for
comp.lang.c http://c-faq.com/
comp.lang.c++ http://www.parashift.com/c++-faq-lite/
alt.comp.lang.learn.c-c++
http://www.club.cc.cmu.edu/~ajo/docs/FAQ-acllc.html
 
K

Keith Thompson

Mark McIntyre said:
Er, I mean mathematically identical, as in the same.

In your example program, the two values you compared would be
identical if computed with mathematically perfect accuracy. The
values (trivial quibble: they're of type double, not float) are not
identical, because the operations that produced them are not perfectly
accurate.

I'm certain I'm not telling you anything you don't know perfectly
well; I'm trying to explain why I disagreed with your referring to
them as "mathematically identical floats". I suppose you could call
them mathematically identical real numbers, to which the double values
produced at runtime are imperfect approximations.

[...]
Note that other than as a demo of the issue, I consider this
meaningless since one is printing greater than the precision of the
type.

Yes, it was a demo of the issue.
Did you confuse me with the OP?

That was meant to be a generic "you", not specifically referring to
you.
Its not, but I did spend several years using numerical analysis as
part of my post-doc research so I hope I'm fairly familiar with FP....

I believe that the majority of programmers underestimate the
complexity of floating-point arithmetic. I certainly don't exclude
myself from that. If you're an exception to that, I congratulate you.
 
C

christian.bau

float t, check;
...
if ( t > check)
    printf("%f is greater than %f\n" t, check);
else
    printf("%f is less than %f\n" t, check);

Don't compare for equality, won't work for float (or double), rather compare
the difference against some delta that suits your needs for accuracy

You are doing nothing but perpetuating stupid myths. In a correct C
implementation (and there are a few incorrect ones around), t < check
is true if t is less than check, t == check is true if t is equal to
check, t > check is true if t is greater than check. In an IEEE 754
compatible implementation, exactly one of these three will be true
unless one of t or check is Not-A-Number, in which case neither
condition is true.

Once you accept this fact, you then need to accept that floating-point
operations are usually not exact, and two calculations that you think
should give the same result will often give results that are close
together but not necessarily the same. For example, if you write

double x = sqrt (2.0) * sqrt (3.0);
double y = sqrt (6.0);

you can be sure that the results will be close together, but which of
x < y, x == y and x > y is true is hard to predict. But the point is:
If they are equal, then x == y is true, and if they are not equal then
x == y is false. The == operator will give you the correct result,
except that the correct result may not be what you thought the exact
result is.
 
U

user923005

You are doing nothing but perpetuating stupid myths. In a correct C
implementation (and there are a few incorrect ones around), t < check
is true if t is less than check, t == check is true if t is equal to
check, t > check is true if t is greater than check. In an IEEE 754
compatible implementation, exactly one of these three will be true
unless one of t or check is Not-A-Number, in which case neither
condition is true.

Once you accept this fact, you then need to accept that floating-point
operations are usually not exact, and two calculations that you think
should give the same result will often give results that are close
together but not necessarily the same. For example, if you write

double x = sqrt (2.0) * sqrt (3.0);
double y = sqrt (6.0);

you can be sure that the results will be close together, but which of
x < y, x == y and x > y is true is hard to predict. But the point is:
If they are equal, then x == y is true, and if they are not equal then
x == y is false. The == operator will give you the correct result,
except that the correct result may not be what you thought the exact
result is.

While what you said is correct, I seen nothing wrong in the post of
"Joachim Schmitz".

In my opinion, any relational operation against floating point values
(either float or double) is a mistake unless it takes into account the
probability of floating point errors in calculation/conversion/
storage. If it were simply two constants stored in the program, we
can solve it by inspection and there is no need even to process the
test. Consider:

#include <stdio.h>
#include <stdlib.h>

int main(void)
{
double da,
db;
float fa,
fb;
double num = (double) rand();
double den = (double) rand() + 1.0;
da = num / den;
num = (double) rand();
den = (double) rand() + 1.0;
db = num / den;
num = (double) rand();
den = (double) rand() + 1.0;
fa = (float) (num / den);
num = (double) rand();
den = (double) rand() + 1.0;
fb = (float) (num / den);

if (da == db)
(void)puts("equal");
else if (da <= db)
(void)puts("less than or equal");
else if (da >= db)
(void)puts("greater than or equal");
else if (da < db)
(void)puts("less than");
else if (da > db)
(void)puts("greater than");
else if (da != db)
(void)puts("not equal");
else
(void)puts("probably a NAN is involved.");

if (fa == fb)
(void)puts("equal");
else if (fa <= fb)
(void)puts("less than or equal");
else if (fa >= fb)
(void)puts("greater than or equal");
else if (fa < fb)
(void)puts("less than");
else if (fa > fb)
(void)puts("greater than");
else if (fa != fb)
(void)puts("not equal");
else
(void)puts("probably a NAN is involved.");

return 0;
}
/*
C:\tmp>splint foo.c
Splint 3.1.1 --- 12 Mar 2007

foo.c: (in function main)
foo.c(23,9): Dangerous comparison involving double types: da == db
Two real (float, double, or long double) values are compared
directly using a
C primitive. This may produce unexpected results since floating
point
representations are inexact. Instead, compare the difference to
FLT_EPSILON
or DBL_EPSILON. (Use -realcompare to inhibit warning)
foo.c(29,14): Dangerous comparison involving double types: da < db
foo.c(31,14): Dangerous comparison involving double types: da > db
foo.c(33,14): Dangerous comparison involving double types: da != db
foo.c(38,9): Dangerous comparison involving float types: fa == fb
foo.c(44,14): Dangerous comparison involving float types: fa < fb
foo.c(46,14): Dangerous comparison involving float types: fa > fb
foo.c(48,14): Dangerous comparison involving float types: fa != fb

Finished checking --- 8 code warnings
*/

It is even possible for a number to not compare equal to itself. If
it were stored in an extended precision register, and then an IRQ
fires the extended precision bits could get shaved off when it is
stored the second time. At any rate, comparisons involving floating
point are very unlikely to work as desired unless you are very, very
careful.

Consider also:

C:\tmp>bar
num:
3
den:
3
num:
5
den:
5
num:
7
den:
7
num:
11
den:
11
greater than or equal
not equal

C:\tmp>type bar.c
#include <stdio.h>
#include <stdlib.h>

static char cNum[256];
static char cDen[256];

int main(void)
{
double da,
db;
float fa,
fb;
double num;
double den;
(void)puts("num:");
num = atof(fgets(cNum, sizeof cNum, stdin));
(void)puts("den:");
den = atof(fgets(cNum, sizeof cDen, stdin));
da = 1.0 / den;
da *= num;
(void)puts("num:");
num = atof(fgets(cNum, sizeof cNum, stdin));
(void)puts("den:");
den = atof(fgets(cNum, sizeof cDen, stdin));
da = 1.0 / den;
da *= num;
(void)puts("num:");
num = atof(fgets(cNum, sizeof cNum, stdin));
(void)puts("den:");
den = atof(fgets(cNum, sizeof cDen, stdin));
fa = 1.0 / den;
fa *= num;
(void)puts("num:");
num = atof(fgets(cNum, sizeof cNum, stdin));
(void)puts("den:");
den = atof(fgets(cNum, sizeof cDen, stdin));
fa = 1.0 / den;
fa *= num;

if (da == db)
(void) puts("equal");
else if (da <= db)
(void) puts("less than or equal");
else if (da >= db)
(void) puts("greater than or equal");
else if (da < db)
(void) puts("less than");
else if (da > db)
(void) puts("greater than");
else if (da != db)
(void) puts("not equal");
else
(void) puts("probably a NAN is involved.");

if (fa == fb)
(void) puts("equal");
else if (fa <= fb)
(void) puts("less than or equal");
else if (fa >= fb)
(void) puts("greater than or equal");
else if (fa < fb)
(void) puts("less than");
else if (fa > fb)
(void) puts("greater than");
else if (fa != fb)
(void) puts("not equal");
else
(void) puts("probably a NAN is involved.");

return 0;
}
 

Ask a Question

Want to reply to this thread or ask your own question?

You'll need to choose a username for the site, which only take a couple of moments. After that, you can post your question and our members will help you out.

Ask a Question

Members online

Forum statistics

Threads
473,731
Messages
2,569,432
Members
44,832
Latest member
GlennSmall

Latest Threads

Top