diff options
Diffstat (limited to 'Modules/mathmodule.c')
-rw-r--r-- | Modules/mathmodule.c | 87 |
1 files changed, 74 insertions, 13 deletions
diff --git a/Modules/mathmodule.c b/Modules/mathmodule.c index 2c4cc73..d693da1 100644 --- a/Modules/mathmodule.c +++ b/Modules/mathmodule.c @@ -55,11 +55,6 @@ raised for division by zero and mod by zero. #include "Python.h" #include "_math.h" -#ifdef _OSF_SOURCE -/* OSF1 5.1 doesn't make this available with XOPEN_SOURCE_EXTENDED defined */ -extern double copysign(double, double); -#endif - /* sin(pi*x), giving accurate results for all finite x (especially x integral or close to an integer). This is here for use in the @@ -244,7 +239,8 @@ m_tgamma(double x) } if (x == 0.0) { errno = EDOM; - return 1.0/x; /* tgamma(+-0.0) = +-inf, divide-by-zero */ + /* tgamma(+-0.0) = +-inf, divide-by-zero */ + return copysign(Py_HUGE_VAL, x); } /* integer arguments */ @@ -582,6 +578,61 @@ m_log(double x) } } +/* + log2: log to base 2. + + Uses an algorithm that should: + + (a) produce exact results for powers of 2, and + (b) give a monotonic log2 (for positive finite floats), + assuming that the system log is monotonic. +*/ + +static double +m_log2(double x) +{ + if (!Py_IS_FINITE(x)) { + if (Py_IS_NAN(x)) + return x; /* log2(nan) = nan */ + else if (x > 0.0) + return x; /* log2(+inf) = +inf */ + else { + errno = EDOM; + return Py_NAN; /* log2(-inf) = nan, invalid-operation */ + } + } + + if (x > 0.0) { +#ifdef HAVE_LOG2 + return log2(x); +#else + double m; + int e; + m = frexp(x, &e); + /* We want log2(m * 2**e) == log(m) / log(2) + e. Care is needed when + * x is just greater than 1.0: in that case e is 1, log(m) is negative, + * and we get significant cancellation error from the addition of + * log(m) / log(2) to e. The slight rewrite of the expression below + * avoids this problem. + */ + if (x >= 1.0) { + return log(2.0 * m) / log(2.0) + (e - 1); + } + else { + return log(m) / log(2.0) + e; + } +#endif + } + else if (x == 0.0) { + errno = EDOM; + return -Py_HUGE_VAL; /* log2(0) = -inf, divide-by-zero */ + } + else { + errno = EDOM; + return Py_NAN; /* log2(-inf) = nan, invalid-operation */ + } +} + static double m_log10(double x) { @@ -836,10 +887,10 @@ FUNC1(atanh, m_atanh, 0, "atanh(x)\n\nReturn the hyperbolic arc tangent (measured in radians) of x.") static PyObject * math_ceil(PyObject *self, PyObject *number) { - static PyObject *ceil_str = NULL; + _Py_IDENTIFIER(__ceil__); PyObject *method, *result; - method = _PyObject_LookupSpecial(number, "__ceil__", &ceil_str); + method = _PyObject_LookupSpecial(number, &PyId___ceil__); if (method == NULL) { if (PyErr_Occurred()) return NULL; @@ -874,10 +925,10 @@ FUNC1(fabs, fabs, 0, "fabs(x)\n\nReturn the absolute value of the float x.") static PyObject * math_floor(PyObject *self, PyObject *number) { - static PyObject *floor_str = NULL; + _Py_IDENTIFIER(__floor__); PyObject *method, *result; - method = _PyObject_LookupSpecial(number, "__floor__", &floor_str); + method = _PyObject_LookupSpecial(number, &PyId___floor__); if (method == NULL) { if (PyErr_Occurred()) return NULL; @@ -1384,7 +1435,7 @@ math_factorial(PyObject *self, PyObject *arg) } /* use lookup table if x is small */ - if (x < (long)(sizeof(SmallFactorials)/sizeof(SmallFactorials[0]))) + if (x < (long)Py_ARRAY_LENGTH(SmallFactorials)) return PyLong_FromUnsignedLong(SmallFactorials[x]); /* else express in the form odd_part * 2**two_valuation, and compute as @@ -1411,7 +1462,7 @@ PyDoc_STRVAR(math_factorial_doc, static PyObject * math_trunc(PyObject *self, PyObject *number) { - static PyObject *trunc_str = NULL; + _Py_IDENTIFIER(__trunc__); PyObject *trunc, *result; if (Py_TYPE(number)->tp_dict == NULL) { @@ -1419,7 +1470,7 @@ math_trunc(PyObject *self, PyObject *number) return NULL; } - trunc = _PyObject_LookupSpecial(number, "__trunc__", &trunc_str); + trunc = _PyObject_LookupSpecial(number, &PyId___trunc__); if (trunc == NULL) { if (!PyErr_Occurred()) PyErr_Format(PyExc_TypeError, @@ -1628,6 +1679,15 @@ Return the logarithm of x to the given base.\n\ If the base not specified, returns the natural logarithm (base e) of x."); static PyObject * +math_log2(PyObject *self, PyObject *arg) +{ + return loghelper(arg, m_log2, "log2"); +} + +PyDoc_STRVAR(math_log2_doc, +"log2(x)\n\nReturn the base 2 logarithm of x."); + +static PyObject * math_log10(PyObject *self, PyObject *arg) { return loghelper(arg, m_log10, "log10"); @@ -1899,6 +1959,7 @@ static PyMethodDef math_methods[] = { {"log", math_log, METH_VARARGS, math_log_doc}, {"log1p", math_log1p, METH_O, math_log1p_doc}, {"log10", math_log10, METH_O, math_log10_doc}, + {"log2", math_log2, METH_O, math_log2_doc}, {"modf", math_modf, METH_O, math_modf_doc}, {"pow", math_pow, METH_VARARGS, math_pow_doc}, {"radians", math_radians, METH_O, math_radians_doc}, |