Home | History | Annotate | Download | only in test
      1 from __future__ import division
      2 # When true division is the default, get rid of this and add it to
      3 # test_long.py instead.  In the meantime, it's too obscure to try to
      4 # trick just part of test_long into using future division.
      5 
      6 import sys
      7 import random
      8 import math
      9 import unittest
     10 from test.test_support import run_unittest
     11 
     12 # decorator for skipping tests on non-IEEE 754 platforms
     13 requires_IEEE_754 = unittest.skipUnless(
     14     float.__getformat__("double").startswith("IEEE"),
     15     "test requires IEEE 754 doubles")
     16 
     17 DBL_MAX = sys.float_info.max
     18 DBL_MAX_EXP = sys.float_info.max_exp
     19 DBL_MIN_EXP = sys.float_info.min_exp
     20 DBL_MANT_DIG = sys.float_info.mant_dig
     21 DBL_MIN_OVERFLOW = 2**DBL_MAX_EXP - 2**(DBL_MAX_EXP - DBL_MANT_DIG - 1)
     22 
     23 # pure Python version of correctly-rounded true division
     24 def truediv(a, b):
     25     """Correctly-rounded true division for integers."""
     26     negative = a^b < 0
     27     a, b = abs(a), abs(b)
     28 
     29     # exceptions:  division by zero, overflow
     30     if not b:
     31         raise ZeroDivisionError("division by zero")
     32     if a >= DBL_MIN_OVERFLOW * b:
     33         raise OverflowError("int/int too large to represent as a float")
     34 
     35    # find integer d satisfying 2**(d - 1) <= a/b < 2**d
     36     d = a.bit_length() - b.bit_length()
     37     if d >= 0 and a >= 2**d * b or d < 0 and a * 2**-d >= b:
     38         d += 1
     39 
     40     # compute 2**-exp * a / b for suitable exp
     41     exp = max(d, DBL_MIN_EXP) - DBL_MANT_DIG
     42     a, b = a << max(-exp, 0), b << max(exp, 0)
     43     q, r = divmod(a, b)
     44 
     45     # round-half-to-even: fractional part is r/b, which is > 0.5 iff
     46     # 2*r > b, and == 0.5 iff 2*r == b.
     47     if 2*r > b or 2*r == b and q % 2 == 1:
     48         q += 1
     49 
     50     result = math.ldexp(float(q), exp)
     51     return -result if negative else result
     52 
     53 class TrueDivisionTests(unittest.TestCase):
     54     def test(self):
     55         huge = 1L << 40000
     56         mhuge = -huge
     57         self.assertEqual(huge / huge, 1.0)
     58         self.assertEqual(mhuge / mhuge, 1.0)
     59         self.assertEqual(huge / mhuge, -1.0)
     60         self.assertEqual(mhuge / huge, -1.0)
     61         self.assertEqual(1 / huge, 0.0)
     62         self.assertEqual(1L / huge, 0.0)
     63         self.assertEqual(1 / mhuge, 0.0)
     64         self.assertEqual(1L / mhuge, 0.0)
     65         self.assertEqual((666 * huge + (huge >> 1)) / huge, 666.5)
     66         self.assertEqual((666 * mhuge + (mhuge >> 1)) / mhuge, 666.5)
     67         self.assertEqual((666 * huge + (huge >> 1)) / mhuge, -666.5)
     68         self.assertEqual((666 * mhuge + (mhuge >> 1)) / huge, -666.5)
     69         self.assertEqual(huge / (huge << 1), 0.5)
     70         self.assertEqual((1000000 * huge) / huge, 1000000)
     71 
     72         namespace = {'huge': huge, 'mhuge': mhuge}
     73 
     74         for overflow in ["float(huge)", "float(mhuge)",
     75                          "huge / 1", "huge / 2L", "huge / -1", "huge / -2L",
     76                          "mhuge / 100", "mhuge / 100L"]:
     77             # If the "eval" does not happen in this module,
     78             # true division is not enabled
     79             with self.assertRaises(OverflowError):
     80                 eval(overflow, namespace)
     81 
     82         for underflow in ["1 / huge", "2L / huge", "-1 / huge", "-2L / huge",
     83                          "100 / mhuge", "100L / mhuge"]:
     84             result = eval(underflow, namespace)
     85             self.assertEqual(result, 0.0, 'expected underflow to 0 '
     86                              'from {!r}'.format(underflow))
     87 
     88         for zero in ["huge / 0", "huge / 0L", "mhuge / 0", "mhuge / 0L"]:
     89             with self.assertRaises(ZeroDivisionError):
     90                 eval(zero, namespace)
     91 
     92     def check_truediv(self, a, b, skip_small=True):
     93         """Verify that the result of a/b is correctly rounded, by
     94         comparing it with a pure Python implementation of correctly
     95         rounded division.  b should be nonzero."""
     96 
     97         a, b = long(a), long(b)
     98 
     99         # skip check for small a and b: in this case, the current
    100         # implementation converts the arguments to float directly and
    101         # then applies a float division.  This can give doubly-rounded
    102         # results on x87-using machines (particularly 32-bit Linux).
    103         if skip_small and max(abs(a), abs(b)) < 2**DBL_MANT_DIG:
    104             return
    105 
    106         try:
    107             # use repr so that we can distinguish between -0.0 and 0.0
    108             expected = repr(truediv(a, b))
    109         except OverflowError:
    110             expected = 'overflow'
    111         except ZeroDivisionError:
    112             expected = 'zerodivision'
    113 
    114         try:
    115             got = repr(a / b)
    116         except OverflowError:
    117             got = 'overflow'
    118         except ZeroDivisionError:
    119             got = 'zerodivision'
    120 
    121         self.assertEqual(expected, got, "Incorrectly rounded division {}/{}: "
    122                          "expected {}, got {}".format(a, b, expected, got))
    123 
    124     @requires_IEEE_754
    125     def test_correctly_rounded_true_division(self):
    126         # more stringent tests than those above, checking that the
    127         # result of true division of ints is always correctly rounded.
    128         # This test should probably be considered CPython-specific.
    129 
    130         # Exercise all the code paths not involving Gb-sized ints.
    131         # ... divisions involving zero
    132         self.check_truediv(123, 0)
    133         self.check_truediv(-456, 0)
    134         self.check_truediv(0, 3)
    135         self.check_truediv(0, -3)
    136         self.check_truediv(0, 0)
    137         # ... overflow or underflow by large margin
    138         self.check_truediv(671 * 12345 * 2**DBL_MAX_EXP, 12345)
    139         self.check_truediv(12345, 345678 * 2**(DBL_MANT_DIG - DBL_MIN_EXP))
    140         # ... a much larger or smaller than b
    141         self.check_truediv(12345*2**100, 98765)
    142         self.check_truediv(12345*2**30, 98765*7**81)
    143         # ... a / b near a boundary: one of 1, 2**DBL_MANT_DIG, 2**DBL_MIN_EXP,
    144         #                 2**DBL_MAX_EXP, 2**(DBL_MIN_EXP-DBL_MANT_DIG)
    145         bases = (0, DBL_MANT_DIG, DBL_MIN_EXP,
    146                  DBL_MAX_EXP, DBL_MIN_EXP - DBL_MANT_DIG)
    147         for base in bases:
    148             for exp in range(base - 15, base + 15):
    149                 self.check_truediv(75312*2**max(exp, 0), 69187*2**max(-exp, 0))
    150                 self.check_truediv(69187*2**max(exp, 0), 75312*2**max(-exp, 0))
    151 
    152         # overflow corner case
    153         for m in [1, 2, 7, 17, 12345, 7**100,
    154                   -1, -2, -5, -23, -67891, -41**50]:
    155             for n in range(-10, 10):
    156                 self.check_truediv(m*DBL_MIN_OVERFLOW + n, m)
    157                 self.check_truediv(m*DBL_MIN_OVERFLOW + n, -m)
    158 
    159         # check detection of inexactness in shifting stage
    160         for n in range(250):
    161             # (2**DBL_MANT_DIG+1)/(2**DBL_MANT_DIG) lies halfway
    162             # between two representable floats, and would usually be
    163             # rounded down under round-half-to-even.  The tiniest of
    164             # additions to the numerator should cause it to be rounded
    165             # up instead.
    166             self.check_truediv((2**DBL_MANT_DIG + 1)*12345*2**200 + 2**n,
    167                            2**DBL_MANT_DIG*12345)
    168 
    169         # 1/2731 is one of the smallest division cases that's subject
    170         # to double rounding on IEEE 754 machines working internally with
    171         # 64-bit precision.  On such machines, the next check would fail,
    172         # were it not explicitly skipped in check_truediv.
    173         self.check_truediv(1, 2731)
    174 
    175         # a particularly bad case for the old algorithm:  gives an
    176         # error of close to 3.5 ulps.
    177         self.check_truediv(295147931372582273023, 295147932265116303360)
    178         for i in range(1000):
    179             self.check_truediv(10**(i+1), 10**i)
    180             self.check_truediv(10**i, 10**(i+1))
    181 
    182         # test round-half-to-even behaviour, normal result
    183         for m in [1, 2, 4, 7, 8, 16, 17, 32, 12345, 7**100,
    184                   -1, -2, -5, -23, -67891, -41**50]:
    185             for n in range(-10, 10):
    186                 self.check_truediv(2**DBL_MANT_DIG*m + n, m)
    187 
    188         # test round-half-to-even, subnormal result
    189         for n in range(-20, 20):
    190             self.check_truediv(n, 2**1076)
    191 
    192         # largeish random divisions: a/b where |a| <= |b| <=
    193         # 2*|a|; |ans| is between 0.5 and 1.0, so error should
    194         # always be bounded by 2**-54 with equality possible only
    195         # if the least significant bit of q=ans*2**53 is zero.
    196         for M in [10**10, 10**100, 10**1000]:
    197             for i in range(1000):
    198                 a = random.randrange(1, M)
    199                 b = random.randrange(a, 2*a+1)
    200                 self.check_truediv(a, b)
    201                 self.check_truediv(-a, b)
    202                 self.check_truediv(a, -b)
    203                 self.check_truediv(-a, -b)
    204 
    205         # and some (genuinely) random tests
    206         for _ in range(10000):
    207             a_bits = random.randrange(1000)
    208             b_bits = random.randrange(1, 1000)
    209             x = random.randrange(2**a_bits)
    210             y = random.randrange(1, 2**b_bits)
    211             self.check_truediv(x, y)
    212             self.check_truediv(x, -y)
    213             self.check_truediv(-x, y)
    214             self.check_truediv(-x, -y)
    215 
    216 
    217 def test_main():
    218     run_unittest(TrueDivisionTests)
    219 
    220 if __name__ == "__main__":
    221     test_main()
    222