Bug #1

GMT_equal_double and GMT_equal_float fail for values close to zero

Added by Florian over 5 years ago. Updated almost 4 years ago.

Status:ClosedStart date:2012-02-23
Priority:NormalDue date:
Assignee:Florian% Done:

100%

Category:-
Target version:-
Affected version:5.x-dev Platform:

Description

Consider this pseudocode example:

float someFloat = 67329.2348f; // arbitrarily chosen
float nextFloat = NearbyFloat(someFloat, 1); // exactly one ULP away from 'someFloat'
float diff = nextFloat - someFloat; // = .0078125000; diff is 1,006,632,960 ULPs away from zero!

Btw, FLT_EPSILON is 872,415,232 ULPs away from zero.


Related issues

Related to Bug #2: Buffer overflow in function GMT_equal_double Closed 2012-02-23

History

#1 Updated by Florian over 5 years ago

  • Status changed from New to Resolved
  • % Done changed from 0 to 100
  • Affected version changed from all to 5.x-dev

Here is another typical catastrophic cancellation:

float degree = 360.0;
float lessDegree = nearbyFloat(degree, -1); // = 359.9999695
float diff = degree - lessDegree; // = 0.0000305
Now here comes the example why this is dangerous.
#define GMT_CONV_LIMIT 1.0e-8
#define GMT_IS_ZERO(x) (fabs (x) < GMT_CONV_LIMIT)
// test if lessDegree ~= degree (very common gmt code):
GMT_IS_ZERO (degree - lessDegree) // returns false even though lessDegree is only 1 ULP away from degree!
This is probably not what we want. Though, I cheated in this example because I used float instead of double, this example shows that such code can lead to unexpected results. In gmt it worked so far because GMT_CONV_LIMIT is very very large in comparison to DBL_EPSILON. GMT_IS_ZERO---if used in such a way---would also fail for large doubles (10 million or so). Here are some new functions that could be used instead:
/* common_math.c */
 * int32_abs                   Like abs but for int32_t
 * int64_abs                   Like abs but for int64_t
 * floatAlmostEqualUlps        Compare 2 floats with ULPs (not safe for values near zero)
 * doubleAlmostEqualUlps       Same as floatAlmostEqualUlps for double
 * floatAlmostEqualUlpsAndAbs  Compare 2 floats with an absolute epsilon check (values near zero), then based on ULPs
 * doubleAlmostEqualUlpsAndAbs Same as floatAlmostEqualUlpsAndAbs for double
 * floatAlmostEqual            Convenience macro for floatAlmostEqualUlps
 * doubleAlmostEqual           Convenience macro for doubleAlmostEqualUlps
 * floatAlmostEqualZero        Convenience macro for floatAlmostEqualUlpsAndAbs
 * doubleAlmostEqualZero       Convenience macro for doubleAlmostEqualUlpsAndAbs
They are mostly inspired by these C++ functions and don't have the overflow bug #2. *AlmostEqualUlpsAndAbs and *AlmostEqualZero also can handle values close to zero, which GMT_equal_double could not. So instead of

GMT_IS_ZERO(lon-360.0)

it is safer to use

doubleAlmostEqual (lon, 360)

or if both parameters might be zero

doubleAlmostEqualZero (lon1, lon2)

This has been fixed in r9839.

#2 Updated by Remko over 5 years ago

  • Status changed from Resolved to Closed

#3 Updated by Florian almost 4 years ago

  • Target version deleted (Candidate for next minor release)

Also available in: Atom PDF