gh-140443: Use fma in loghelper to improve accuracy of log for very large integers (#140469)

* gh-140443:use fma in loghelper to improve accuracy of log for very large integers

Use fused multiply-add in log_helper() for huge ints.

Saving a rounding here is remarkably effective. Across some millions
of randomized test cases with ints up to a billion bits, on Windows
and using log10, the ULP error distribution was dramatically
flattened, and its range was nearly cut in half. In fact, the largest
error Tim saw was under 0.6 ULP.

---------

Co-authored-by: abhi210 <27881020+Abhi210@users.noreply.github.com>
Co-authored-by: Stan Ulbrych <89152624+StanFromIreland@users.noreply.github.com>
This commit is contained in:
Abhishek Tiwari
2025-10-23 22:35:12 +05:30
committed by GitHub
parent 918a9ac9f4
commit f0291c3f2d
3 changed files with 7 additions and 1 deletions

View File

@@ -1921,6 +1921,7 @@ Tim Tisdall
Jason Tishler Jason Tishler
Christian Tismer Christian Tismer
Jim Tittsler Jim Tittsler
Abhishek Tiwari
Frank J. Tobin Frank J. Tobin
James Tocknell James Tocknell
Bennett Todd Bennett Todd

View File

@@ -0,0 +1,5 @@
The logarithm functions (such as :func:`math.log10` and :func:`math.log`) may now produce
slightly different results for extremely large integers that cannot be
converted to floats without overflow. These results are generally more
accurate, with reduced worst-case error and a tighter overall error
distribution.

View File

@@ -2309,7 +2309,7 @@ loghelper(PyObject* arg, double (*func)(double))
assert(e >= 0); assert(e >= 0);
assert(!PyErr_Occurred()); assert(!PyErr_Occurred());
/* Value is ~= x * 2**e, so the log ~= log(x) + log(2) * e. */ /* Value is ~= x * 2**e, so the log ~= log(x) + log(2) * e. */
result = func(x) + func(2.0) * e; result = fma(func(2.0), (double)e, func(x));
} }
else else
/* Successfully converted x to a double. */ /* Successfully converted x to a double. */