float_pow: Put *all* of the burden on the libm pow in normal
cases.
powu: Deleted.
This started with a nonsensical error msg:
>>> x = -1.
>>> import sys
>>> x**(-sys.maxint-1L)
Traceback (most recent call last):
File "<stdin>", line 1, in ?
ValueError: negative number cannot be raised to a fractional power
>>>
The special-casing in float_pow was simply wrong in this case (there's
not even anything peculiar about these inputs), and I don't see any point
to it in *any* case: a decent libm pow should have worst-case error under
1 ULP, so in particular should deliver the exact result whenever the exact
result is representable (else its error is at least 1 ULP). Thus our
special fiddling for integral values "shouldn't" buy anything in accuracy,
and, to the contrary, repeated multiplication is less accurate than a
decent pow when the true result isn't exactly representable. So just
letting pow() do its job here (we may not be able to trust libm x-platform
in exceptional cases, but these are normal cases).
diff --git a/Objects/floatobject.c b/Objects/floatobject.c
index 044d1d3..69aede4 100644
--- a/Objects/floatobject.c
+++ b/Objects/floatobject.c
@@ -469,25 +469,10 @@
return Py_BuildValue("(dd)", floordiv, mod);
}
-static double powu(double x, long n)
-{
- double r = 1.;
- double p = x;
- long mask = 1;
- while (mask > 0 && n >= mask) {
- if (n & mask)
- r *= p;
- mask <<= 1;
- p *= p;
- }
- return r;
-}
-
static PyObject *
float_pow(PyObject *v, PyObject *w, PyObject *z)
{
double iv, iw, ix;
- long intw;
/* XXX Doesn't handle overflows if z!=None yet; it may never do so :(
* The z parameter is really only going to be useful for integers and
* long integers. Maybe something clever with logarithms could be done.
@@ -495,23 +480,23 @@
*/
CONVERT_TO_DOUBLE(v, iv);
CONVERT_TO_DOUBLE(w, iw);
- intw = (long)iw;
/* Sort out special cases here instead of relying on pow() */
- if (iw == 0) { /* x**0 is 1, even 0**0 */
+ if (iw == 0) { /* v**0 is 1, even 0**0 */
PyFPE_START_PROTECT("pow", return NULL)
if ((PyObject *)z != Py_None) {
double iz;
CONVERT_TO_DOUBLE(z, iz);
- ix=fmod(1.0, iz);
- if (ix!=0 && iz<0) ix+=iz;
+ ix = fmod(1.0, iz);
+ if (ix != 0 && iz < 0)
+ ix += iz;
}
else
ix = 1.0;
PyFPE_END_PROTECT(ix)
return PyFloat_FromDouble(ix);
}
- if (iv == 0.0) {
+ if (iv == 0.0) { /* 0**w is error if w<0, else 1 */
if (iw < 0.0) {
PyErr_SetString(PyExc_ZeroDivisionError,
"0.0 cannot be raised to a negative power");
@@ -519,29 +504,15 @@
}
return PyFloat_FromDouble(0.0);
}
-
- if (iw == intw && intw > LONG_MIN) {
- /* ruled out LONG_MIN because -LONG_MIN isn't representable */
- errno = 0;
- PyFPE_START_PROTECT("pow", return NULL)
- if (intw > 0)
- ix = powu(iv, intw);
- else
- ix = 1./powu(iv, -intw);
- PyFPE_END_PROTECT(ix)
+ if (iv < 0.0 && iw != floor(iw)) {
+ PyErr_SetString(PyExc_ValueError,
+ "negative number cannot be raised to a fractional power");
+ return NULL;
}
- else {
- /* Sort out special cases here instead of relying on pow() */
- if (iv < 0.0) {
- PyErr_SetString(PyExc_ValueError,
- "negative number cannot be raised to a fractional power");
- return NULL;
- }
- errno = 0;
- PyFPE_START_PROTECT("pow", return NULL)
- ix = pow(iv, iw);
- PyFPE_END_PROTECT(ix)
- }
+ errno = 0;
+ PyFPE_START_PROTECT("pow", return NULL)
+ ix = pow(iv, iw);
+ PyFPE_END_PROTECT(ix)
CHECK(ix);
if (errno != 0) {
/* XXX could it be another type of error? */
@@ -552,9 +523,9 @@
double iz;
CONVERT_TO_DOUBLE(z, iz);
PyFPE_START_PROTECT("pow", return 0)
- ix=fmod(ix, iz); /* XXX To Be Rewritten */
- if (ix!=0 && ((iv<0 && iz>0) || (iv>0 && iz<0) )) {
- ix+=iz;
+ ix = fmod(ix, iz); /* XXX To Be Rewritten */
+ if (ix != 0 && ((iv < 0 && iz > 0) || (iv > 0 && iz < 0) )) {
+ ix += iz;
}
PyFPE_END_PROTECT(ix)
}