comparing doubles for equality

Discussion in 'C Programming' started by John Smith, Dec 30, 2006.

  1. John Smith

    John Smith Guest

    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?
    John Smith, Dec 30, 2006
    #1
    1. Advertising

  2. John Smith

    Tim Prince Guest

    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.
    Tim Prince, Dec 30, 2006
    #2
    1. Advertising

  3. In article <I8Clh.533642$1T2.401780@pd7urf2no> John Smith <> writes:
    > 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, +31205924131
    home: bovenover 215, 1025 jn amsterdam, nederland; http://www.cwi.nl/~dik/
    Dik T. Winter, Dec 31, 2006
    #3
  4. John Smith

    Thad Smith Guest

    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
    Thad Smith, Dec 31, 2006
    #4
  5. John Smith

    CBFalconer Guest

    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, Dec 31, 2006
    #5
  6. John Smith

    CBFalconer Guest

    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>
    CBFalconer, Dec 31, 2006
    #6
  7. John Smith

    Thad Smith Guest

    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
    Thad Smith, Dec 31, 2006
    #7
  8. John Smith

    CBFalconer Guest

    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>
    CBFalconer, Dec 31, 2006
    #8
  9. John Smith

    SM Ryan Guest

    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?

    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.
    SM Ryan, Dec 31, 2006
    #9
  10. In article <> 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, +31205924131
    home: bovenover 215, 1025 jn amsterdam, nederland; http://www.cwi.nl/~dik/
    Dik T. Winter, Jan 1, 2007
    #10
  11. John Smith

    Thad Smith Guest

    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
    Thad Smith, Jan 1, 2007
    #11
  12. John Smith

    CBFalconer Guest

    "Dik T. Winter" wrote:
    > 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>
    CBFalconer, Jan 1, 2007
    #12
  13. John Smith

    pete Guest

    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
    pete, Jan 3, 2007
    #13
    1. Advertising

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

It takes just 2 minutes to sign up (and it's free!). Just click the sign up button to choose a username and then you can ask your own questions on the forum.
Similar Threads
  1. {AGUT2} {H}-IWIK
    Replies:
    4
    Views:
    2,839
    Marcelo Pinto
    Sep 12, 2003
  2. nicolas

    Comparing two doubles

    nicolas, Sep 20, 2003, in forum: C++
    Replies:
    3
    Views:
    4,094
    Jacek Dziedzic
    Sep 21, 2003
  3. dan
    Replies:
    1
    Views:
    2,304
    Jack Klein
    Nov 26, 2003
  4. Thomas Kowalski

    Comparing doubles

    Thomas Kowalski, Jul 9, 2007, in forum: C++
    Replies:
    28
    Views:
    1,108
    kwikius
    Jul 11, 2007
  5. vk
    Replies:
    23
    Views:
    1,355
    James Kanze
    Jan 14, 2009
Loading...

Share This Page