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