blob: 0e91f8fd6a48604ca59673b3038821342a30b00c [file] [log] [blame]
Guido van Rossum85a5fbb1990-10-14 12:07:46 +00001/* Math module -- standard C math library functions, pi and e */
2
Christian Heimes6f341092008-04-18 23:13:07 +00003/* 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
7 exceptions.
8
9These are the "spirit of 754" rules:
10
111. If the mathematical result is a real number, but of magnitude too
12large to approximate by a machine float, overflow is signaled and the
13result is an infinity (with the appropriate sign).
14
152. If the mathematical result is a real number, but of magnitude too
16small to approximate by a machine float, underflow is signaled and the
17result is a zero (with the appropriate sign).
18
193. At a singularity (a value x such that the limit of f(y) as y
20approaches x exists and is an infinity), "divide by zero" is signaled
21and the result is an infinity (with the appropriate sign). This is
22complicated a little by that the left-side and right-side limits may
23not be the same; e.g., 1/x approaches +inf or -inf as x approaches 0
24from the positive or negative directions. In that specific case, the
25sign of the zero determines the result of 1/0.
26
274. At a point where a function has no defined result in the extended
28reals (i.e., the reals plus an infinity or two), invalid operation is
29signaled and a NaN is returned.
30
31And these are what Python has historically /tried/ to do (but not
32always successfully, as platform libm behavior varies a lot):
33
34For #1, raise OverflowError.
35
36For #2, return a zero (with the appropriate sign if that happens by
37accident ;-)).
38
39For #3 and #4, raise ValueError. It may have made sense to raise
40Python's ZeroDivisionError in #3, but historically that's only been
41raised for division by zero and mod by zero.
42
43*/
44
45/*
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
52 returned.
53 */
54
Barry Warsaw8b43b191996-12-09 22:32:36 +000055#include "Python.h"
Michael W. Hudson9ef852c2005-04-06 13:05:18 +000056#include "longintrepr.h" /* just for SHIFT */
Guido van Rossum85a5fbb1990-10-14 12:07:46 +000057
Neal Norwitz5f95a792008-01-25 08:04:16 +000058#ifdef _OSF_SOURCE
59/* OSF1 5.1 doesn't make this available with XOPEN_SOURCE_EXTENDED defined */
60extern double copysign(double, double);
61#endif
62
Tim Peters1d120612000-10-12 06:10:25 +000063/* 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.
66 */
67static int
68is_error(double x)
Guido van Rossum8832b621991-12-16 15:44:24 +000069{
Tim Peters1d120612000-10-12 06:10:25 +000070 int result = 1; /* presumption of guilt */
Tim Peters2bf405a2000-10-12 19:42:00 +000071 assert(errno); /* non-zero errno is a precondition for calling */
Guido van Rossum8832b621991-12-16 15:44:24 +000072 if (errno == EDOM)
Barry Warsaw8b43b191996-12-09 22:32:36 +000073 PyErr_SetString(PyExc_ValueError, "math domain error");
Tim Petersa40c7932001-09-05 22:36:56 +000074
Tim Peters1d120612000-10-12 06:10:25 +000075 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
Tim Petersa40c7932001-09-05 22:36:56 +000079 * 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).
Tim Peters1d120612000-10-12 06:10:25 +000085 */
86 if (x)
Tim Petersfe71f812001-08-07 22:10:00 +000087 PyErr_SetString(PyExc_OverflowError,
Tim Peters1d120612000-10-12 06:10:25 +000088 "math range error");
89 else
90 result = 0;
91 }
Guido van Rossum8832b621991-12-16 15:44:24 +000092 else
Barry Warsaw8b43b191996-12-09 22:32:36 +000093 /* Unexpected math error */
94 PyErr_SetFromErrno(PyExc_ValueError);
Tim Peters1d120612000-10-12 06:10:25 +000095 return result;
Guido van Rossum8832b621991-12-16 15:44:24 +000096}
97
Christian Heimes6f341092008-04-18 23:13:07 +000098/*
99 math_1 is used to wrap a libm function f that takes a double
100 arguments and returns a double.
101
102 The error reporting follows these rules, which are designed to do
103 the right thing on C89/C99 platforms and IEEE 754/non IEEE 754
104 platforms.
105
106 - a NaN result from non-NaN inputs causes ValueError to be raised
107 - an infinite result from finite inputs causes OverflowError to be
108 raised if can_overflow is 1, or raises ValueError if can_overflow
109 is 0.
110 - if the result is finite and errno == EDOM then ValueError is
111 raised
112 - if the result is finite and nonzero and errno == ERANGE then
113 OverflowError is raised
114
115 The last rule is used to catch overflow on platforms which follow
116 C89 but for which HUGE_VAL is not an infinity.
117
118 For the majority of one-argument functions these rules are enough
119 to ensure that Python's functions behave as specified in 'Annex F'
120 of the C99 standard, with the 'invalid' and 'divide-by-zero'
121 floating-point exceptions mapping to Python's ValueError and the
122 'overflow' floating-point exception mapping to OverflowError.
123 math_1 only works for functions that don't have singularities *and*
124 the possibility of overflow; fortunately, that covers everything we
125 care about right now.
126*/
127
Barry Warsaw8b43b191996-12-09 22:32:36 +0000128static PyObject *
Christian Heimes6f341092008-04-18 23:13:07 +0000129math_1(PyObject *arg, double (*func) (double), int can_overflow)
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000130{
Christian Heimes6f341092008-04-18 23:13:07 +0000131 double x, r;
132 x = PyFloat_AsDouble(arg);
Neal Norwitz45e230a2006-11-19 21:26:53 +0000133 if (x == -1.0 && PyErr_Occurred())
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000134 return NULL;
135 errno = 0;
Christian Heimes6f341092008-04-18 23:13:07 +0000136 PyFPE_START_PROTECT("in math_1", return 0);
137 r = (*func)(x);
138 PyFPE_END_PROTECT(r);
139 if (Py_IS_NAN(r)) {
140 if (!Py_IS_NAN(x))
141 errno = EDOM;
142 else
143 errno = 0;
144 }
145 else if (Py_IS_INFINITY(r)) {
146 if (Py_IS_FINITE(x))
147 errno = can_overflow ? ERANGE : EDOM;
148 else
149 errno = 0;
150 }
151 if (errno && is_error(r))
Tim Peters1d120612000-10-12 06:10:25 +0000152 return NULL;
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000153 else
Christian Heimes6f341092008-04-18 23:13:07 +0000154 return PyFloat_FromDouble(r);
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000155}
156
Christian Heimes6f341092008-04-18 23:13:07 +0000157/*
158 math_2 is used to wrap a libm function f that takes two double
159 arguments and returns a double.
160
161 The error reporting follows these rules, which are designed to do
162 the right thing on C89/C99 platforms and IEEE 754/non IEEE 754
163 platforms.
164
165 - a NaN result from non-NaN inputs causes ValueError to be raised
166 - an infinite result from finite inputs causes OverflowError to be
167 raised.
168 - if the result is finite and errno == EDOM then ValueError is
169 raised
170 - if the result is finite and nonzero and errno == ERANGE then
171 OverflowError is raised
172
173 The last rule is used to catch overflow on platforms which follow
174 C89 but for which HUGE_VAL is not an infinity.
175
176 For most two-argument functions (copysign, fmod, hypot, atan2)
177 these rules are enough to ensure that Python's functions behave as
178 specified in 'Annex F' of the C99 standard, with the 'invalid' and
179 'divide-by-zero' floating-point exceptions mapping to Python's
180 ValueError and the 'overflow' floating-point exception mapping to
181 OverflowError.
182*/
183
Barry Warsaw8b43b191996-12-09 22:32:36 +0000184static PyObject *
Neal Norwitz45e230a2006-11-19 21:26:53 +0000185math_2(PyObject *args, double (*func) (double, double), char *funcname)
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000186{
Neal Norwitz45e230a2006-11-19 21:26:53 +0000187 PyObject *ox, *oy;
Christian Heimes6f341092008-04-18 23:13:07 +0000188 double x, y, r;
Neal Norwitz45e230a2006-11-19 21:26:53 +0000189 if (! PyArg_UnpackTuple(args, funcname, 2, 2, &ox, &oy))
190 return NULL;
191 x = PyFloat_AsDouble(ox);
192 y = PyFloat_AsDouble(oy);
193 if ((x == -1.0 || y == -1.0) && PyErr_Occurred())
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000194 return NULL;
195 errno = 0;
Christian Heimes6f341092008-04-18 23:13:07 +0000196 PyFPE_START_PROTECT("in math_2", return 0);
197 r = (*func)(x, y);
198 PyFPE_END_PROTECT(r);
199 if (Py_IS_NAN(r)) {
200 if (!Py_IS_NAN(x) && !Py_IS_NAN(y))
201 errno = EDOM;
202 else
203 errno = 0;
204 }
205 else if (Py_IS_INFINITY(r)) {
206 if (Py_IS_FINITE(x) && Py_IS_FINITE(y))
207 errno = ERANGE;
208 else
209 errno = 0;
210 }
211 if (errno && is_error(r))
Tim Peters1d120612000-10-12 06:10:25 +0000212 return NULL;
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000213 else
Christian Heimes6f341092008-04-18 23:13:07 +0000214 return PyFloat_FromDouble(r);
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000215}
216
Christian Heimes6f341092008-04-18 23:13:07 +0000217#define FUNC1(funcname, func, can_overflow, docstring) \
Fred Drake40c48682000-07-03 18:11:56 +0000218 static PyObject * math_##funcname(PyObject *self, PyObject *args) { \
Christian Heimes6f341092008-04-18 23:13:07 +0000219 return math_1(args, func, can_overflow); \
Guido van Rossumc6e22901998-12-04 19:26:43 +0000220 }\
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000221 PyDoc_STRVAR(math_##funcname##_doc, docstring);
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000222
Fred Drake40c48682000-07-03 18:11:56 +0000223#define FUNC2(funcname, func, docstring) \
224 static PyObject * math_##funcname(PyObject *self, PyObject *args) { \
Neal Norwitz45e230a2006-11-19 21:26:53 +0000225 return math_2(args, func, #funcname); \
Guido van Rossumc6e22901998-12-04 19:26:43 +0000226 }\
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000227 PyDoc_STRVAR(math_##funcname##_doc, docstring);
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000228
Christian Heimes6f341092008-04-18 23:13:07 +0000229FUNC1(acos, acos, 0,
Tim Petersfe71f812001-08-07 22:10:00 +0000230 "acos(x)\n\nReturn the arc cosine (measured in radians) of x.")
Christian Heimes6f341092008-04-18 23:13:07 +0000231FUNC1(acosh, acosh, 0,
232 "acosh(x)\n\nReturn the hyperbolic arc cosine (measured in radians) of x.")
233FUNC1(asin, asin, 0,
Tim Petersfe71f812001-08-07 22:10:00 +0000234 "asin(x)\n\nReturn the arc sine (measured in radians) of x.")
Christian Heimes6f341092008-04-18 23:13:07 +0000235FUNC1(asinh, asinh, 0,
236 "asinh(x)\n\nReturn the hyperbolic arc sine (measured in radians) of x.")
237FUNC1(atan, atan, 0,
Tim Petersfe71f812001-08-07 22:10:00 +0000238 "atan(x)\n\nReturn the arc tangent (measured in radians) of x.")
Martin v. Löwis387c5472001-09-06 08:16:17 +0000239FUNC2(atan2, atan2,
Tim Petersfe71f812001-08-07 22:10:00 +0000240 "atan2(y, x)\n\nReturn the arc tangent (measured in radians) of y/x.\n"
241 "Unlike atan(y/x), the signs of both x and y are considered.")
Christian Heimes6f341092008-04-18 23:13:07 +0000242FUNC1(atanh, atanh, 0,
243 "atanh(x)\n\nReturn the hyperbolic arc tangent (measured in radians) of x.")
244FUNC1(ceil, ceil, 0,
Jeffrey Yasskin9871d8f2008-01-05 08:47:13 +0000245 "ceil(x)\n\nReturn the ceiling of x as a float.\n"
246 "This is the smallest integral value >= x.")
Christian Heimeseebb79c2008-01-03 22:32:26 +0000247FUNC2(copysign, copysign,
Christian Heimes6f341092008-04-18 23:13:07 +0000248 "copysign(x,y)\n\nReturn x with the sign of y.")
249FUNC1(cos, cos, 0,
250 "cos(x)\n\nReturn the cosine of x (measured in radians).")
251FUNC1(cosh, cosh, 1,
252 "cosh(x)\n\nReturn the hyperbolic cosine of x.")
253FUNC1(exp, exp, 1,
Guido van Rossumc6e22901998-12-04 19:26:43 +0000254 "exp(x)\n\nReturn e raised to the power of x.")
Christian Heimes6f341092008-04-18 23:13:07 +0000255FUNC1(fabs, fabs, 0,
Tim Petersfe71f812001-08-07 22:10:00 +0000256 "fabs(x)\n\nReturn the absolute value of the float x.")
Christian Heimes6f341092008-04-18 23:13:07 +0000257FUNC1(floor, floor, 0,
Jeffrey Yasskin9871d8f2008-01-05 08:47:13 +0000258 "floor(x)\n\nReturn the floor of x as a float.\n"
259 "This is the largest integral value <= x.")
Christian Heimes6f341092008-04-18 23:13:07 +0000260FUNC1(log1p, log1p, 1,
261 "log1p(x)\n\nReturn the natural logarithm of 1+x (base e).\n\
262 The result is computed in a way which is accurate for x near zero.")
263FUNC1(sin, sin, 0,
Tim Petersfe71f812001-08-07 22:10:00 +0000264 "sin(x)\n\nReturn the sine of x (measured in radians).")
Christian Heimes6f341092008-04-18 23:13:07 +0000265FUNC1(sinh, sinh, 1,
Guido van Rossumc6e22901998-12-04 19:26:43 +0000266 "sinh(x)\n\nReturn the hyperbolic sine of x.")
Christian Heimes6f341092008-04-18 23:13:07 +0000267FUNC1(sqrt, sqrt, 0,
Guido van Rossumc6e22901998-12-04 19:26:43 +0000268 "sqrt(x)\n\nReturn the square root of x.")
Christian Heimes6f341092008-04-18 23:13:07 +0000269FUNC1(tan, tan, 0,
Tim Petersfe71f812001-08-07 22:10:00 +0000270 "tan(x)\n\nReturn the tangent of x (measured in radians).")
Christian Heimes6f341092008-04-18 23:13:07 +0000271FUNC1(tanh, tanh, 0,
Guido van Rossumc6e22901998-12-04 19:26:43 +0000272 "tanh(x)\n\nReturn the hyperbolic tangent of x.")
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000273
Barry Warsaw8b43b191996-12-09 22:32:36 +0000274static PyObject *
Jeffrey Yasskinca2b69f2008-02-01 06:22:46 +0000275math_trunc(PyObject *self, PyObject *number)
276{
Jeffrey Yasskinca2b69f2008-02-01 06:22:46 +0000277 return PyObject_CallMethod(number, "__trunc__", NULL);
278}
279
280PyDoc_STRVAR(math_trunc_doc,
281"trunc(x:Real) -> Integral\n"
282"\n"
Raymond Hettingerfe424f72008-02-02 05:24:44 +0000283"Truncates x to the nearest Integral toward 0. Uses the __trunc__ magic method.");
Jeffrey Yasskinca2b69f2008-02-01 06:22:46 +0000284
285static PyObject *
Neal Norwitz45e230a2006-11-19 21:26:53 +0000286math_frexp(PyObject *self, PyObject *arg)
Guido van Rossumd18ad581991-10-24 14:57:21 +0000287{
Guido van Rossumd18ad581991-10-24 14:57:21 +0000288 int i;
Neal Norwitz45e230a2006-11-19 21:26:53 +0000289 double x = PyFloat_AsDouble(arg);
290 if (x == -1.0 && PyErr_Occurred())
Guido van Rossumd18ad581991-10-24 14:57:21 +0000291 return NULL;
Christian Heimes6f341092008-04-18 23:13:07 +0000292 /* deal with special cases directly, to sidestep platform
293 differences */
294 if (Py_IS_NAN(x) || Py_IS_INFINITY(x) || !x) {
295 i = 0;
296 }
297 else {
298 PyFPE_START_PROTECT("in math_frexp", return 0);
299 x = frexp(x, &i);
300 PyFPE_END_PROTECT(x);
301 }
302 return Py_BuildValue("(di)", x, i);
Guido van Rossumd18ad581991-10-24 14:57:21 +0000303}
304
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000305PyDoc_STRVAR(math_frexp_doc,
Tim Peters63c94532001-09-04 23:17:42 +0000306"frexp(x)\n"
307"\n"
308"Return the mantissa and exponent of x, as pair (m, e).\n"
309"m is a float and e is an int, such that x = m * 2.**e.\n"
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000310"If x is 0, m and e are both 0. Else 0.5 <= abs(m) < 1.0.");
Guido van Rossumc6e22901998-12-04 19:26:43 +0000311
Barry Warsaw8b43b191996-12-09 22:32:36 +0000312static PyObject *
Fred Drake40c48682000-07-03 18:11:56 +0000313math_ldexp(PyObject *self, PyObject *args)
Guido van Rossumd18ad581991-10-24 14:57:21 +0000314{
Christian Heimes6f341092008-04-18 23:13:07 +0000315 double x, r;
Guido van Rossumc5545052000-05-08 14:29:38 +0000316 int exp;
Fred Drake40c48682000-07-03 18:11:56 +0000317 if (! PyArg_ParseTuple(args, "di:ldexp", &x, &exp))
Guido van Rossumd18ad581991-10-24 14:57:21 +0000318 return NULL;
319 errno = 0;
Christian Heimes6f341092008-04-18 23:13:07 +0000320 PyFPE_START_PROTECT("in math_ldexp", return 0)
321 r = ldexp(x, exp);
322 PyFPE_END_PROTECT(r)
323 if (Py_IS_FINITE(x) && Py_IS_INFINITY(r))
324 errno = ERANGE;
325 /* Windows MSVC8 sets errno = EDOM on ldexp(NaN, i);
326 we unset it to avoid raising a ValueError here. */
327 if (errno == EDOM)
328 errno = 0;
329 if (errno && is_error(r))
Tim Peters1d120612000-10-12 06:10:25 +0000330 return NULL;
Guido van Rossumd18ad581991-10-24 14:57:21 +0000331 else
Christian Heimes6f341092008-04-18 23:13:07 +0000332 return PyFloat_FromDouble(r);
Guido van Rossumd18ad581991-10-24 14:57:21 +0000333}
334
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000335PyDoc_STRVAR(math_ldexp_doc,
336"ldexp(x, i) -> x * (2**i)");
Guido van Rossumc6e22901998-12-04 19:26:43 +0000337
Barry Warsaw8b43b191996-12-09 22:32:36 +0000338static PyObject *
Neal Norwitz45e230a2006-11-19 21:26:53 +0000339math_modf(PyObject *self, PyObject *arg)
Guido van Rossumd18ad581991-10-24 14:57:21 +0000340{
Neal Norwitz45e230a2006-11-19 21:26:53 +0000341 double y, x = PyFloat_AsDouble(arg);
342 if (x == -1.0 && PyErr_Occurred())
Guido van Rossumd18ad581991-10-24 14:57:21 +0000343 return NULL;
Mark Dickinsonb2f70902008-04-20 01:39:24 +0000344 /* some platforms don't do the right thing for NaNs and
345 infinities, so we take care of special cases directly. */
346 if (!Py_IS_FINITE(x)) {
347 if (Py_IS_INFINITY(x))
348 return Py_BuildValue("(dd)", copysign(0., x), x);
349 else if (Py_IS_NAN(x))
350 return Py_BuildValue("(dd)", x, x);
351 }
352
Guido van Rossumd18ad581991-10-24 14:57:21 +0000353 errno = 0;
Christian Heimes6f341092008-04-18 23:13:07 +0000354 PyFPE_START_PROTECT("in math_modf", return 0);
Guido van Rossumd18ad581991-10-24 14:57:21 +0000355 x = modf(x, &y);
Christian Heimes6f341092008-04-18 23:13:07 +0000356 PyFPE_END_PROTECT(x);
357 return Py_BuildValue("(dd)", x, y);
Guido van Rossumd18ad581991-10-24 14:57:21 +0000358}
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000359
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000360PyDoc_STRVAR(math_modf_doc,
Tim Peters63c94532001-09-04 23:17:42 +0000361"modf(x)\n"
362"\n"
363"Return the fractional and integer parts of x. Both results carry the sign\n"
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000364"of x. The integer part is returned as a real.");
Guido van Rossumc6e22901998-12-04 19:26:43 +0000365
Tim Peters78526162001-09-05 00:53:45 +0000366/* A decent logarithm is easy to compute even for huge longs, but libm can't
367 do that by itself -- loghelper can. func is log or log10, and name is
368 "log" or "log10". Note that overflow isn't possible: a long can contain
369 no more than INT_MAX * SHIFT bits, so has value certainly less than
370 2**(2**64 * 2**16) == 2**2**80, and log2 of that is 2**80, which is
371 small enough to fit in an IEEE single. log and log10 are even smaller.
372*/
373
374static PyObject*
Neal Norwitz45e230a2006-11-19 21:26:53 +0000375loghelper(PyObject* arg, double (*func)(double), char *funcname)
Tim Peters78526162001-09-05 00:53:45 +0000376{
Tim Peters78526162001-09-05 00:53:45 +0000377 /* If it is long, do it ourselves. */
378 if (PyLong_Check(arg)) {
379 double x;
380 int e;
381 x = _PyLong_AsScaledDouble(arg, &e);
382 if (x <= 0.0) {
383 PyErr_SetString(PyExc_ValueError,
384 "math domain error");
385 return NULL;
386 }
Christian Heimes543cabc2008-01-25 14:54:23 +0000387 /* Value is ~= x * 2**(e*PyLong_SHIFT), so the log ~=
388 log(x) + log(2) * e * PyLong_SHIFT.
389 CAUTION: e*PyLong_SHIFT may overflow using int arithmetic,
Tim Peters78526162001-09-05 00:53:45 +0000390 so force use of double. */
Christian Heimes543cabc2008-01-25 14:54:23 +0000391 x = func(x) + (e * (double)PyLong_SHIFT) * func(2.0);
Tim Peters78526162001-09-05 00:53:45 +0000392 return PyFloat_FromDouble(x);
393 }
394
395 /* Else let libm handle it by itself. */
Christian Heimes6f341092008-04-18 23:13:07 +0000396 return math_1(arg, func, 0);
Tim Peters78526162001-09-05 00:53:45 +0000397}
398
399static PyObject *
400math_log(PyObject *self, PyObject *args)
401{
Raymond Hettinger866964c2002-12-14 19:51:34 +0000402 PyObject *arg;
403 PyObject *base = NULL;
404 PyObject *num, *den;
405 PyObject *ans;
Raymond Hettinger866964c2002-12-14 19:51:34 +0000406
Raymond Hettingerea3fdf42002-12-29 16:33:45 +0000407 if (!PyArg_UnpackTuple(args, "log", 1, 2, &arg, &base))
Raymond Hettinger866964c2002-12-14 19:51:34 +0000408 return NULL;
Raymond Hettinger866964c2002-12-14 19:51:34 +0000409
Neal Norwitz45e230a2006-11-19 21:26:53 +0000410 num = loghelper(arg, log, "log");
411 if (num == NULL || base == NULL)
412 return num;
Raymond Hettinger866964c2002-12-14 19:51:34 +0000413
Neal Norwitz45e230a2006-11-19 21:26:53 +0000414 den = loghelper(base, log, "log");
Raymond Hettinger866964c2002-12-14 19:51:34 +0000415 if (den == NULL) {
416 Py_DECREF(num);
417 return NULL;
418 }
419
420 ans = PyNumber_Divide(num, den);
421 Py_DECREF(num);
422 Py_DECREF(den);
423 return ans;
Tim Peters78526162001-09-05 00:53:45 +0000424}
425
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000426PyDoc_STRVAR(math_log_doc,
Raymond Hettinger866964c2002-12-14 19:51:34 +0000427"log(x[, base]) -> the logarithm of x to the given base.\n\
428If the base not specified, returns the natural logarithm (base e) of x.");
Tim Peters78526162001-09-05 00:53:45 +0000429
430static PyObject *
Neal Norwitz45e230a2006-11-19 21:26:53 +0000431math_log10(PyObject *self, PyObject *arg)
Tim Peters78526162001-09-05 00:53:45 +0000432{
Neal Norwitz45e230a2006-11-19 21:26:53 +0000433 return loghelper(arg, log10, "log10");
Tim Peters78526162001-09-05 00:53:45 +0000434}
435
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000436PyDoc_STRVAR(math_log10_doc,
437"log10(x) -> the base 10 logarithm of x.");
Tim Peters78526162001-09-05 00:53:45 +0000438
Christian Heimes6f341092008-04-18 23:13:07 +0000439static PyObject *
440math_fmod(PyObject *self, PyObject *args)
441{
442 PyObject *ox, *oy;
443 double r, x, y;
444 if (! PyArg_UnpackTuple(args, "fmod", 2, 2, &ox, &oy))
445 return NULL;
446 x = PyFloat_AsDouble(ox);
447 y = PyFloat_AsDouble(oy);
448 if ((x == -1.0 || y == -1.0) && PyErr_Occurred())
449 return NULL;
450 /* fmod(x, +/-Inf) returns x for finite x. */
451 if (Py_IS_INFINITY(y) && Py_IS_FINITE(x))
452 return PyFloat_FromDouble(x);
453 errno = 0;
454 PyFPE_START_PROTECT("in math_fmod", return 0);
455 r = fmod(x, y);
456 PyFPE_END_PROTECT(r);
457 if (Py_IS_NAN(r)) {
458 if (!Py_IS_NAN(x) && !Py_IS_NAN(y))
459 errno = EDOM;
460 else
461 errno = 0;
462 }
463 if (errno && is_error(r))
464 return NULL;
465 else
466 return PyFloat_FromDouble(r);
467}
468
469PyDoc_STRVAR(math_fmod_doc,
470"fmod(x,y)\n\nReturn fmod(x, y), according to platform C."
471" x % y may differ.");
472
473static PyObject *
474math_hypot(PyObject *self, PyObject *args)
475{
476 PyObject *ox, *oy;
477 double r, x, y;
478 if (! PyArg_UnpackTuple(args, "hypot", 2, 2, &ox, &oy))
479 return NULL;
480 x = PyFloat_AsDouble(ox);
481 y = PyFloat_AsDouble(oy);
482 if ((x == -1.0 || y == -1.0) && PyErr_Occurred())
483 return NULL;
484 /* hypot(x, +/-Inf) returns Inf, even if x is a NaN. */
485 if (Py_IS_INFINITY(x))
486 return PyFloat_FromDouble(fabs(x));
487 if (Py_IS_INFINITY(y))
488 return PyFloat_FromDouble(fabs(y));
489 errno = 0;
490 PyFPE_START_PROTECT("in math_hypot", return 0);
491 r = hypot(x, y);
492 PyFPE_END_PROTECT(r);
493 if (Py_IS_NAN(r)) {
494 if (!Py_IS_NAN(x) && !Py_IS_NAN(y))
495 errno = EDOM;
496 else
497 errno = 0;
498 }
499 else if (Py_IS_INFINITY(r)) {
500 if (Py_IS_FINITE(x) && Py_IS_FINITE(y))
501 errno = ERANGE;
502 else
503 errno = 0;
504 }
505 if (errno && is_error(r))
506 return NULL;
507 else
508 return PyFloat_FromDouble(r);
509}
510
511PyDoc_STRVAR(math_hypot_doc,
512"hypot(x,y)\n\nReturn the Euclidean distance, sqrt(x*x + y*y).");
513
514/* pow can't use math_2, but needs its own wrapper: the problem is
515 that an infinite result can arise either as a result of overflow
516 (in which case OverflowError should be raised) or as a result of
517 e.g. 0.**-5. (for which ValueError needs to be raised.)
518*/
519
520static PyObject *
521math_pow(PyObject *self, PyObject *args)
522{
523 PyObject *ox, *oy;
524 double r, x, y;
Mark Dickinsone941d972008-04-19 18:51:48 +0000525 int y_is_odd;
Christian Heimes6f341092008-04-18 23:13:07 +0000526
527 if (! PyArg_UnpackTuple(args, "pow", 2, 2, &ox, &oy))
528 return NULL;
529 x = PyFloat_AsDouble(ox);
530 y = PyFloat_AsDouble(oy);
531 if ((x == -1.0 || y == -1.0) && PyErr_Occurred())
532 return NULL;
Mark Dickinsona1293eb2008-04-19 19:41:52 +0000533
534 /* deal directly with various special cases, to cope with problems on
535 various platforms whose semantics don't exactly match C99 */
536
537 /* 1**x, x**0, and (-1)**(+-infinity) return 1., even if x is NaN or
538 an infinity. */
539 if (x == 1. || y == 0. || (x == -1. && Py_IS_INFINITY(y)))
Christian Heimes6f341092008-04-18 23:13:07 +0000540 return PyFloat_FromDouble(1.);
Mark Dickinsona1293eb2008-04-19 19:41:52 +0000541 /* otherwise, return a NaN if either input was a NaN */
542 if (Py_IS_NAN(x))
543 return PyFloat_FromDouble(x);
544 if (Py_IS_NAN(y))
545 return PyFloat_FromDouble(y);
Mark Dickinsone941d972008-04-19 18:51:48 +0000546 /* inf ** (nonzero, non-NaN) is one of +-0, +-infinity */
547 if (Py_IS_INFINITY(x) && !Py_IS_NAN(y)) {
548 y_is_odd = Py_IS_FINITE(y) && fmod(fabs(y), 2.0) == 1.0;
549 if (y > 0.)
550 r = y_is_odd ? x : fabs(x);
551 else
552 r = y_is_odd ? copysign(0., x) : 0.;
553 return PyFloat_FromDouble(r);
554 }
555
Christian Heimes6f341092008-04-18 23:13:07 +0000556 errno = 0;
557 PyFPE_START_PROTECT("in math_pow", return 0);
558 r = pow(x, y);
559 PyFPE_END_PROTECT(r);
560 if (Py_IS_NAN(r)) {
Mark Dickinsona1293eb2008-04-19 19:41:52 +0000561 errno = EDOM;
Christian Heimes6f341092008-04-18 23:13:07 +0000562 }
563 /* an infinite result arises either from:
564
565 (A) (+/-0.)**negative,
566 (B) overflow of x**y with both x and y finite (and x nonzero)
567 (C) (+/-inf)**positive, or
568 (D) x**inf with |x| > 1, or x**-inf with |x| < 1.
569
570 In case (A) we want ValueError to be raised. In case (B)
571 OverflowError should be raised. In cases (C) and (D) the infinite
572 result should be returned.
573 */
574 else if (Py_IS_INFINITY(r)) {
575 if (x == 0.)
576 errno = EDOM;
577 else if (Py_IS_FINITE(x) && Py_IS_FINITE(y))
578 errno = ERANGE;
579 else
580 errno = 0;
581 }
582
583 if (errno && is_error(r))
584 return NULL;
585 else
586 return PyFloat_FromDouble(r);
587}
588
589PyDoc_STRVAR(math_pow_doc,
590"pow(x,y)\n\nReturn x**y (x to the power of y).");
591
Christian Heimese2ca4242008-01-03 20:23:15 +0000592static const double degToRad = Py_MATH_PI / 180.0;
593static const double radToDeg = 180.0 / Py_MATH_PI;
Raymond Hettingerd6f22672002-05-13 03:56:10 +0000594
595static PyObject *
Neal Norwitz45e230a2006-11-19 21:26:53 +0000596math_degrees(PyObject *self, PyObject *arg)
Raymond Hettingerd6f22672002-05-13 03:56:10 +0000597{
Neal Norwitz45e230a2006-11-19 21:26:53 +0000598 double x = PyFloat_AsDouble(arg);
599 if (x == -1.0 && PyErr_Occurred())
Raymond Hettingerd6f22672002-05-13 03:56:10 +0000600 return NULL;
Christian Heimese2ca4242008-01-03 20:23:15 +0000601 return PyFloat_FromDouble(x * radToDeg);
Raymond Hettingerd6f22672002-05-13 03:56:10 +0000602}
603
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000604PyDoc_STRVAR(math_degrees_doc,
605"degrees(x) -> converts angle x from radians to degrees");
Raymond Hettingerd6f22672002-05-13 03:56:10 +0000606
607static PyObject *
Neal Norwitz45e230a2006-11-19 21:26:53 +0000608math_radians(PyObject *self, PyObject *arg)
Raymond Hettingerd6f22672002-05-13 03:56:10 +0000609{
Neal Norwitz45e230a2006-11-19 21:26:53 +0000610 double x = PyFloat_AsDouble(arg);
611 if (x == -1.0 && PyErr_Occurred())
Raymond Hettingerd6f22672002-05-13 03:56:10 +0000612 return NULL;
613 return PyFloat_FromDouble(x * degToRad);
614}
615
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000616PyDoc_STRVAR(math_radians_doc,
617"radians(x) -> converts angle x from degrees to radians");
Tim Peters78526162001-09-05 00:53:45 +0000618
Christian Heimese2ca4242008-01-03 20:23:15 +0000619static PyObject *
620math_isnan(PyObject *self, PyObject *arg)
621{
622 double x = PyFloat_AsDouble(arg);
623 if (x == -1.0 && PyErr_Occurred())
624 return NULL;
625 return PyBool_FromLong((long)Py_IS_NAN(x));
626}
627
628PyDoc_STRVAR(math_isnan_doc,
629"isnan(x) -> bool\n\
630Checks if float x is not a number (NaN)");
631
632static PyObject *
633math_isinf(PyObject *self, PyObject *arg)
634{
635 double x = PyFloat_AsDouble(arg);
636 if (x == -1.0 && PyErr_Occurred())
637 return NULL;
638 return PyBool_FromLong((long)Py_IS_INFINITY(x));
639}
640
641PyDoc_STRVAR(math_isinf_doc,
642"isinf(x) -> bool\n\
643Checks if float x is infinite (positive or negative)");
644
Barry Warsaw8b43b191996-12-09 22:32:36 +0000645static PyMethodDef math_methods[] = {
Neal Norwitz45e230a2006-11-19 21:26:53 +0000646 {"acos", math_acos, METH_O, math_acos_doc},
Christian Heimes6f341092008-04-18 23:13:07 +0000647 {"acosh", math_acosh, METH_O, math_acosh_doc},
Neal Norwitz45e230a2006-11-19 21:26:53 +0000648 {"asin", math_asin, METH_O, math_asin_doc},
Christian Heimes6f341092008-04-18 23:13:07 +0000649 {"asinh", math_asinh, METH_O, math_asinh_doc},
Neal Norwitz45e230a2006-11-19 21:26:53 +0000650 {"atan", math_atan, METH_O, math_atan_doc},
Fred Drake40c48682000-07-03 18:11:56 +0000651 {"atan2", math_atan2, METH_VARARGS, math_atan2_doc},
Christian Heimes6f341092008-04-18 23:13:07 +0000652 {"atanh", math_atanh, METH_O, math_atanh_doc},
Neal Norwitz45e230a2006-11-19 21:26:53 +0000653 {"ceil", math_ceil, METH_O, math_ceil_doc},
Christian Heimeseebb79c2008-01-03 22:32:26 +0000654 {"copysign", math_copysign, METH_VARARGS, math_copysign_doc},
Neal Norwitz45e230a2006-11-19 21:26:53 +0000655 {"cos", math_cos, METH_O, math_cos_doc},
656 {"cosh", math_cosh, METH_O, math_cosh_doc},
657 {"degrees", math_degrees, METH_O, math_degrees_doc},
658 {"exp", math_exp, METH_O, math_exp_doc},
659 {"fabs", math_fabs, METH_O, math_fabs_doc},
660 {"floor", math_floor, METH_O, math_floor_doc},
Fred Drake40c48682000-07-03 18:11:56 +0000661 {"fmod", math_fmod, METH_VARARGS, math_fmod_doc},
Neal Norwitz45e230a2006-11-19 21:26:53 +0000662 {"frexp", math_frexp, METH_O, math_frexp_doc},
Fred Drake40c48682000-07-03 18:11:56 +0000663 {"hypot", math_hypot, METH_VARARGS, math_hypot_doc},
Christian Heimese2ca4242008-01-03 20:23:15 +0000664 {"isinf", math_isinf, METH_O, math_isinf_doc},
665 {"isnan", math_isnan, METH_O, math_isnan_doc},
Fred Drake40c48682000-07-03 18:11:56 +0000666 {"ldexp", math_ldexp, METH_VARARGS, math_ldexp_doc},
667 {"log", math_log, METH_VARARGS, math_log_doc},
Christian Heimes6f341092008-04-18 23:13:07 +0000668 {"log1p", math_log1p, METH_O, math_log1p_doc},
Neal Norwitz45e230a2006-11-19 21:26:53 +0000669 {"log10", math_log10, METH_O, math_log10_doc},
670 {"modf", math_modf, METH_O, math_modf_doc},
Fred Drake40c48682000-07-03 18:11:56 +0000671 {"pow", math_pow, METH_VARARGS, math_pow_doc},
Neal Norwitz45e230a2006-11-19 21:26:53 +0000672 {"radians", math_radians, METH_O, math_radians_doc},
673 {"sin", math_sin, METH_O, math_sin_doc},
674 {"sinh", math_sinh, METH_O, math_sinh_doc},
675 {"sqrt", math_sqrt, METH_O, math_sqrt_doc},
676 {"tan", math_tan, METH_O, math_tan_doc},
677 {"tanh", math_tanh, METH_O, math_tanh_doc},
Jeffrey Yasskinca2b69f2008-02-01 06:22:46 +0000678 {"trunc", math_trunc, METH_O, math_trunc_doc},
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000679 {NULL, NULL} /* sentinel */
680};
681
Guido van Rossumc6e22901998-12-04 19:26:43 +0000682
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000683PyDoc_STRVAR(module_doc,
Tim Peters63c94532001-09-04 23:17:42 +0000684"This module is always available. It provides access to the\n"
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +0000685"mathematical functions defined by the C standard.");
Guido van Rossumc6e22901998-12-04 19:26:43 +0000686
Mark Hammondfe51c6d2002-08-02 02:27:13 +0000687PyMODINIT_FUNC
Thomas Woutersf3f33dc2000-07-21 06:00:07 +0000688initmath(void)
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000689{
Christian Heimes6f341092008-04-18 23:13:07 +0000690 PyObject *m;
Tim Petersfe71f812001-08-07 22:10:00 +0000691
Guido van Rossumc6e22901998-12-04 19:26:43 +0000692 m = Py_InitModule3("math", math_methods, module_doc);
Neal Norwitz1ac754f2006-01-19 06:09:39 +0000693 if (m == NULL)
694 goto finally;
Barry Warsawfc93f751996-12-17 00:47:03 +0000695
Christian Heimes6f341092008-04-18 23:13:07 +0000696 PyModule_AddObject(m, "pi", PyFloat_FromDouble(Py_MATH_PI));
697 PyModule_AddObject(m, "e", PyFloat_FromDouble(Py_MATH_E));
Barry Warsawfc93f751996-12-17 00:47:03 +0000698
Christian Heimes6f341092008-04-18 23:13:07 +0000699 finally:
Barry Warsaw9bfd2bf2000-09-01 09:01:32 +0000700 return;
Guido van Rossum85a5fbb1990-10-14 12:07:46 +0000701}