pythonfloating-pointprecisionpari

the pow() function isn't accurate

I read Calculating very large exponents in python , and Ignacio Vazquez-Abrams's answer was to use the pow() function

So , I run the following command in Python 3.9.6 :

print(int(pow(1.5,96)))

and gives me 80308380747696832

Now trying the same using PARI/GP , I get :

? (1.5)^96
%3 = 80308380747696837.016746608902131688101
?

So no one should rely on pow anymore ?

Solution

  • Python floats are limited by IEEE 754 logic, but you can get an exact Answer by converting to integer math (which has effectively unlimited precision in Python and while not always possible, can be used for this case) or with SymPy, which does symbolic math and can also use mpmath (or gmpy) for more precise evaluation

    >>> a = 15
>>> b = 10
>>> a**96 // b**96  # integer space equivalent distributing (1.5*10/10)**96
80308380747696837

    https://docs.sympy.org/latest/explanation/gotchas.html#evaluating-expressions-with-floats-and-rationals

    >>> from sympy import *
>>> Rational("1.5").evalf()
1.50000000000000
>>> Rational("1.5").evalf(16)
1.500000000000000
>>> Rational("1.5").evalf(50)
1.5000000000000000000000000000000000000000000000000

    >>> pow(Rational("1.5"), 96)   # exact unevaluated fraction from symbolics
6362685441135942358474828762538534230890216321/79228162514264337593543950336
>>> int(pow(Rational("1.5"), 96).evalf())     # default precision
80308380747696832
>>> int(pow(Rational("1.5"), 96).evalf(100))  # increased precision
80308380747696837

    Adjusting the precision changes the lower bits through several values (but note that these are extremely small relative to the size of the number)

    >>> int(pow(Rational("1.5"), 96).evalf(15))  # default
80308380747696832
>>> int(pow(Rational("1.5"), 96).evalf(16))
80308380747696836
>>> int(pow(Rational("1.5"), 96).evalf(17))
80308380747696837

    PARI uses higher-precision math by default (very likely via GMP, though it's an optional build dependency) and likely doing additional work internally to keep values in a good range while evaluating them rather than following IEE 754
    (note, I'd originally never heard of "pari" and read "perl", but all floating-point math precision problems fundamentally exist in this space)

    This is generally a case where a "mathematically correct" and a "computed" result can be at odds due to differences in what correct means to each space

    Very large and very small floating-point values should also generally not be compared literally, but compared with tolerance, and care should be taken to avoid compounding errors, preferring a dedicated library to increase precision when required How should I do floating point comparison?

    Further reading