1 /* Math module -- standard C math library functions, pi and e */
3 /* Here are some comments from Tim Peters, extracted from the
4 discussion attached to http://bugs.python.org/issue1640. They
5 describe the general aims of the math module with respect to
6 special values, IEEE-754 floating-point exceptions, and Python
9 These are the "spirit of 754" rules:
11 1. If the mathematical result is a real number, but of magnitude too
12 large to approximate by a machine float, overflow is signaled and the
13 result is an infinity (with the appropriate sign).
15 2. If the mathematical result is a real number, but of magnitude too
16 small to approximate by a machine float, underflow is signaled and the
17 result is a zero (with the appropriate sign).
19 3. At a singularity (a value x such that the limit of f(y) as y
20 approaches x exists and is an infinity), "divide by zero" is signaled
21 and the result is an infinity (with the appropriate sign). This is
22 complicated a little by that the left-side and right-side limits may
23 not be the same; e.g., 1/x approaches +inf or -inf as x approaches 0
24 from the positive or negative directions. In that specific case, the
25 sign of the zero determines the result of 1/0.
27 4. At a point where a function has no defined result in the extended
28 reals (i.e., the reals plus an infinity or two), invalid operation is
29 signaled and a NaN is returned.
31 And these are what Python has historically /tried/ to do (but not
32 always successfully, as platform libm behavior varies a lot):
34 For #1, raise OverflowError.
36 For #2, return a zero (with the appropriate sign if that happens by
39 For #3 and #4, raise ValueError. It may have made sense to raise
40 Python's ZeroDivisionError in #3, but historically that's only been
41 raised for division by zero and mod by zero.
46 In general, on an IEEE-754 platform the aim is to follow the C99
47 standard, including Annex 'F', whenever possible. Where the
48 standard recommends raising the 'divide-by-zero' or 'invalid'
49 floating-point exceptions, Python should raise a ValueError. Where
50 the standard recommends raising 'overflow', Python should raise an
51 OverflowError. In all other circumstances a value should be
56 #include "longintrepr.h" /* just for SHIFT */
59 /* OSF1 5.1 doesn't make this available with XOPEN_SOURCE_EXTENDED defined */
60 extern double copysign(double, double);
63 /* Call is_error when errno != 0, and where x is the result libm
64 * returned. is_error will usually set up an exception and return
65 * true (1), but may return false (0) without setting up an exception.
70 int result
= 1; /* presumption of guilt */
71 assert(errno
); /* non-zero errno is a precondition for calling */
73 PyErr_SetString(PyExc_ValueError
, "math domain error");
75 else if (errno
== ERANGE
) {
76 /* ANSI C generally requires libm functions to set ERANGE
77 * on overflow, but also generally *allows* them to set
78 * ERANGE on underflow too. There's no consistency about
79 * the latter across platforms.
80 * Alas, C99 never requires that errno be set.
81 * Here we suppress the underflow errors (libm functions
82 * should return a zero on underflow, and +- HUGE_VAL on
83 * overflow, so testing the result for zero suffices to
84 * distinguish the cases).
87 PyErr_SetString(PyExc_OverflowError
,
93 /* Unexpected math error */
94 PyErr_SetFromErrno(PyExc_ValueError
);
99 wrapper for atan2 that deals directly with special cases before
100 delegating to the platform libm for the remaining cases. This
101 is necessary to get consistent behaviour across platforms.
102 Windows, FreeBSD and alpha Tru64 are amongst platforms that don't
107 m_atan2(double y
, double x
)
109 if (Py_IS_NAN(x
) || Py_IS_NAN(y
))
111 if (Py_IS_INFINITY(y
)) {
112 if (Py_IS_INFINITY(x
)) {
113 if (copysign(1., x
) == 1.)
114 /* atan2(+-inf, +inf) == +-pi/4 */
115 return copysign(0.25*Py_MATH_PI
, y
);
117 /* atan2(+-inf, -inf) == +-pi*3/4 */
118 return copysign(0.75*Py_MATH_PI
, y
);
120 /* atan2(+-inf, x) == +-pi/2 for finite x */
121 return copysign(0.5*Py_MATH_PI
, y
);
123 if (Py_IS_INFINITY(x
) || y
== 0.) {
124 if (copysign(1., x
) == 1.)
125 /* atan2(+-y, +inf) = atan2(+-0, +x) = +-0. */
126 return copysign(0., y
);
128 /* atan2(+-y, -inf) = atan2(+-0., -x) = +-pi. */
129 return copysign(Py_MATH_PI
, y
);
135 math_1 is used to wrap a libm function f that takes a double
136 arguments and returns a double.
138 The error reporting follows these rules, which are designed to do
139 the right thing on C89/C99 platforms and IEEE 754/non IEEE 754
142 - a NaN result from non-NaN inputs causes ValueError to be raised
143 - an infinite result from finite inputs causes OverflowError to be
144 raised if can_overflow is 1, or raises ValueError if can_overflow
146 - if the result is finite and errno == EDOM then ValueError is
148 - if the result is finite and nonzero and errno == ERANGE then
149 OverflowError is raised
151 The last rule is used to catch overflow on platforms which follow
152 C89 but for which HUGE_VAL is not an infinity.
154 For the majority of one-argument functions these rules are enough
155 to ensure that Python's functions behave as specified in 'Annex F'
156 of the C99 standard, with the 'invalid' and 'divide-by-zero'
157 floating-point exceptions mapping to Python's ValueError and the
158 'overflow' floating-point exception mapping to OverflowError.
159 math_1 only works for functions that don't have singularities *and*
160 the possibility of overflow; fortunately, that covers everything we
161 care about right now.
165 math_1(PyObject
*arg
, double (*func
) (double), int can_overflow
)
168 x
= PyFloat_AsDouble(arg
);
169 if (x
== -1.0 && PyErr_Occurred())
172 PyFPE_START_PROTECT("in math_1", return 0);
174 PyFPE_END_PROTECT(r
);
181 else if (Py_IS_INFINITY(r
)) {
183 errno
= can_overflow
? ERANGE
: EDOM
;
187 if (errno
&& is_error(r
))
190 return PyFloat_FromDouble(r
);
194 math_2 is used to wrap a libm function f that takes two double
195 arguments and returns a double.
197 The error reporting follows these rules, which are designed to do
198 the right thing on C89/C99 platforms and IEEE 754/non IEEE 754
201 - a NaN result from non-NaN inputs causes ValueError to be raised
202 - an infinite result from finite inputs causes OverflowError to be
204 - if the result is finite and errno == EDOM then ValueError is
206 - if the result is finite and nonzero and errno == ERANGE then
207 OverflowError is raised
209 The last rule is used to catch overflow on platforms which follow
210 C89 but for which HUGE_VAL is not an infinity.
212 For most two-argument functions (copysign, fmod, hypot, atan2)
213 these rules are enough to ensure that Python's functions behave as
214 specified in 'Annex F' of the C99 standard, with the 'invalid' and
215 'divide-by-zero' floating-point exceptions mapping to Python's
216 ValueError and the 'overflow' floating-point exception mapping to
221 math_2(PyObject
*args
, double (*func
) (double, double), char *funcname
)
225 if (! PyArg_UnpackTuple(args
, funcname
, 2, 2, &ox
, &oy
))
227 x
= PyFloat_AsDouble(ox
);
228 y
= PyFloat_AsDouble(oy
);
229 if ((x
== -1.0 || y
== -1.0) && PyErr_Occurred())
232 PyFPE_START_PROTECT("in math_2", return 0);
234 PyFPE_END_PROTECT(r
);
236 if (!Py_IS_NAN(x
) && !Py_IS_NAN(y
))
241 else if (Py_IS_INFINITY(r
)) {
242 if (Py_IS_FINITE(x
) && Py_IS_FINITE(y
))
247 if (errno
&& is_error(r
))
250 return PyFloat_FromDouble(r
);
253 #define FUNC1(funcname, func, can_overflow, docstring) \
254 static PyObject * math_##funcname(PyObject *self, PyObject *args) { \
255 return math_1(args, func, can_overflow); \
257 PyDoc_STRVAR(math_##funcname##_doc, docstring);
259 #define FUNC2(funcname, func, docstring) \
260 static PyObject * math_##funcname(PyObject *self, PyObject *args) { \
261 return math_2(args, func, #funcname); \
263 PyDoc_STRVAR(math_##funcname##_doc, docstring);
266 "acos(x)\n\nReturn the arc cosine (measured in radians) of x.")
267 FUNC1(acosh
, acosh
, 0,
268 "acosh(x)\n\nReturn the hyperbolic arc cosine (measured in radians) of x.")
270 "asin(x)\n\nReturn the arc sine (measured in radians) of x.")
271 FUNC1(asinh
, asinh
, 0,
272 "asinh(x)\n\nReturn the hyperbolic arc sine (measured in radians) of x.")
274 "atan(x)\n\nReturn the arc tangent (measured in radians) of x.")
275 FUNC2(atan2
, m_atan2
,
276 "atan2(y, x)\n\nReturn the arc tangent (measured in radians) of y/x.\n"
277 "Unlike atan(y/x), the signs of both x and y are considered.")
278 FUNC1(atanh
, atanh
, 0,
279 "atanh(x)\n\nReturn the hyperbolic arc tangent (measured in radians) of x.")
281 "ceil(x)\n\nReturn the ceiling of x as a float.\n"
282 "This is the smallest integral value >= x.")
283 FUNC2(copysign
, copysign
,
284 "copysign(x,y)\n\nReturn x with the sign of y.")
286 "cos(x)\n\nReturn the cosine of x (measured in radians).")
288 "cosh(x)\n\nReturn the hyperbolic cosine of x.")
290 "exp(x)\n\nReturn e raised to the power of x.")
292 "fabs(x)\n\nReturn the absolute value of the float x.")
293 FUNC1(floor
, floor
, 0,
294 "floor(x)\n\nReturn the floor of x as a float.\n"
295 "This is the largest integral value <= x.")
296 FUNC1(log1p
, log1p
, 1,
297 "log1p(x)\n\nReturn the natural logarithm of 1+x (base e).\n\
298 The result is computed in a way which is accurate for x near zero.")
300 "sin(x)\n\nReturn the sine of x (measured in radians).")
302 "sinh(x)\n\nReturn the hyperbolic sine of x.")
304 "sqrt(x)\n\nReturn the square root of x.")
306 "tan(x)\n\nReturn the tangent of x (measured in radians).")
308 "tanh(x)\n\nReturn the hyperbolic tangent of x.")
310 /* Precision summation function as msum() by Raymond Hettinger in
311 <http://aspn.activestate.com/ASPN/Cookbook/Python/Recipe/393090>,
312 enhanced with the exact partials sum and roundoff from Mark
313 Dickinson's post at <http://bugs.python.org/file10357/msum4.py>.
314 See those links for more details, proofs and other references.
316 Note 1: IEEE 754R floating point semantics are assumed,
317 but the current implementation does not re-establish special
318 value semantics across iterations (i.e. handling -Inf + Inf).
320 Note 2: No provision is made for intermediate overflow handling;
321 therefore, sum([1e+308, 1e-308, 1e+308]) returns 1e+308 while
322 sum([1e+308, 1e+308, 1e-308]) raises an OverflowError due to the
323 overflow of the first partial sum.
325 Note 3: The itermediate values lo, yr, and hi are declared volatile so
326 aggressive compilers won't algebraicly reduce lo to always be exactly 0.0.
327 Also, the volatile declaration forces the values to be stored in memory as
328 regular doubles instead of extended long precision (80-bit) values. This
329 prevents double rounding because any addition or substraction of two doubles
330 can be resolved exactly into double-sized hi and lo values. As long as the
331 hi value gets forced into a double before yr and lo are computed, the extra
332 bits in downstream extended precision operations (x87 for example) will be
333 exactly zero and therefore can be losslessly stored back into a double,
334 thereby preventing double rounding.
336 Note 4: A similar implementation is in Modules/cmathmodule.c.
337 Be sure to update both when making changes.
339 Note 5: The signature of math.sum() differs from __builtin__.sum()
340 because the start argument doesn't make sense in the context of
341 accurate summation. Since the partials table is collapsed before
342 returning a result, sum(seq2, start=sum(seq1)) may not equal the
343 accurate result returned by sum(itertools.chain(seq1, seq2)).
346 #define NUM_PARTIALS 32 /* initial partials array size, on stack */
348 /* Extend the partials array p[] by doubling its size. */
349 static int /* non-zero on error */
350 _sum_realloc(double **p_ptr
, Py_ssize_t n
,
351 double *ps
, Py_ssize_t
*m_ptr
)
354 Py_ssize_t m
= *m_ptr
;
357 if (n
< m
&& m
< (PY_SSIZE_T_MAX
/ sizeof(double))) {
360 v
= PyMem_Malloc(sizeof(double) * m
);
362 memcpy(v
, ps
, sizeof(double) * n
);
365 v
= PyMem_Realloc(p
, sizeof(double) * m
);
367 if (v
== NULL
) { /* size overflow or no memory */
368 PyErr_SetString(PyExc_MemoryError
, "math sum partials");
371 *p_ptr
= (double*) v
;
376 /* Full precision summation of a sequence of floats.
379 partials = [] # sorted, non-overlapping partial sums
392 return sum_exact(partials)
394 Rounded x+y stored in hi with the roundoff stored in lo. Together hi+lo
395 are exactly equal to x+y. The inner loop applies hi/lo summation to each
396 partial so that the list of partial sums remains exact.
398 Sum_exact() adds the partial sums exactly and correctly rounds the final
399 result (using the round-half-to-even rule). The items in partials remain
400 non-zero, non-special, non-overlapping and strictly increasing in
401 magnitude, but possibly not all having the same sign.
403 Depends on IEEE 754 arithmetic guarantees and half-even rounding.
407 math_sum(PyObject
*self
, PyObject
*seq
)
409 PyObject
*item
, *iter
, *sum
= NULL
;
410 Py_ssize_t i
, j
, n
= 0, m
= NUM_PARTIALS
;
411 double x
, y
, t
, ps
[NUM_PARTIALS
], *p
= ps
;
412 volatile double hi
, yr
, lo
;
414 iter
= PyObject_GetIter(seq
);
418 PyFPE_START_PROTECT("sum", Py_DECREF(iter
); return NULL
)
420 for(;;) { /* for x in iterable */
421 assert(0 <= n
&& n
<= m
);
422 assert((m
== NUM_PARTIALS
&& p
== ps
) ||
423 (m
> NUM_PARTIALS
&& p
!= NULL
));
425 item
= PyIter_Next(iter
);
427 if (PyErr_Occurred())
431 x
= PyFloat_AsDouble(item
);
433 if (PyErr_Occurred())
436 for (i
= j
= 0; j
< n
; j
++) { /* for y in partials */
438 if (fabs(x
) < fabs(y
)) {
449 n
= i
; /* ps[i:] = [x] */
451 /* If non-finite, reset partials, effectively
452 adding subsequent items without roundoff
453 and yielding correct non-finite results,
454 provided IEEE 754 rules are observed */
455 if (! Py_IS_FINITE(x
))
457 else if (n
>= m
&& _sum_realloc(&p
, n
, ps
, &m
))
466 if (Py_IS_FINITE(hi
)) {
467 /* sum_exact(ps, hi) from the top, stop when the sum becomes inexact. */
471 assert(fabs(y
) < fabs(x
));
478 /* Make half-even rounding work across multiple partials. Needed
479 so that sum([1e-16, 1, 1e16]) will round-up the last digit to
480 two instead of down to zero (the 1e-16 makes the 1 slightly
481 closer to two). With a potential 1 ULP rounding error fixed-up,
482 math.sum() can guarantee commutativity. */
483 if (n
> 0 && ((lo
< 0.0 && p
[n
-1] < 0.0) ||
484 (lo
> 0.0 && p
[n
-1] > 0.0))) {
492 else { /* raise exception corresponding to a special value */
493 errno
= Py_IS_NAN(hi
) ? EDOM
: ERANGE
;
498 sum
= PyFloat_FromDouble(hi
);
501 PyFPE_END_PROTECT(hi
)
510 PyDoc_STRVAR(math_sum_doc
,
512 Return an accurate floating point sum of values in the iterable.\n\
513 Assumes IEEE-754 floating point arithmetic.");
516 math_trunc(PyObject
*self
, PyObject
*number
)
518 return PyObject_CallMethod(number
, "__trunc__", NULL
);
521 PyDoc_STRVAR(math_trunc_doc
,
522 "trunc(x:Real) -> Integral\n"
524 "Truncates x to the nearest Integral toward 0. Uses the __trunc__ magic method.");
527 math_frexp(PyObject
*self
, PyObject
*arg
)
530 double x
= PyFloat_AsDouble(arg
);
531 if (x
== -1.0 && PyErr_Occurred())
533 /* deal with special cases directly, to sidestep platform
535 if (Py_IS_NAN(x
) || Py_IS_INFINITY(x
) || !x
) {
539 PyFPE_START_PROTECT("in math_frexp", return 0);
541 PyFPE_END_PROTECT(x
);
543 return Py_BuildValue("(di)", x
, i
);
546 PyDoc_STRVAR(math_frexp_doc
,
549 "Return the mantissa and exponent of x, as pair (m, e).\n"
550 "m is a float and e is an int, such that x = m * 2.**e.\n"
551 "If x is 0, m and e are both 0. Else 0.5 <= abs(m) < 1.0.");
554 math_ldexp(PyObject
*self
, PyObject
*args
)
559 if (! PyArg_ParseTuple(args
, "dO:ldexp", &x
, &oexp
))
562 if (PyLong_Check(oexp
)) {
563 /* on overflow, replace exponent with either LONG_MAX
564 or LONG_MIN, depending on the sign. */
565 exp
= PyLong_AsLong(oexp
);
566 if (exp
== -1 && PyErr_Occurred()) {
567 if (PyErr_ExceptionMatches(PyExc_OverflowError
)) {
568 if (Py_SIZE(oexp
) < 0) {
577 /* propagate any unexpected exception */
582 else if (PyInt_Check(oexp
)) {
583 exp
= PyInt_AS_LONG(oexp
);
586 PyErr_SetString(PyExc_TypeError
,
587 "Expected an int or long as second argument "
592 if (x
== 0. || !Py_IS_FINITE(x
)) {
593 /* NaNs, zeros and infinities are returned unchanged */
596 } else if (exp
> INT_MAX
) {
598 r
= copysign(Py_HUGE_VAL
, x
);
600 } else if (exp
< INT_MIN
) {
601 /* underflow to +-0 */
606 PyFPE_START_PROTECT("in math_ldexp", return 0);
607 r
= ldexp(x
, (int)exp
);
608 PyFPE_END_PROTECT(r
);
609 if (Py_IS_INFINITY(r
))
613 if (errno
&& is_error(r
))
615 return PyFloat_FromDouble(r
);
618 PyDoc_STRVAR(math_ldexp_doc
,
619 "ldexp(x, i) -> x * (2**i)");
622 math_modf(PyObject
*self
, PyObject
*arg
)
624 double y
, x
= PyFloat_AsDouble(arg
);
625 if (x
== -1.0 && PyErr_Occurred())
627 /* some platforms don't do the right thing for NaNs and
628 infinities, so we take care of special cases directly. */
629 if (!Py_IS_FINITE(x
)) {
630 if (Py_IS_INFINITY(x
))
631 return Py_BuildValue("(dd)", copysign(0., x
), x
);
632 else if (Py_IS_NAN(x
))
633 return Py_BuildValue("(dd)", x
, x
);
637 PyFPE_START_PROTECT("in math_modf", return 0);
639 PyFPE_END_PROTECT(x
);
640 return Py_BuildValue("(dd)", x
, y
);
643 PyDoc_STRVAR(math_modf_doc
,
646 "Return the fractional and integer parts of x. Both results carry the sign\n"
647 "of x. The integer part is returned as a real.");
649 /* A decent logarithm is easy to compute even for huge longs, but libm can't
650 do that by itself -- loghelper can. func is log or log10, and name is
651 "log" or "log10". Note that overflow isn't possible: a long can contain
652 no more than INT_MAX * SHIFT bits, so has value certainly less than
653 2**(2**64 * 2**16) == 2**2**80, and log2 of that is 2**80, which is
654 small enough to fit in an IEEE single. log and log10 are even smaller.
658 loghelper(PyObject
* arg
, double (*func
)(double), char *funcname
)
660 /* If it is long, do it ourselves. */
661 if (PyLong_Check(arg
)) {
664 x
= _PyLong_AsScaledDouble(arg
, &e
);
666 PyErr_SetString(PyExc_ValueError
,
667 "math domain error");
670 /* Value is ~= x * 2**(e*PyLong_SHIFT), so the log ~=
671 log(x) + log(2) * e * PyLong_SHIFT.
672 CAUTION: e*PyLong_SHIFT may overflow using int arithmetic,
673 so force use of double. */
674 x
= func(x
) + (e
* (double)PyLong_SHIFT
) * func(2.0);
675 return PyFloat_FromDouble(x
);
678 /* Else let libm handle it by itself. */
679 return math_1(arg
, func
, 0);
683 math_log(PyObject
*self
, PyObject
*args
)
686 PyObject
*base
= NULL
;
690 if (!PyArg_UnpackTuple(args
, "log", 1, 2, &arg
, &base
))
693 num
= loghelper(arg
, log
, "log");
694 if (num
== NULL
|| base
== NULL
)
697 den
= loghelper(base
, log
, "log");
703 ans
= PyNumber_Divide(num
, den
);
709 PyDoc_STRVAR(math_log_doc
,
710 "log(x[, base]) -> the logarithm of x to the given base.\n\
711 If the base not specified, returns the natural logarithm (base e) of x.");
714 math_log10(PyObject
*self
, PyObject
*arg
)
716 return loghelper(arg
, log10
, "log10");
719 PyDoc_STRVAR(math_log10_doc
,
720 "log10(x) -> the base 10 logarithm of x.");
723 math_fmod(PyObject
*self
, PyObject
*args
)
727 if (! PyArg_UnpackTuple(args
, "fmod", 2, 2, &ox
, &oy
))
729 x
= PyFloat_AsDouble(ox
);
730 y
= PyFloat_AsDouble(oy
);
731 if ((x
== -1.0 || y
== -1.0) && PyErr_Occurred())
733 /* fmod(x, +/-Inf) returns x for finite x. */
734 if (Py_IS_INFINITY(y
) && Py_IS_FINITE(x
))
735 return PyFloat_FromDouble(x
);
737 PyFPE_START_PROTECT("in math_fmod", return 0);
739 PyFPE_END_PROTECT(r
);
741 if (!Py_IS_NAN(x
) && !Py_IS_NAN(y
))
746 if (errno
&& is_error(r
))
749 return PyFloat_FromDouble(r
);
752 PyDoc_STRVAR(math_fmod_doc
,
753 "fmod(x,y)\n\nReturn fmod(x, y), according to platform C."
754 " x % y may differ.");
757 math_hypot(PyObject
*self
, PyObject
*args
)
761 if (! PyArg_UnpackTuple(args
, "hypot", 2, 2, &ox
, &oy
))
763 x
= PyFloat_AsDouble(ox
);
764 y
= PyFloat_AsDouble(oy
);
765 if ((x
== -1.0 || y
== -1.0) && PyErr_Occurred())
767 /* hypot(x, +/-Inf) returns Inf, even if x is a NaN. */
768 if (Py_IS_INFINITY(x
))
769 return PyFloat_FromDouble(fabs(x
));
770 if (Py_IS_INFINITY(y
))
771 return PyFloat_FromDouble(fabs(y
));
773 PyFPE_START_PROTECT("in math_hypot", return 0);
775 PyFPE_END_PROTECT(r
);
777 if (!Py_IS_NAN(x
) && !Py_IS_NAN(y
))
782 else if (Py_IS_INFINITY(r
)) {
783 if (Py_IS_FINITE(x
) && Py_IS_FINITE(y
))
788 if (errno
&& is_error(r
))
791 return PyFloat_FromDouble(r
);
794 PyDoc_STRVAR(math_hypot_doc
,
795 "hypot(x,y)\n\nReturn the Euclidean distance, sqrt(x*x + y*y).");
797 /* pow can't use math_2, but needs its own wrapper: the problem is
798 that an infinite result can arise either as a result of overflow
799 (in which case OverflowError should be raised) or as a result of
800 e.g. 0.**-5. (for which ValueError needs to be raised.)
804 math_pow(PyObject
*self
, PyObject
*args
)
810 if (! PyArg_UnpackTuple(args
, "pow", 2, 2, &ox
, &oy
))
812 x
= PyFloat_AsDouble(ox
);
813 y
= PyFloat_AsDouble(oy
);
814 if ((x
== -1.0 || y
== -1.0) && PyErr_Occurred())
817 /* deal directly with IEEE specials, to cope with problems on various
818 platforms whose semantics don't exactly match C99 */
819 r
= 0.; /* silence compiler warning */
820 if (!Py_IS_FINITE(x
) || !Py_IS_FINITE(y
)) {
823 r
= y
== 0. ? 1. : x
; /* NaN**0 = 1 */
824 else if (Py_IS_NAN(y
))
825 r
= x
== 1. ? 1. : y
; /* 1**NaN = 1 */
826 else if (Py_IS_INFINITY(x
)) {
827 odd_y
= Py_IS_FINITE(y
) && fmod(fabs(y
), 2.0) == 1.0;
829 r
= odd_y
? x
: fabs(x
);
833 r
= odd_y
? copysign(0., x
) : 0.;
835 else if (Py_IS_INFINITY(y
)) {
838 else if (y
> 0. && fabs(x
) > 1.0)
840 else if (y
< 0. && fabs(x
) < 1.0) {
841 r
= -y
; /* result is +inf */
842 if (x
== 0.) /* 0**-inf: divide-by-zero */
850 /* let libm handle finite**finite */
852 PyFPE_START_PROTECT("in math_pow", return 0);
854 PyFPE_END_PROTECT(r
);
855 /* a NaN result should arise only from (-ve)**(finite
856 non-integer); in this case we want to raise ValueError. */
857 if (!Py_IS_FINITE(r
)) {
862 an infinite result here arises either from:
863 (A) (+/-0.)**negative (-> divide-by-zero)
864 (B) overflow of x**y with x and y finite
866 else if (Py_IS_INFINITY(r
)) {
875 if (errno
&& is_error(r
))
878 return PyFloat_FromDouble(r
);
881 PyDoc_STRVAR(math_pow_doc
,
882 "pow(x,y)\n\nReturn x**y (x to the power of y).");
884 static const double degToRad
= Py_MATH_PI
/ 180.0;
885 static const double radToDeg
= 180.0 / Py_MATH_PI
;
888 math_degrees(PyObject
*self
, PyObject
*arg
)
890 double x
= PyFloat_AsDouble(arg
);
891 if (x
== -1.0 && PyErr_Occurred())
893 return PyFloat_FromDouble(x
* radToDeg
);
896 PyDoc_STRVAR(math_degrees_doc
,
897 "degrees(x) -> converts angle x from radians to degrees");
900 math_radians(PyObject
*self
, PyObject
*arg
)
902 double x
= PyFloat_AsDouble(arg
);
903 if (x
== -1.0 && PyErr_Occurred())
905 return PyFloat_FromDouble(x
* degToRad
);
908 PyDoc_STRVAR(math_radians_doc
,
909 "radians(x) -> converts angle x from degrees to radians");
912 math_isnan(PyObject
*self
, PyObject
*arg
)
914 double x
= PyFloat_AsDouble(arg
);
915 if (x
== -1.0 && PyErr_Occurred())
917 return PyBool_FromLong((long)Py_IS_NAN(x
));
920 PyDoc_STRVAR(math_isnan_doc
,
922 Checks if float x is not a number (NaN)");
925 math_isinf(PyObject
*self
, PyObject
*arg
)
927 double x
= PyFloat_AsDouble(arg
);
928 if (x
== -1.0 && PyErr_Occurred())
930 return PyBool_FromLong((long)Py_IS_INFINITY(x
));
933 PyDoc_STRVAR(math_isinf_doc
,
935 Checks if float x is infinite (positive or negative)");
937 static PyMethodDef math_methods
[] = {
938 {"acos", math_acos
, METH_O
, math_acos_doc
},
939 {"acosh", math_acosh
, METH_O
, math_acosh_doc
},
940 {"asin", math_asin
, METH_O
, math_asin_doc
},
941 {"asinh", math_asinh
, METH_O
, math_asinh_doc
},
942 {"atan", math_atan
, METH_O
, math_atan_doc
},
943 {"atan2", math_atan2
, METH_VARARGS
, math_atan2_doc
},
944 {"atanh", math_atanh
, METH_O
, math_atanh_doc
},
945 {"ceil", math_ceil
, METH_O
, math_ceil_doc
},
946 {"copysign", math_copysign
, METH_VARARGS
, math_copysign_doc
},
947 {"cos", math_cos
, METH_O
, math_cos_doc
},
948 {"cosh", math_cosh
, METH_O
, math_cosh_doc
},
949 {"degrees", math_degrees
, METH_O
, math_degrees_doc
},
950 {"exp", math_exp
, METH_O
, math_exp_doc
},
951 {"fabs", math_fabs
, METH_O
, math_fabs_doc
},
952 {"floor", math_floor
, METH_O
, math_floor_doc
},
953 {"fmod", math_fmod
, METH_VARARGS
, math_fmod_doc
},
954 {"frexp", math_frexp
, METH_O
, math_frexp_doc
},
955 {"hypot", math_hypot
, METH_VARARGS
, math_hypot_doc
},
956 {"isinf", math_isinf
, METH_O
, math_isinf_doc
},
957 {"isnan", math_isnan
, METH_O
, math_isnan_doc
},
958 {"ldexp", math_ldexp
, METH_VARARGS
, math_ldexp_doc
},
959 {"log", math_log
, METH_VARARGS
, math_log_doc
},
960 {"log1p", math_log1p
, METH_O
, math_log1p_doc
},
961 {"log10", math_log10
, METH_O
, math_log10_doc
},
962 {"modf", math_modf
, METH_O
, math_modf_doc
},
963 {"pow", math_pow
, METH_VARARGS
, math_pow_doc
},
964 {"radians", math_radians
, METH_O
, math_radians_doc
},
965 {"sin", math_sin
, METH_O
, math_sin_doc
},
966 {"sinh", math_sinh
, METH_O
, math_sinh_doc
},
967 {"sqrt", math_sqrt
, METH_O
, math_sqrt_doc
},
968 {"sum", math_sum
, METH_O
, math_sum_doc
},
969 {"tan", math_tan
, METH_O
, math_tan_doc
},
970 {"tanh", math_tanh
, METH_O
, math_tanh_doc
},
971 {"trunc", math_trunc
, METH_O
, math_trunc_doc
},
972 {NULL
, NULL
} /* sentinel */
976 PyDoc_STRVAR(module_doc
,
977 "This module is always available. It provides access to the\n"
978 "mathematical functions defined by the C standard.");
985 m
= Py_InitModule3("math", math_methods
, module_doc
);
989 PyModule_AddObject(m
, "pi", PyFloat_FromDouble(Py_MATH_PI
));
990 PyModule_AddObject(m
, "e", PyFloat_FromDouble(Py_MATH_E
));