CODEWITHSUNDEEP
pythonfloating-point
devtk
devtk

Reputation: 2169

consistency of floating point division

It is well known that floating point values cannot exactly represent every decimal value. Thus, the floating point value for 1/3 is not exactly 1/3. For this reason, it is generally inadvisable to compare floating point values directly.

However, in this application, I am trying to determine if two fractions a / b and c / d are equivalent. If they are, then there exist integers e and f such that a * e = c * f and b * e = d * f. Assume a, b, c, d, e, f are all positive integers which can be represented exactly in a floating point value.

In practice, simply comparing a/b to c/d worked, but is it guaranteed to work? Is there something within Python and/or IEEE-754 that would guarantee such a scheme to work?

Example code (that shows this scheme works for a reasonable number of values):

a = 1.0
b = 3.0
for c in xrange(1, 999):
    assert a / b == (a * c) / (b * c)

If this is not guaranteed, is there a counterexample with values a, b, c, d such that a/b = c/d (in math) but Python fails to compare a / b == c / d? Again, these are all positive integers which Python can exactly represent in its floating point value.

Upvotes: 3

Views: 528

Answers (3)

Eric Postpischil
Eric Postpischil

Reputation: 222273

By specification of division in IEEE-754, if a/b = c/d, then a/b == c/d, provided that the operations in a/b == c/d are evaluated as specified in IEEE-754 using one floating-point format. (Python does not make this guarantee; it inherits the floating-point behavior of whatever platform it is implemented on, and IEEE-754 behavior is not universally guaranteed.) This is because IEEE-754 specifies that the result is the exact mathematical result rounded to the nearest representable value (according to the applicable rounding mode). Since the mathematical results are identical for a/b and c/d, the computed results are identical.

However, a/b == c/d does not imply a/b = c/d. A counterexample is that 1 / (0x1p53-1) and 1 / (0x1p53-2) both produce 1.11022302462515678694266454965700950366517665087069677287701097156968899071216583251953125 • 10−16 (in hexadecimal floating-point, 0x1.0000000000001p-53).

Upvotes: 2

user2357112
user2357112

Reputation: 280182

IEEE 754 division is specified to behave as if you computed the exact result of the division and then rounded the exact result. If a/b and c/d have the same value in exact arithmetic, then since IEEE 754 rounding is consistent and deterministic, a/b and c/d must have the same result in IEEE 754 floating point.

This only holds if a/b and c/d really do have the same value in exact arithmetic, though. Rounding error in a, b, c, or d can throw this off. Also, the converse doesn't hold; if a/b == c/d in floating point, that doesn't mean a/b and c/d are equal in exact arithmetic.

Rather than dealing with floating-point rounding, why not stick to integers?

a*d == b*c

can be done in unrounded integer arithmetic. Alternatively, why not use an exact rational type?

from fractions import Fraction

Fraction(a, b) == Fraction(c, d)

Upvotes: 5

jpp
jpp

Reputation: 164613

Well, I tried a few random numbers and found this. But as @user2357112 points out, this more due to a deficiency of the test than a good example.

Leaving this here to demonstrate the fact.

a = 653543435456556.0
b = 3.0

for c in range(1, 999):
    assert a / b == (a * c) / (b * c), (c, a / b - (a * c) / (b * c))

---------------------------------------------------------------------------
AssertionError                            Traceback (most recent call last)
<ipython-input-128-405615b37738> in <module>()
      2 b = 3.0
      3 for c in range(1, 999):
----> 4     assert a / b == (a * c) / (b * c), (c, a / b - (a * c) / (b * c))
      5 

AssertionError: (57, -0.03125)

Upvotes: 1

Related Questions