diff -r dca18f0ec280 Lib/statistics.py --- a/Lib/statistics.py Sat Oct 01 03:11:04 2016 +0000 +++ b/Lib/statistics.py Sat Oct 01 13:14:00 2016 +0100 @@ -478,6 +478,55 @@ assert type(_nth_root) is type(lambda: None) +def _frexp_gen(n): + """Version of math.frexp that works for floats and arbitrary integers. + + This version of frexp avoids the OverflowError that occurs for math.frexp + applied to integers outside the range of a float. + + Examples + -------- + >>> _frexp_gen(10) + (0.625, 4) + >>> _frexp_gen(16) + (0.5, 5) + >>> _frexp_gen(-10) + (-0.625, 4) + >>> _frexp_gen(2**1023) + (0.5, 1024) + >>> _frexp_gen(2**1024) + (0.5, 1025) + >>> _frexp_gen(2**1024 - 1) + (0.5, 1025) + >>> _frexp_gen(2**1024 - 2**970) # round-ties-to-even check + (0.5, 1025) + >>> _frexp_gen(2**1024 + 2**971) # round-ties-to-even check + (0.5, 1025) + >>> _frexp_gen(2**1024 + 3 * 2**971) # round-ties-to-even check + (0.5000000000000002, 1025) + >>> _frexp_gen(2**1024 + 5 * 2**971) # round-ties-to-even check + (0.5000000000000002, 1025) + >>> _frexp_gen(2**1024 + 7 * 2**971) # round-ties-to-even check + (0.5000000000000004, 1025) + >>> _frexp_gen(0) + (0.0, 0) + >>> _frexp_gen(-2**1024) + (-0.5, 1025) + """ + try: + return math.frexp(n) + except OverflowError: + # Does correct rounding assuming IEEE 754 binary64 format, but + # shouldn't give badly wrong results for other formats. For abs(n) in + # the range [2**53, 2**1024 - 2**970), the code below should give + # identical results to math.frexp(n). For the purposes of _product, + # the abs(m) == 1.0 check can be omitted - it's safe to simply return + # (m, e). + e = n.bit_length() + m = ((n >> e - 54) - (-n >> e - 54)) / (1 << 55) + return (m / 2.0, e + 1) if abs(m) == 1.0 else (m, e) + + def _product(values): """Return product of values as (exponent, mantissa).""" errmsg = 'mixed Decimal and float is not supported' @@ -499,11 +548,11 @@ # # x1*x2 = 2**p1*m1 * 2**p2*m2 = 2**(p1+p2)*(m1*m2) # - mant, scale = 1, 0 #math.frexp(prod) # FIXME + mant, scale = _frexp_gen(prod) for y in chain([x], values): if isinstance(y, Decimal): raise TypeError(errmsg) - m1, e1 = math.frexp(y) + m1, e1 = _frexp_gen(y) m2, e2 = math.frexp(mant) scale += (e1 + e2) mant = m1*m2