Modules/mathmodule.c

   1 /* Math module -- standard C math library functions, pi and e */
   2
   3 #include "Python.h"
   4 #include "longintrepr.h" /* just for SHIFT */
   5
   6 #ifndef _MSC_VER
   7 #ifndef __STDC__
   8 extern double fmod (double, double);
   9 extern double frexp (double, int *);
  10 extern double ldexp (double, int);
  11 extern double modf (double, double *);
  12 #endif /* __STDC__ */
  13 #endif /* _MSC_VER */
  14
  15 /* Call is_error when errno != 0, and where x is the result libm
  16  * returned.  is_error will usually set up an exception and return
  17  * true (1), but may return false (0) without setting up an exception.
  18  */
  19 static int
  20 is_error(double x)
  21 {
  22         int result = 1; /* presumption of guilt */
  23         assert(errno);  /* non-zero errno is a precondition for calling */
  24         if (errno == EDOM)
  25                 PyErr_SetString(PyExc_ValueError, "math domain error");
  26
  27         else if (errno == ERANGE) {
  28                 /* ANSI C generally requires libm functions to set ERANGE
  29                  * on overflow, but also generally *allows* them to set
  30                  * ERANGE on underflow too.  There's no consistency about
  31                  * the latter across platforms.
  32                  * Alas, C99 never requires that errno be set.
  33                  * Here we suppress the underflow errors (libm functions
  34                  * should return a zero on underflow, and +- HUGE_VAL on
  35                  * overflow, so testing the result for zero suffices to
  36                  * distinguish the cases).
  37                  */
  38                 if (x)
  39                         PyErr_SetString(PyExc_OverflowError,
  40                                         "math range error");
  41                 else
  42                         result = 0;
  43         }
  44         else
  45                 /* Unexpected math error */
  46                 PyErr_SetFromErrno(PyExc_ValueError);
  47         return result;
  48 }
  49
  50 static PyObject *
  51 math_1(PyObject *arg, double (*func) (double))
  52 {
  53         double x = PyFloat_AsDouble(arg);
  54         if (x == -1.0 && PyErr_Occurred())
  55                 return NULL;
  56         errno = 0;
  57         PyFPE_START_PROTECT("in math_1", return 0)
  58         x = (*func)(x);
  59         PyFPE_END_PROTECT(x)
  60         Py_SET_ERRNO_ON_MATH_ERROR(x);
  61         if (errno && is_error(x))
  62                 return NULL;
  63         else
  64                 return PyFloat_FromDouble(x);
  65 }
  66
  67 static PyObject *
  68 math_2(PyObject *args, double (*func) (double, double), char *funcname)
  69 {
  70         PyObject *ox, *oy;
  71         double x, y;
  72         if (! PyArg_UnpackTuple(args, funcname, 2, 2, &ox, &oy))
  73                 return NULL;
  74         x = PyFloat_AsDouble(ox);
  75         y = PyFloat_AsDouble(oy);
  76         if ((x == -1.0 || y == -1.0) && PyErr_Occurred())
  77                 return NULL;
  78         errno = 0;
  79         PyFPE_START_PROTECT("in math_2", return 0)
  80         x = (*func)(x, y);
  81         PyFPE_END_PROTECT(x)
  82         Py_SET_ERRNO_ON_MATH_ERROR(x);
  83         if (errno && is_error(x))
  84                 return NULL;
  85         else
  86                 return PyFloat_FromDouble(x);
  87 }
  88
  89 #define FUNC1(funcname, func, docstring) \
  90         static PyObject * math_##funcname(PyObject *self, PyObject *args) { \
  91                 return math_1(args, func); \
  92         }\
  93         PyDoc_STRVAR(math_##funcname##_doc, docstring);
  94
  95 #define FUNC2(funcname, func, docstring) \
  96         static PyObject * math_##funcname(PyObject *self, PyObject *args) { \
  97                 return math_2(args, func, #funcname); \
  98         }\
  99         PyDoc_STRVAR(math_##funcname##_doc, docstring);
 100
 101 FUNC1(acos, acos,
 102       "acos(x)\n\nReturn the arc cosine (measured in radians) of x.")
 103 FUNC1(asin, asin,
 104       "asin(x)\n\nReturn the arc sine (measured in radians) of x.")
 105 FUNC1(atan, atan,
 106       "atan(x)\n\nReturn the arc tangent (measured in radians) of x.")
 107 FUNC2(atan2, atan2,
 108       "atan2(y, x)\n\nReturn the arc tangent (measured in radians) of y/x.\n"
 109       "Unlike atan(y/x), the signs of both x and y are considered.")
 110 FUNC1(ceil, ceil,
 111       "ceil(x)\n\nReturn the ceiling of x as a float.\n"
 112       "This is the smallest integral value >= x.")
 113 FUNC1(cos, cos,
 114       "cos(x)\n\nReturn the cosine of x (measured in radians).")
 115 FUNC1(cosh, cosh,
 116       "cosh(x)\n\nReturn the hyperbolic cosine of x.")
 117 FUNC1(exp, exp,
 118       "exp(x)\n\nReturn e raised to the power of x.")
 119 FUNC1(fabs, fabs,
 120       "fabs(x)\n\nReturn the absolute value of the float x.")
 121 FUNC1(floor, floor,
 122       "floor(x)\n\nReturn the floor of x as a float.\n"
 123       "This is the largest integral value <= x.")
 124 FUNC2(fmod, fmod,
 125       "fmod(x,y)\n\nReturn fmod(x, y), according to platform C."
 126       "  x % y may differ.")
 127 FUNC2(hypot, hypot,
 128       "hypot(x,y)\n\nReturn the Euclidean distance, sqrt(x*x + y*y).")
 129 FUNC2(pow, pow,
 130       "pow(x,y)\n\nReturn x**y (x to the power of y).")
 131 FUNC1(sin, sin,
 132       "sin(x)\n\nReturn the sine of x (measured in radians).")
 133 FUNC1(sinh, sinh,
 134       "sinh(x)\n\nReturn the hyperbolic sine of x.")
 135 FUNC1(sqrt, sqrt,
 136       "sqrt(x)\n\nReturn the square root of x.")
 137 FUNC1(tan, tan,
 138       "tan(x)\n\nReturn the tangent of x (measured in radians).")
 139 FUNC1(tanh, tanh,
 140       "tanh(x)\n\nReturn the hyperbolic tangent of x.")
 141
 142 static PyObject *
 143 math_frexp(PyObject *self, PyObject *arg)
 144 {
 145         int i;
 146         double x = PyFloat_AsDouble(arg);
 147         if (x == -1.0 && PyErr_Occurred())
 148                 return NULL;
 149         errno = 0;
 150         x = frexp(x, &i);
 151         Py_SET_ERRNO_ON_MATH_ERROR(x);
 152         if (errno && is_error(x))
 153                 return NULL;
 154         else
 155                 return Py_BuildValue("(di)", x, i);
 156 }
 157
 158 PyDoc_STRVAR(math_frexp_doc,
 159 "frexp(x)\n"
 160 "\n"
 161 "Return the mantissa and exponent of x, as pair (m, e).\n"
 162 "m is a float and e is an int, such that x = m * 2.**e.\n"
 163 "If x is 0, m and e are both 0.  Else 0.5 <= abs(m) < 1.0.");
 164
 165 static PyObject *
 166 math_ldexp(PyObject *self, PyObject *args)
 167 {
 168         double x;
 169         int exp;
 170         if (! PyArg_ParseTuple(args, "di:ldexp", &x, &exp))
 171                 return NULL;
 172         errno = 0;
 173         PyFPE_START_PROTECT("ldexp", return 0)
 174         x = ldexp(x, exp);
 175         PyFPE_END_PROTECT(x)
 176         Py_SET_ERRNO_ON_MATH_ERROR(x);
 177         if (errno && is_error(x))
 178                 return NULL;
 179         else
 180                 return PyFloat_FromDouble(x);
 181 }
 182
 183 PyDoc_STRVAR(math_ldexp_doc,
 184 "ldexp(x, i) -> x * (2**i)");
 185
 186 static PyObject *
 187 math_modf(PyObject *self, PyObject *arg)
 188 {
 189         double y, x = PyFloat_AsDouble(arg);
 190         if (x == -1.0 && PyErr_Occurred())
 191                 return NULL;
 192         errno = 0;
 193         x = modf(x, &y);
 194         Py_SET_ERRNO_ON_MATH_ERROR(x);
 195         if (errno && is_error(x))
 196                 return NULL;
 197         else
 198                 return Py_BuildValue("(dd)", x, y);
 199 }
 200
 201 PyDoc_STRVAR(math_modf_doc,
 202 "modf(x)\n"
 203 "\n"
 204 "Return the fractional and integer parts of x.  Both results carry the sign\n"
 205 "of x.  The integer part is returned as a real.");
 206
 207 /* A decent logarithm is easy to compute even for huge longs, but libm can't
 208    do that by itself -- loghelper can.  func is log or log10, and name is
 209    "log" or "log10".  Note that overflow isn't possible:  a long can contain
 210    no more than INT_MAX * SHIFT bits, so has value certainly less than
 211    2**(2**64 * 2**16) == 2**2**80, and log2 of that is 2**80, which is
 212    small enough to fit in an IEEE single.  log and log10 are even smaller.
 213 */
 214
 215 static PyObject*
 216 loghelper(PyObject* arg, double (*func)(double), char *funcname)
 217 {
 218         /* If it is long, do it ourselves. */
 219         if (PyLong_Check(arg)) {
 220                 double x;
 221                 int e;
 222                 x = _PyLong_AsScaledDouble(arg, &e);
 223                 if (x <= 0.0) {
 224                         PyErr_SetString(PyExc_ValueError,
 225                                         "math domain error");
 226                         return NULL;
 227                 }
 228                 /* Value is ~= x * 2**(e*SHIFT), so the log ~=
 229                    log(x) + log(2) * e * SHIFT.
 230                    CAUTION:  e*SHIFT may overflow using int arithmetic,
 231                    so force use of double. */
 232                 x = func(x) + (e * (double)SHIFT) * func(2.0);
 233                 return PyFloat_FromDouble(x);
 234         }
 235
 236         /* Else let libm handle it by itself. */
 237         return math_1(arg, func);
 238 }
 239
 240 static PyObject *
 241 math_log(PyObject *self, PyObject *args)
 242 {
 243         PyObject *arg;
 244         PyObject *base = NULL;
 245         PyObject *num, *den;
 246         PyObject *ans;
 247
 248         if (!PyArg_UnpackTuple(args, "log", 1, 2, &arg, &base))
 249                 return NULL;
 250
 251         num = loghelper(arg, log, "log");
 252         if (num == NULL || base == NULL)
 253                 return num;
 254
 255         den = loghelper(base, log, "log");
 256         if (den == NULL) {
 257                 Py_DECREF(num);
 258                 return NULL;
 259         }
 260
 261         ans = PyNumber_Divide(num, den);
 262         Py_DECREF(num);
 263         Py_DECREF(den);
 264         return ans;
 265 }
 266
 267 PyDoc_STRVAR(math_log_doc,
 268 "log(x[, base]) -> the logarithm of x to the given base.\n\
 269 If the base not specified, returns the natural logarithm (base e) of x.");
 270
 271 static PyObject *
 272 math_log10(PyObject *self, PyObject *arg)
 273 {
 274         return loghelper(arg, log10, "log10");
 275 }
 276
 277 PyDoc_STRVAR(math_log10_doc,
 278 "log10(x) -> the base 10 logarithm of x.");
 279
 280 /* XXX(nnorwitz): Should we use the platform M_PI or something more accurate
 281    like: 3.14159265358979323846264338327950288 */
 282 static const double degToRad = 3.141592653589793238462643383 / 180.0;
 283
 284 static PyObject *
 285 math_degrees(PyObject *self, PyObject *arg)
 286 {
 287         double x = PyFloat_AsDouble(arg);
 288         if (x == -1.0 && PyErr_Occurred())
 289                 return NULL;
 290         return PyFloat_FromDouble(x / degToRad);
 291 }
 292
 293 PyDoc_STRVAR(math_degrees_doc,
 294 "degrees(x) -> converts angle x from radians to degrees");
 295
 296 static PyObject *
 297 math_radians(PyObject *self, PyObject *arg)
 298 {
 299         double x = PyFloat_AsDouble(arg);
 300         if (x == -1.0 && PyErr_Occurred())
 301                 return NULL;
 302         return PyFloat_FromDouble(x * degToRad);
 303 }
 304
 305 PyDoc_STRVAR(math_radians_doc,
 306 "radians(x) -> converts angle x from degrees to radians");
 307
 308 static PyMethodDef math_methods[] = {
 309         {"acos",        math_acos,      METH_O,         math_acos_doc},
 310         {"asin",        math_asin,      METH_O,         math_asin_doc},
 311         {"atan",        math_atan,      METH_O,         math_atan_doc},
 312         {"atan2",       math_atan2,     METH_VARARGS,   math_atan2_doc},
 313         {"ceil",        math_ceil,      METH_O,         math_ceil_doc},
 314         {"cos",         math_cos,       METH_O,         math_cos_doc},
 315         {"cosh",        math_cosh,      METH_O,         math_cosh_doc},
 316         {"degrees",     math_degrees,   METH_O,         math_degrees_doc},
 317         {"exp",         math_exp,       METH_O,         math_exp_doc},
 318         {"fabs",        math_fabs,      METH_O,         math_fabs_doc},
 319         {"floor",       math_floor,     METH_O,         math_floor_doc},
 320         {"fmod",        math_fmod,      METH_VARARGS,   math_fmod_doc},
 321         {"frexp",       math_frexp,     METH_O,         math_frexp_doc},
 322         {"hypot",       math_hypot,     METH_VARARGS,   math_hypot_doc},
 323         {"ldexp",       math_ldexp,     METH_VARARGS,   math_ldexp_doc},
 324         {"log",         math_log,       METH_VARARGS,   math_log_doc},
 325         {"log10",       math_log10,     METH_O,         math_log10_doc},
 326         {"modf",        math_modf,      METH_O,         math_modf_doc},
 327         {"pow",         math_pow,       METH_VARARGS,   math_pow_doc},
 328         {"radians",     math_radians,   METH_O,         math_radians_doc},
 329         {"sin",         math_sin,       METH_O,         math_sin_doc},
 330         {"sinh",        math_sinh,      METH_O,         math_sinh_doc},
 331         {"sqrt",        math_sqrt,      METH_O,         math_sqrt_doc},
 332         {"tan",         math_tan,       METH_O,         math_tan_doc},
 333         {"tanh",        math_tanh,      METH_O,         math_tanh_doc},
 334         {NULL,          NULL}           /* sentinel */
 335 };
 336
 337
 338 PyDoc_STRVAR(module_doc,
 339 "This module is always available.  It provides access to the\n"
 340 "mathematical functions defined by the C standard.");
 341
 342 PyMODINIT_FUNC
 343 initmath(void)
 344 {
 345         PyObject *m, *d, *v;
 346
 347         m = Py_InitModule3("math", math_methods, module_doc);
 348         if (m == NULL)
 349                 goto finally;
 350         d = PyModule_GetDict(m);
 351         if (d == NULL)
 352                 goto finally;
 353
 354         if (!(v = PyFloat_FromDouble(atan(1.0) * 4.0)))
 355                 goto finally;
 356         if (PyDict_SetItemString(d, "pi", v) < 0)
 357                 goto finally;
 358         Py_DECREF(v);
 359
 360         if (!(v = PyFloat_FromDouble(exp(1.0))))
 361                 goto finally;
 362         if (PyDict_SetItemString(d, "e", v) < 0)
 363                 goto finally;
 364         Py_DECREF(v);
 365
 366   finally:
 367         return;
 368 }