Floating Point Equality and Tolerances

Floating point equality and tolerances

This blogpost contains an example, fairly foolproof implementation, and detailed theory behind it
http://randomascii.wordpress.com/2012/02/25/comparing-floating-point-numbers-2012-edition/
it is also one of a series, so you can always read more.
In short: use ULP for most numbers, use epsilon for numbers near zero, but there are still caveats. If you want to be sure about your floating point math i recommend reading whole series.

Floating point comparison without any tolerance

Please note that there is no such thing as exactness in floating points due to the fundamental limitation of the binary code. You cannot express decimal in binary exactly, instead halves, quarters (and other ratios with denominator being a power of two) are used to express values that are 'close enough' within some acceptable tolerance. You cannot remove this tolerance and compare exact values, because otherwise you cannot express decimals in binary.

Useful link.

UPDATE:
Just found on the stackoverflow this great response, so that I don't have to write it all over. Read this and maybe you will consider using BigDecimal for the task at hand.

How should I do floating point comparison?

Comparing for greater/smaller is not really a problem unless you're working right at the edge of the float/double precision limit.

For a "fuzzy equals" comparison, this (Java code, should be easy to adapt) is what I came up with for The Floating-Point Guide after a lot of work and taking into account lots of criticism:

public static boolean nearlyEqual(float a, float b, float epsilon) {
    final float absA = Math.abs(a);
    final float absB = Math.abs(b);
    final float diff = Math.abs(a - b);

    if (a == b) { // shortcut, handles infinities
        return true;
    } else if (a == 0 || b == 0 || diff < Float.MIN_NORMAL) {
        // a or b is zero or both are extremely close to it
        // relative error is less meaningful here
        return diff < (epsilon * Float.MIN_NORMAL);
    } else { // use relative error
        return diff / (absA + absB) < epsilon;
    }
}

It comes with a test suite. You should immediately dismiss any solution that doesn't, because it is virtually guaranteed to fail in some edge cases like having one value 0, two very small values opposite of zero, or infinities.

An alternative (see link above for more details) is to convert the floats' bit patterns to integer and accept everything within a fixed integer distance.

In any case, there probably isn't any solution that is perfect for all applications. Ideally, you'd develop/adapt your own with a test suite covering your actual use cases.

What is the best way to compare floats for almost-equality in Python?

Python 3.5 adds the math.isclose and cmath.isclose functions as described in PEP 485.

If you're using an earlier version of Python, the equivalent function is given in the documentation.

def isclose(a, b, rel_tol=1e-09, abs_tol=0.0):
    return abs(a-b) <= max(rel_tol * max(abs(a), abs(b)), abs_tol)

rel_tol is a relative tolerance, it is multiplied by the greater of the magnitudes of the two arguments; as the values get larger, so does the allowed difference between them while still considering them equal.

abs_tol is an absolute tolerance that is applied as-is in all cases. If the difference is less than either of those tolerances, the values are considered equal.

Floating point Arithmetics

1. Efficient way to compare two floating values.

A simple double a,b; if (a == b) is an efficient way to compare two floating values. Yet as OP noticed, this may not meet the overall coding goal. Better ways depend on the context of the compare, something not supplied by OP. See far below.

2. How to add a floating value to another one. Example. Add 0.1111 to 94.4345 to get the exact value as 94.5456

Floating values as source code have effective unlimited range and precision such as 1.23456789012345678901234567890e1234567. Conversion of this text to a double is limited typically to one of 2⁶⁴ different values. The closest is selected, but that may not be an exact match.

Neither 0.1111, 94.4345, 94.5456 can be representably exactly as a typical double.

OP has choices:

1.) Use another type other than double, float. Various libraries offer decimal floating point types.

2) Limit code to rare platforms that support double to a base 10 form such that FLT_RADIX == 10.

3) Write your own code to handle user input like "0.1111" into a structure/string and perform the needed operations.

4) Treat user input as strings and the convert to some integer type, again with supported routines to read/compute/and write.

5) Accept that floating point operations are not mathematically exact and handle round-off error.

double a = 0.1111;
printf("a:   %.*e\n", DBL_DECIMAL_DIG -1 , a);
double b = 94.4345;
printf("b:   %.*e\n", DBL_DECIMAL_DIG -1 , b);
double sum = a + b;
printf("sum: %.*e\n", DBL_DECIMAL_DIG -1 , sum);
printf("%.4f\n", sum);

Output

a:   1.1110000000000000e-01
b:   9.4434500000000000e+01
sum: 9.4545599999999993e+01
94.5456  // Desired textual output based on a rounded `sum` to the nearest 0.0001

---

More on #1

If an exact compare is not sought but some sort of "are the two values close enough?", a definition of "close enough" is needed - of which there are many.

The following "close enough" compares the distance by examining the ULP of the two numbers. It is a linear difference when the values are in the same power-of-two and becomes logarithmic other wise. Of course, change of sign is an issue.

float example:

Consider all finite float ordered from most negative to most positive. The following, somewhat-portable code, returns an integer for each float with that same order.

uint32_t sequence_f(float x) {
  union {
    float f;
    uint32_t u32;
  } u;
  assert(sizeof(float) == sizeof(uint32_t));
  u.f = x;
  if (u.u32 & 0x80000000) {
    u.u32 ^= 0x80000000;
    return 0x80000000 - u.u32;
  }
  return u.u3
}

Now, to determine if two float are "close enough", simple compare two integers.

static bool close_enough(float x, float y, uint32_t ULP_delta) {
  uint32_t ullx = sequence_f(x);
  uint32_t ully = sequence_f(y);
  if (ullx > ully) return (ullx - ully) <= ULP_delta;
  return (ully - ullx) <= ULP_delta;
}

Compare floats in R

.Machine$double.eps is the difference between 1 and the smallest representable value greater than 1. The difference between 0.1 and the smallest representable value greater than 0.1 is smaller than .Machine$double.eps and the difference between 100 and the smallest representable value greater than 100 is larger than .Machine$double.eps. Have a look at: What is the correct/standard way to check if difference is smaller than machine precision?.

.Machine$double.eps is

.Machine$double.eps
[1] 2.220446e-16

When you make the calculations the intern sotred values will be approximately:

print(1.0 + .Machine$double.eps, 20)
#[1] 1.000000000000000222
print(1.0 - .Machine$double.eps, 20)
#[1] 0.99999999999999977796
print(0.9 + .Machine$double.eps, 20)
#[1] 0.90000000000000024425
print(2.0 + .Machine$double.eps, 20)
#[1] 2

Using tolerance = .Machine$double.eps all.equal returns TRUE or FALSE depending if the difference of the intern stored values is lager or not than tolerance.

To compare 2 numbers in R if the are intern stored equal use ==.

---

Related Topics