racket
racket copied to clipboard
Inaccurate (expt flonum exact-int)
What version of Racket are you using? Welcome to Racket v8.16.0.2 [cs].
What program did you run? Please include a short example that triggers the bug
Welcome to Racket v8.16.0.2 [cs].
> (expt 1.000000000000001 (expt 2 56))
5.540619075645279e+34
> (expt -1.000000000000001 (expt 2 56))
5.540619075645279e+34
> (expt 1.000000000000001 (+ 1 (expt 2 56)))
5.5406190756452855e+34
> (expt -1.000000000000001 (+ 1 (expt 2 56)))
-5.5406190756452855e+34
What should have happened?
Roughly (from Gambit):
> (expt 1.000000000000001 (expt 2 56))
5.540622384393264e34
> (expt -1.000000000000001 (expt 2 56))
5.540622384393264e34
> (expt 1.000000000000001 (+ 1 (expt 2 56)))
5.540622384393264e34
> (expt -1.000000000000001 (+ 1 (expt 2 56)))
-5.540622384393264e34
If you got an error message, please include it here.
Please include any other relevant details e.g., the operating system used or how you are running the code.
The system C library power function should compute this very accurately, (log 1.000000000000001) is easily computed to machine accuracy, multiplication by (expt 2 56) is exact, and then calculating exp is also very accurate:
> (exp (* (log 1.000000000000001) (expt 2 56)))
5.540622384393274e34
> (flexpt 1.000000000000001 (inexact (expt 2 56)))
5.540622384393264e34
The iterative method used in Racket's library is very inaccurate for this example (getting only five correct digits in the result). One can compute the correctly rounded results (here with my computable reals library, but also with MPFR, etc.):
> (computable->inexact (computable-expt (->computable 1.000000000000001) (expt 2 56)))
5.540622384393264e34
> (computable->inexact (computable-expt (->computable 1.000000000000001) (+ 1 (expt 2 56))))
5.54062238439327e34
For (expt flonum exact-int) I would recommend computing (flonum-expt (flabs flonum) (->inexact exact-int)) and then changing the sign of the result if flonum is negative and exact-int is odd.