[Python-Dev] Floor division

Tim Peters tim.peters at gmail.com
Sat Jan 20 02:33:23 CET 2007


[Raymond Hettinger]
> I bumped into an oddity today:
>
>     6.0 // 0.001 != math.floor(6.0 / 0.001)
>
> In looking at Objects/floatobject.c, I was surprised to find that
> float_floor_division() is implemented in terms of float_divmod().  Does anyone
> know why it takes such a circuitous path?  I had expected something simpler and
> faster:
>
>      return PyFloat_FromDouble(floor(a/b));

To preserve, so far as is possible with floats, that

(*)    a == (a//b)*b + a%b

In this case the binary double closest to decimal 0.001 is

0.001000000000000000020816681711721685132943093776702880859375

which is a little bit larger than 1/1000.  Therefore the mathematical
value of a/b is a little bit smaller than 6/(1/1000) = 6000, and the
true floor of the mathematically correct result is 5999.

a % b is always computed exactly (mathematical result == machine
result) under Python's definition whenever a and b have the same sign
(under C99's and the `decimal` module's definition it's always exact
regardless of signs), and getting the exact value for a%b implies the
exact value for a//b must also be returned (when possible) in order to
preserve identity (*) above.

IOW, since

>>> 6.0 % 0.001
0.00099999999999987512

it would be inconsistent to return 6000 for 6.0 // 0.001:

>>> 6.0 - 6000 * 0.001   # this isn't close to the value of a%b
0.0
>>> 6.0 - 5999 * 0.001   # but this is close
0.00099999999999944578

Note that two rounding errors occur in computing a - N*b via binary
doubles.  If there were no rounding errors, we'd have

    6 % b ==  6.0 - 5999 * b

exactly where

b =
0.001000000000000000020816681711721685132943093776702880859375

is the value actually stored in the machine for 0.001:

>>> import decimal
>>> decimal.getcontext().prec = 1000 # way more than needed
>>> import math
>>> # Set b to exact decimal value of binary double closest to 0.001.
>>> m, e = math.frexp(0.001)
>>> b = decimal.Decimal(int(m*2**53)) / decimal.Decimal(1 << (53-e))
>>> # Then 6%b is exactly equal to 6 - 5999*b
>>> 6 % b == 6 - 5999*b
True
>>> # Confirm that all decimal calculations were exact.
>>> decimal.getcontext().flags[decimal.Inexact]
False
>>> # Confirm that floor(6/b) is 5999.
>>> int(6/b)
5999
>>> print 6/b
5999.99999999999987509990972966989180234686226...

All that said, (*) doesn't make a lot of sense for floats to begin
with (for that matter, Python's definition of mod alone doesn't make
much sense for floats when the signs differ -- e.g.,

>>> -1 % 1e100
1e+100
>>> decimal.Decimal(-1) % decimal.Decimal("1e100")
Decimal("-1")

).


More information about the Python-Dev mailing list