blob: cba42a7257a08378355ba20abe6b3dd1929bd738 [file] [log] [blame]
Guido van Rossum71aa32f1996-01-12 01:34:57 +00001/* Complex math module */
2
3/* much code borrowed from mathmodule.c */
4
Roger E. Masse24070ca1996-12-09 22:59:53 +00005#include "Python.h"
Mark Dickinsonf3718592009-12-21 15:27:41 +00006#include "_math.h"
Christian Heimes53876d92008-04-19 00:31:39 +00007/* we need DBL_MAX, DBL_MIN, DBL_EPSILON, DBL_MANT_DIG and FLT_RADIX from
8 float.h. We assume that FLT_RADIX is either 2 or 16. */
9#include <float.h>
Guido van Rossum71aa32f1996-01-12 01:34:57 +000010
Brett Cannonb0fc4902014-10-14 17:37:02 -040011#include "clinic/cmathmodule.c.h"
12/*[clinic input]
Brett Cannonb0fc4902014-10-14 17:37:02 -040013module cmath
14[clinic start generated code]*/
Serhiy Storchaka1009bf12015-04-03 23:53:51 +030015/*[clinic end generated code: output=da39a3ee5e6b4b0d input=308d6839f4a46333]*/
Brett Cannonb0fc4902014-10-14 17:37:02 -040016
17/*[python input]
18class Py_complex_protected_converter(Py_complex_converter):
19 def modify(self):
20 return 'errno = 0; PyFPE_START_PROTECT("complex function", goto exit);'
21
22
23class Py_complex_protected_return_converter(CReturnConverter):
24 type = "Py_complex"
25
26 def render(self, function, data):
27 self.declare(data)
28 data.return_conversion.append("""
29PyFPE_END_PROTECT(_return_value);
Serhiy Storchakaebe95fd2016-06-09 16:02:15 +030030if (errno == EDOM) {
Brett Cannonb0fc4902014-10-14 17:37:02 -040031 PyErr_SetString(PyExc_ValueError, "math domain error");
32 goto exit;
Serhiy Storchakaebe95fd2016-06-09 16:02:15 +030033}
34else if (errno == ERANGE) {
Brett Cannonb0fc4902014-10-14 17:37:02 -040035 PyErr_SetString(PyExc_OverflowError, "math range error");
36 goto exit;
Serhiy Storchakaebe95fd2016-06-09 16:02:15 +030037}
38else {
Brett Cannonb0fc4902014-10-14 17:37:02 -040039 return_value = PyComplex_FromCComplex(_return_value);
Serhiy Storchakaebe95fd2016-06-09 16:02:15 +030040}
Brett Cannonb0fc4902014-10-14 17:37:02 -040041""".strip())
42[python start generated code]*/
Serhiy Storchakaebe95fd2016-06-09 16:02:15 +030043/*[python end generated code: output=da39a3ee5e6b4b0d input=345daa075b1028e7]*/
Brett Cannonb0fc4902014-10-14 17:37:02 -040044
Christian Heimes53876d92008-04-19 00:31:39 +000045#if (FLT_RADIX != 2 && FLT_RADIX != 16)
46#error "Modules/cmathmodule.c expects FLT_RADIX to be 2 or 16"
Guido van Rossum71aa32f1996-01-12 01:34:57 +000047#endif
48
Christian Heimes53876d92008-04-19 00:31:39 +000049#ifndef M_LN2
50#define M_LN2 (0.6931471805599453094) /* natural log of 2 */
51#endif
Guido van Rossum71aa32f1996-01-12 01:34:57 +000052
Christian Heimes53876d92008-04-19 00:31:39 +000053#ifndef M_LN10
54#define M_LN10 (2.302585092994045684) /* natural log of 10 */
55#endif
56
57/*
58 CM_LARGE_DOUBLE is used to avoid spurious overflow in the sqrt, log,
59 inverse trig and inverse hyperbolic trig functions. Its log is used in the
Ezio Melotti13925002011-03-16 11:05:33 +020060 evaluation of exp, cos, cosh, sin, sinh, tan, and tanh to avoid unnecessary
Christian Heimes53876d92008-04-19 00:31:39 +000061 overflow.
62 */
63
64#define CM_LARGE_DOUBLE (DBL_MAX/4.)
65#define CM_SQRT_LARGE_DOUBLE (sqrt(CM_LARGE_DOUBLE))
66#define CM_LOG_LARGE_DOUBLE (log(CM_LARGE_DOUBLE))
67#define CM_SQRT_DBL_MIN (sqrt(DBL_MIN))
68
Antoine Pitrouf95a1b32010-05-09 15:52:27 +000069/*
Christian Heimes53876d92008-04-19 00:31:39 +000070 CM_SCALE_UP is an odd integer chosen such that multiplication by
71 2**CM_SCALE_UP is sufficient to turn a subnormal into a normal.
72 CM_SCALE_DOWN is (-(CM_SCALE_UP+1)/2). These scalings are used to compute
73 square roots accurately when the real and imaginary parts of the argument
74 are subnormal.
75*/
76
77#if FLT_RADIX==2
78#define CM_SCALE_UP (2*(DBL_MANT_DIG/2) + 1)
79#elif FLT_RADIX==16
80#define CM_SCALE_UP (4*DBL_MANT_DIG+1)
81#endif
82#define CM_SCALE_DOWN (-(CM_SCALE_UP+1)/2)
Guido van Rossum71aa32f1996-01-12 01:34:57 +000083
84/* forward declarations */
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +030085static Py_complex cmath_asinh_impl(PyObject *, Py_complex);
86static Py_complex cmath_atanh_impl(PyObject *, Py_complex);
87static Py_complex cmath_cosh_impl(PyObject *, Py_complex);
88static Py_complex cmath_sinh_impl(PyObject *, Py_complex);
89static Py_complex cmath_sqrt_impl(PyObject *, Py_complex);
90static Py_complex cmath_tanh_impl(PyObject *, Py_complex);
Raymond Hettingerb67ad7e2004-06-14 07:40:10 +000091static PyObject * math_error(void);
Guido van Rossum71aa32f1996-01-12 01:34:57 +000092
Christian Heimes53876d92008-04-19 00:31:39 +000093/* Code to deal with special values (infinities, NaNs, etc.). */
94
95/* special_type takes a double and returns an integer code indicating
96 the type of the double as follows:
97*/
98
99enum special_types {
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000100 ST_NINF, /* 0, negative infinity */
101 ST_NEG, /* 1, negative finite number (nonzero) */
102 ST_NZERO, /* 2, -0. */
103 ST_PZERO, /* 3, +0. */
104 ST_POS, /* 4, positive finite number (nonzero) */
105 ST_PINF, /* 5, positive infinity */
106 ST_NAN /* 6, Not a Number */
Christian Heimes53876d92008-04-19 00:31:39 +0000107};
108
109static enum special_types
110special_type(double d)
111{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000112 if (Py_IS_FINITE(d)) {
113 if (d != 0) {
114 if (copysign(1., d) == 1.)
115 return ST_POS;
116 else
117 return ST_NEG;
118 }
119 else {
120 if (copysign(1., d) == 1.)
121 return ST_PZERO;
122 else
123 return ST_NZERO;
124 }
125 }
126 if (Py_IS_NAN(d))
127 return ST_NAN;
128 if (copysign(1., d) == 1.)
129 return ST_PINF;
130 else
131 return ST_NINF;
Christian Heimes53876d92008-04-19 00:31:39 +0000132}
133
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000134#define SPECIAL_VALUE(z, table) \
135 if (!Py_IS_FINITE((z).real) || !Py_IS_FINITE((z).imag)) { \
136 errno = 0; \
137 return table[special_type((z).real)] \
138 [special_type((z).imag)]; \
139 }
Christian Heimes53876d92008-04-19 00:31:39 +0000140
141#define P Py_MATH_PI
142#define P14 0.25*Py_MATH_PI
143#define P12 0.5*Py_MATH_PI
144#define P34 0.75*Py_MATH_PI
Christian Heimesa342c012008-04-20 21:01:16 +0000145#define INF Py_HUGE_VAL
146#define N Py_NAN
Christian Heimes53876d92008-04-19 00:31:39 +0000147#define U -9.5426319407711027e33 /* unlikely value, used as placeholder */
148
149/* First, the C functions that do the real work. Each of the c_*
150 functions computes and returns the C99 Annex G recommended result
151 and also sets errno as follows: errno = 0 if no floating-point
152 exception is associated with the result; errno = EDOM if C99 Annex
153 G recommends raising divide-by-zero or invalid for this result; and
154 errno = ERANGE where the overflow floating-point signal should be
155 raised.
156*/
157
Christian Heimesa342c012008-04-20 21:01:16 +0000158static Py_complex acos_special_values[7][7];
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000159
Brett Cannonb0fc4902014-10-14 17:37:02 -0400160/*[clinic input]
161cmath.acos -> Py_complex_protected
162
163 z: Py_complex_protected
164 /
165
166Return the arc cosine of z.
167[clinic start generated code]*/
168
Tim Peters14e26402001-02-20 20:15:19 +0000169static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300170cmath_acos_impl(PyObject *module, Py_complex z)
171/*[clinic end generated code: output=40bd42853fd460ae input=bd6cbd78ae851927]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000172{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000173 Py_complex s1, s2, r;
Christian Heimes53876d92008-04-19 00:31:39 +0000174
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000175 SPECIAL_VALUE(z, acos_special_values);
Christian Heimes53876d92008-04-19 00:31:39 +0000176
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000177 if (fabs(z.real) > CM_LARGE_DOUBLE || fabs(z.imag) > CM_LARGE_DOUBLE) {
178 /* avoid unnecessary overflow for large arguments */
179 r.real = atan2(fabs(z.imag), z.real);
180 /* split into cases to make sure that the branch cut has the
181 correct continuity on systems with unsigned zeros */
182 if (z.real < 0.) {
183 r.imag = -copysign(log(hypot(z.real/2., z.imag/2.)) +
184 M_LN2*2., z.imag);
185 } else {
186 r.imag = copysign(log(hypot(z.real/2., z.imag/2.)) +
187 M_LN2*2., -z.imag);
188 }
189 } else {
190 s1.real = 1.-z.real;
191 s1.imag = -z.imag;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400192 s1 = cmath_sqrt_impl(module, s1);
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000193 s2.real = 1.+z.real;
194 s2.imag = z.imag;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400195 s2 = cmath_sqrt_impl(module, s2);
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000196 r.real = 2.*atan2(s1.real, s2.real);
197 r.imag = m_asinh(s2.real*s1.imag - s2.imag*s1.real);
198 }
199 errno = 0;
200 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000201}
202
Guido van Rossumc6e22901998-12-04 19:26:43 +0000203
Christian Heimesa342c012008-04-20 21:01:16 +0000204static Py_complex acosh_special_values[7][7];
Christian Heimes53876d92008-04-19 00:31:39 +0000205
Brett Cannonb0fc4902014-10-14 17:37:02 -0400206/*[clinic input]
207cmath.acosh = cmath.acos
208
Mark Dickinsoncc8617b2015-01-11 13:22:44 +0000209Return the inverse hyperbolic cosine of z.
Brett Cannonb0fc4902014-10-14 17:37:02 -0400210[clinic start generated code]*/
211
Tim Peters14e26402001-02-20 20:15:19 +0000212static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300213cmath_acosh_impl(PyObject *module, Py_complex z)
214/*[clinic end generated code: output=3e2454d4fcf404ca input=3f61bee7d703e53c]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000215{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000216 Py_complex s1, s2, r;
Christian Heimes53876d92008-04-19 00:31:39 +0000217
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000218 SPECIAL_VALUE(z, acosh_special_values);
Christian Heimes53876d92008-04-19 00:31:39 +0000219
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000220 if (fabs(z.real) > CM_LARGE_DOUBLE || fabs(z.imag) > CM_LARGE_DOUBLE) {
221 /* avoid unnecessary overflow for large arguments */
222 r.real = log(hypot(z.real/2., z.imag/2.)) + M_LN2*2.;
223 r.imag = atan2(z.imag, z.real);
224 } else {
225 s1.real = z.real - 1.;
226 s1.imag = z.imag;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400227 s1 = cmath_sqrt_impl(module, s1);
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000228 s2.real = z.real + 1.;
229 s2.imag = z.imag;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400230 s2 = cmath_sqrt_impl(module, s2);
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000231 r.real = m_asinh(s1.real*s2.real + s1.imag*s2.imag);
232 r.imag = 2.*atan2(s1.imag, s2.real);
233 }
234 errno = 0;
235 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000236}
237
Brett Cannonb0fc4902014-10-14 17:37:02 -0400238/*[clinic input]
239cmath.asin = cmath.acos
Guido van Rossumc6e22901998-12-04 19:26:43 +0000240
Brett Cannonb0fc4902014-10-14 17:37:02 -0400241Return the arc sine of z.
242[clinic start generated code]*/
Guido van Rossumc6e22901998-12-04 19:26:43 +0000243
Tim Peters14e26402001-02-20 20:15:19 +0000244static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300245cmath_asin_impl(PyObject *module, Py_complex z)
246/*[clinic end generated code: output=3b264cd1b16bf4e1 input=be0bf0cfdd5239c5]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000247{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000248 /* asin(z) = -i asinh(iz) */
249 Py_complex s, r;
250 s.real = -z.imag;
251 s.imag = z.real;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400252 s = cmath_asinh_impl(module, s);
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000253 r.real = s.imag;
254 r.imag = -s.real;
255 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000256}
257
Guido van Rossumc6e22901998-12-04 19:26:43 +0000258
Christian Heimesa342c012008-04-20 21:01:16 +0000259static Py_complex asinh_special_values[7][7];
Christian Heimes53876d92008-04-19 00:31:39 +0000260
Brett Cannonb0fc4902014-10-14 17:37:02 -0400261/*[clinic input]
262cmath.asinh = cmath.acos
263
Mark Dickinsoncc8617b2015-01-11 13:22:44 +0000264Return the inverse hyperbolic sine of z.
Brett Cannonb0fc4902014-10-14 17:37:02 -0400265[clinic start generated code]*/
266
Tim Peters14e26402001-02-20 20:15:19 +0000267static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300268cmath_asinh_impl(PyObject *module, Py_complex z)
269/*[clinic end generated code: output=733d8107841a7599 input=5c09448fcfc89a79]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000270{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000271 Py_complex s1, s2, r;
Christian Heimes53876d92008-04-19 00:31:39 +0000272
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000273 SPECIAL_VALUE(z, asinh_special_values);
Christian Heimes53876d92008-04-19 00:31:39 +0000274
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000275 if (fabs(z.real) > CM_LARGE_DOUBLE || fabs(z.imag) > CM_LARGE_DOUBLE) {
276 if (z.imag >= 0.) {
277 r.real = copysign(log(hypot(z.real/2., z.imag/2.)) +
278 M_LN2*2., z.real);
279 } else {
280 r.real = -copysign(log(hypot(z.real/2., z.imag/2.)) +
281 M_LN2*2., -z.real);
282 }
283 r.imag = atan2(z.imag, fabs(z.real));
284 } else {
285 s1.real = 1.+z.imag;
286 s1.imag = -z.real;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400287 s1 = cmath_sqrt_impl(module, s1);
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000288 s2.real = 1.-z.imag;
289 s2.imag = z.real;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400290 s2 = cmath_sqrt_impl(module, s2);
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000291 r.real = m_asinh(s1.real*s2.imag-s2.real*s1.imag);
292 r.imag = atan2(z.imag, s1.real*s2.real-s1.imag*s2.imag);
293 }
294 errno = 0;
295 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000296}
297
Guido van Rossumc6e22901998-12-04 19:26:43 +0000298
Brett Cannonb0fc4902014-10-14 17:37:02 -0400299/*[clinic input]
300cmath.atan = cmath.acos
301
302Return the arc tangent of z.
303[clinic start generated code]*/
Guido van Rossumc6e22901998-12-04 19:26:43 +0000304
Tim Peters14e26402001-02-20 20:15:19 +0000305static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300306cmath_atan_impl(PyObject *module, Py_complex z)
307/*[clinic end generated code: output=b6bfc497058acba4 input=3b21ff7d5eac632a]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000308{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000309 /* atan(z) = -i atanh(iz) */
310 Py_complex s, r;
311 s.real = -z.imag;
312 s.imag = z.real;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400313 s = cmath_atanh_impl(module, s);
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000314 r.real = s.imag;
315 r.imag = -s.real;
316 return r;
Christian Heimes53876d92008-04-19 00:31:39 +0000317}
318
Christian Heimese57950f2008-04-21 13:08:03 +0000319/* Windows screws up atan2 for inf and nan, and alpha Tru64 5.1 doesn't follow
320 C99 for atan2(0., 0.). */
Christian Heimes53876d92008-04-19 00:31:39 +0000321static double
322c_atan2(Py_complex z)
323{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000324 if (Py_IS_NAN(z.real) || Py_IS_NAN(z.imag))
325 return Py_NAN;
326 if (Py_IS_INFINITY(z.imag)) {
327 if (Py_IS_INFINITY(z.real)) {
328 if (copysign(1., z.real) == 1.)
329 /* atan2(+-inf, +inf) == +-pi/4 */
330 return copysign(0.25*Py_MATH_PI, z.imag);
331 else
332 /* atan2(+-inf, -inf) == +-pi*3/4 */
333 return copysign(0.75*Py_MATH_PI, z.imag);
334 }
335 /* atan2(+-inf, x) == +-pi/2 for finite x */
336 return copysign(0.5*Py_MATH_PI, z.imag);
337 }
338 if (Py_IS_INFINITY(z.real) || z.imag == 0.) {
339 if (copysign(1., z.real) == 1.)
340 /* atan2(+-y, +inf) = atan2(+-0, +x) = +-0. */
341 return copysign(0., z.imag);
342 else
343 /* atan2(+-y, -inf) = atan2(+-0., -x) = +-pi. */
344 return copysign(Py_MATH_PI, z.imag);
345 }
346 return atan2(z.imag, z.real);
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000347}
348
Guido van Rossumc6e22901998-12-04 19:26:43 +0000349
Christian Heimesa342c012008-04-20 21:01:16 +0000350static Py_complex atanh_special_values[7][7];
Christian Heimes53876d92008-04-19 00:31:39 +0000351
Brett Cannonb0fc4902014-10-14 17:37:02 -0400352/*[clinic input]
353cmath.atanh = cmath.acos
354
Mark Dickinsoncc8617b2015-01-11 13:22:44 +0000355Return the inverse hyperbolic tangent of z.
Brett Cannonb0fc4902014-10-14 17:37:02 -0400356[clinic start generated code]*/
357
Tim Peters14e26402001-02-20 20:15:19 +0000358static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300359cmath_atanh_impl(PyObject *module, Py_complex z)
360/*[clinic end generated code: output=e83355f93a989c9e input=2b3fdb82fb34487b]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000361{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000362 Py_complex r;
363 double ay, h;
Christian Heimes53876d92008-04-19 00:31:39 +0000364
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000365 SPECIAL_VALUE(z, atanh_special_values);
Christian Heimes53876d92008-04-19 00:31:39 +0000366
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000367 /* Reduce to case where z.real >= 0., using atanh(z) = -atanh(-z). */
368 if (z.real < 0.) {
Brett Cannonb0fc4902014-10-14 17:37:02 -0400369 return _Py_c_neg(cmath_atanh_impl(module, _Py_c_neg(z)));
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000370 }
Christian Heimes53876d92008-04-19 00:31:39 +0000371
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000372 ay = fabs(z.imag);
373 if (z.real > CM_SQRT_LARGE_DOUBLE || ay > CM_SQRT_LARGE_DOUBLE) {
374 /*
375 if abs(z) is large then we use the approximation
376 atanh(z) ~ 1/z +/- i*pi/2 (+/- depending on the sign
377 of z.imag)
378 */
379 h = hypot(z.real/2., z.imag/2.); /* safe from overflow */
380 r.real = z.real/4./h/h;
381 /* the two negations in the next line cancel each other out
382 except when working with unsigned zeros: they're there to
383 ensure that the branch cut has the correct continuity on
384 systems that don't support signed zeros */
385 r.imag = -copysign(Py_MATH_PI/2., -z.imag);
386 errno = 0;
387 } else if (z.real == 1. && ay < CM_SQRT_DBL_MIN) {
388 /* C99 standard says: atanh(1+/-0.) should be inf +/- 0i */
389 if (ay == 0.) {
390 r.real = INF;
391 r.imag = z.imag;
392 errno = EDOM;
393 } else {
394 r.real = -log(sqrt(ay)/sqrt(hypot(ay, 2.)));
395 r.imag = copysign(atan2(2., -ay)/2, z.imag);
396 errno = 0;
397 }
398 } else {
399 r.real = m_log1p(4.*z.real/((1-z.real)*(1-z.real) + ay*ay))/4.;
400 r.imag = -atan2(-2.*z.imag, (1-z.real)*(1+z.real) - ay*ay)/2.;
401 errno = 0;
402 }
403 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000404}
405
Guido van Rossumc6e22901998-12-04 19:26:43 +0000406
Brett Cannonb0fc4902014-10-14 17:37:02 -0400407/*[clinic input]
408cmath.cos = cmath.acos
409
410Return the cosine of z.
411[clinic start generated code]*/
Guido van Rossumc6e22901998-12-04 19:26:43 +0000412
Tim Peters14e26402001-02-20 20:15:19 +0000413static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300414cmath_cos_impl(PyObject *module, Py_complex z)
415/*[clinic end generated code: output=fd64918d5b3186db input=6022e39b77127ac7]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000416{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000417 /* cos(z) = cosh(iz) */
418 Py_complex r;
419 r.real = -z.imag;
420 r.imag = z.real;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400421 r = cmath_cosh_impl(module, r);
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000422 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000423}
424
Guido van Rossumc6e22901998-12-04 19:26:43 +0000425
Christian Heimes53876d92008-04-19 00:31:39 +0000426/* cosh(infinity + i*y) needs to be dealt with specially */
Christian Heimesa342c012008-04-20 21:01:16 +0000427static Py_complex cosh_special_values[7][7];
Christian Heimes53876d92008-04-19 00:31:39 +0000428
Brett Cannonb0fc4902014-10-14 17:37:02 -0400429/*[clinic input]
430cmath.cosh = cmath.acos
431
432Return the hyperbolic cosine of z.
433[clinic start generated code]*/
434
Tim Peters14e26402001-02-20 20:15:19 +0000435static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300436cmath_cosh_impl(PyObject *module, Py_complex z)
437/*[clinic end generated code: output=2e969047da601bdb input=d6b66339e9cc332b]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000438{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000439 Py_complex r;
440 double x_minus_one;
Christian Heimes53876d92008-04-19 00:31:39 +0000441
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000442 /* special treatment for cosh(+/-inf + iy) if y is not a NaN */
443 if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) {
444 if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) &&
445 (z.imag != 0.)) {
446 if (z.real > 0) {
447 r.real = copysign(INF, cos(z.imag));
448 r.imag = copysign(INF, sin(z.imag));
449 }
450 else {
451 r.real = copysign(INF, cos(z.imag));
452 r.imag = -copysign(INF, sin(z.imag));
453 }
454 }
455 else {
456 r = cosh_special_values[special_type(z.real)]
457 [special_type(z.imag)];
458 }
459 /* need to set errno = EDOM if y is +/- infinity and x is not
460 a NaN */
461 if (Py_IS_INFINITY(z.imag) && !Py_IS_NAN(z.real))
462 errno = EDOM;
463 else
464 errno = 0;
465 return r;
466 }
Christian Heimes53876d92008-04-19 00:31:39 +0000467
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000468 if (fabs(z.real) > CM_LOG_LARGE_DOUBLE) {
469 /* deal correctly with cases where cosh(z.real) overflows but
470 cosh(z) does not. */
471 x_minus_one = z.real - copysign(1., z.real);
472 r.real = cos(z.imag) * cosh(x_minus_one) * Py_MATH_E;
473 r.imag = sin(z.imag) * sinh(x_minus_one) * Py_MATH_E;
474 } else {
475 r.real = cos(z.imag) * cosh(z.real);
476 r.imag = sin(z.imag) * sinh(z.real);
477 }
478 /* detect overflow, and set errno accordingly */
479 if (Py_IS_INFINITY(r.real) || Py_IS_INFINITY(r.imag))
480 errno = ERANGE;
481 else
482 errno = 0;
483 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000484}
485
Guido van Rossumc6e22901998-12-04 19:26:43 +0000486
Christian Heimes53876d92008-04-19 00:31:39 +0000487/* exp(infinity + i*y) and exp(-infinity + i*y) need special treatment for
488 finite y */
Christian Heimesa342c012008-04-20 21:01:16 +0000489static Py_complex exp_special_values[7][7];
Christian Heimes53876d92008-04-19 00:31:39 +0000490
Brett Cannonb0fc4902014-10-14 17:37:02 -0400491/*[clinic input]
492cmath.exp = cmath.acos
493
494Return the exponential value e**z.
495[clinic start generated code]*/
496
Tim Peters14e26402001-02-20 20:15:19 +0000497static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300498cmath_exp_impl(PyObject *module, Py_complex z)
499/*[clinic end generated code: output=edcec61fb9dfda6c input=8b9e6cf8a92174c3]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000500{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000501 Py_complex r;
502 double l;
Christian Heimes53876d92008-04-19 00:31:39 +0000503
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000504 if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) {
505 if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag)
506 && (z.imag != 0.)) {
507 if (z.real > 0) {
508 r.real = copysign(INF, cos(z.imag));
509 r.imag = copysign(INF, sin(z.imag));
510 }
511 else {
512 r.real = copysign(0., cos(z.imag));
513 r.imag = copysign(0., sin(z.imag));
514 }
515 }
516 else {
517 r = exp_special_values[special_type(z.real)]
518 [special_type(z.imag)];
519 }
520 /* need to set errno = EDOM if y is +/- infinity and x is not
521 a NaN and not -infinity */
522 if (Py_IS_INFINITY(z.imag) &&
523 (Py_IS_FINITE(z.real) ||
524 (Py_IS_INFINITY(z.real) && z.real > 0)))
525 errno = EDOM;
526 else
527 errno = 0;
528 return r;
529 }
Christian Heimes53876d92008-04-19 00:31:39 +0000530
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000531 if (z.real > CM_LOG_LARGE_DOUBLE) {
532 l = exp(z.real-1.);
533 r.real = l*cos(z.imag)*Py_MATH_E;
534 r.imag = l*sin(z.imag)*Py_MATH_E;
535 } else {
536 l = exp(z.real);
537 r.real = l*cos(z.imag);
538 r.imag = l*sin(z.imag);
539 }
540 /* detect overflow, and set errno accordingly */
541 if (Py_IS_INFINITY(r.real) || Py_IS_INFINITY(r.imag))
542 errno = ERANGE;
543 else
544 errno = 0;
545 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000546}
547
Christian Heimesa342c012008-04-20 21:01:16 +0000548static Py_complex log_special_values[7][7];
Christian Heimes53876d92008-04-19 00:31:39 +0000549
Tim Peters14e26402001-02-20 20:15:19 +0000550static Py_complex
Christian Heimes53876d92008-04-19 00:31:39 +0000551c_log(Py_complex z)
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000552{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000553 /*
554 The usual formula for the real part is log(hypot(z.real, z.imag)).
555 There are four situations where this formula is potentially
556 problematic:
Christian Heimes53876d92008-04-19 00:31:39 +0000557
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000558 (1) the absolute value of z is subnormal. Then hypot is subnormal,
559 so has fewer than the usual number of bits of accuracy, hence may
560 have large relative error. This then gives a large absolute error
561 in the log. This can be solved by rescaling z by a suitable power
562 of 2.
Christian Heimes53876d92008-04-19 00:31:39 +0000563
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000564 (2) the absolute value of z is greater than DBL_MAX (e.g. when both
565 z.real and z.imag are within a factor of 1/sqrt(2) of DBL_MAX)
566 Again, rescaling solves this.
Christian Heimes53876d92008-04-19 00:31:39 +0000567
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000568 (3) the absolute value of z is close to 1. In this case it's
569 difficult to achieve good accuracy, at least in part because a
570 change of 1ulp in the real or imaginary part of z can result in a
571 change of billions of ulps in the correctly rounded answer.
Christian Heimes53876d92008-04-19 00:31:39 +0000572
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000573 (4) z = 0. The simplest thing to do here is to call the
574 floating-point log with an argument of 0, and let its behaviour
575 (returning -infinity, signaling a floating-point exception, setting
576 errno, or whatever) determine that of c_log. So the usual formula
577 is fine here.
Christian Heimes53876d92008-04-19 00:31:39 +0000578
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000579 */
Christian Heimes53876d92008-04-19 00:31:39 +0000580
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000581 Py_complex r;
582 double ax, ay, am, an, h;
Christian Heimes53876d92008-04-19 00:31:39 +0000583
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000584 SPECIAL_VALUE(z, log_special_values);
Christian Heimes53876d92008-04-19 00:31:39 +0000585
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000586 ax = fabs(z.real);
587 ay = fabs(z.imag);
Christian Heimes53876d92008-04-19 00:31:39 +0000588
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000589 if (ax > CM_LARGE_DOUBLE || ay > CM_LARGE_DOUBLE) {
590 r.real = log(hypot(ax/2., ay/2.)) + M_LN2;
591 } else if (ax < DBL_MIN && ay < DBL_MIN) {
592 if (ax > 0. || ay > 0.) {
593 /* catch cases where hypot(ax, ay) is subnormal */
594 r.real = log(hypot(ldexp(ax, DBL_MANT_DIG),
595 ldexp(ay, DBL_MANT_DIG))) - DBL_MANT_DIG*M_LN2;
596 }
597 else {
598 /* log(+/-0. +/- 0i) */
599 r.real = -INF;
600 r.imag = atan2(z.imag, z.real);
601 errno = EDOM;
602 return r;
603 }
604 } else {
605 h = hypot(ax, ay);
606 if (0.71 <= h && h <= 1.73) {
607 am = ax > ay ? ax : ay; /* max(ax, ay) */
608 an = ax > ay ? ay : ax; /* min(ax, ay) */
609 r.real = m_log1p((am-1)*(am+1)+an*an)/2.;
610 } else {
611 r.real = log(h);
612 }
613 }
614 r.imag = atan2(z.imag, z.real);
615 errno = 0;
616 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000617}
618
Guido van Rossumc6e22901998-12-04 19:26:43 +0000619
Brett Cannonb0fc4902014-10-14 17:37:02 -0400620/*[clinic input]
621cmath.log10 = cmath.acos
622
623Return the base-10 logarithm of z.
624[clinic start generated code]*/
625
Tim Peters14e26402001-02-20 20:15:19 +0000626static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300627cmath_log10_impl(PyObject *module, Py_complex z)
628/*[clinic end generated code: output=2922779a7c38cbe1 input=cff5644f73c1519c]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000629{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000630 Py_complex r;
631 int errno_save;
Christian Heimes53876d92008-04-19 00:31:39 +0000632
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000633 r = c_log(z);
634 errno_save = errno; /* just in case the divisions affect errno */
635 r.real = r.real / M_LN10;
636 r.imag = r.imag / M_LN10;
637 errno = errno_save;
638 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000639}
640
Guido van Rossumc6e22901998-12-04 19:26:43 +0000641
Brett Cannonb0fc4902014-10-14 17:37:02 -0400642/*[clinic input]
643cmath.sin = cmath.acos
644
645Return the sine of z.
646[clinic start generated code]*/
Guido van Rossumc6e22901998-12-04 19:26:43 +0000647
Tim Peters14e26402001-02-20 20:15:19 +0000648static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300649cmath_sin_impl(PyObject *module, Py_complex z)
650/*[clinic end generated code: output=980370d2ff0bb5aa input=2d3519842a8b4b85]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000651{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000652 /* sin(z) = -i sin(iz) */
653 Py_complex s, r;
654 s.real = -z.imag;
655 s.imag = z.real;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400656 s = cmath_sinh_impl(module, s);
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000657 r.real = s.imag;
658 r.imag = -s.real;
659 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000660}
661
Guido van Rossumc6e22901998-12-04 19:26:43 +0000662
Christian Heimes53876d92008-04-19 00:31:39 +0000663/* sinh(infinity + i*y) needs to be dealt with specially */
Christian Heimesa342c012008-04-20 21:01:16 +0000664static Py_complex sinh_special_values[7][7];
Christian Heimes53876d92008-04-19 00:31:39 +0000665
Brett Cannonb0fc4902014-10-14 17:37:02 -0400666/*[clinic input]
667cmath.sinh = cmath.acos
668
669Return the hyperbolic sine of z.
670[clinic start generated code]*/
671
Tim Peters14e26402001-02-20 20:15:19 +0000672static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300673cmath_sinh_impl(PyObject *module, Py_complex z)
674/*[clinic end generated code: output=38b0a6cce26f3536 input=d2d3fc8c1ddfd2dd]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000675{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000676 Py_complex r;
677 double x_minus_one;
Christian Heimes53876d92008-04-19 00:31:39 +0000678
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000679 /* special treatment for sinh(+/-inf + iy) if y is finite and
680 nonzero */
681 if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) {
682 if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag)
683 && (z.imag != 0.)) {
684 if (z.real > 0) {
685 r.real = copysign(INF, cos(z.imag));
686 r.imag = copysign(INF, sin(z.imag));
687 }
688 else {
689 r.real = -copysign(INF, cos(z.imag));
690 r.imag = copysign(INF, sin(z.imag));
691 }
692 }
693 else {
694 r = sinh_special_values[special_type(z.real)]
695 [special_type(z.imag)];
696 }
697 /* need to set errno = EDOM if y is +/- infinity and x is not
698 a NaN */
699 if (Py_IS_INFINITY(z.imag) && !Py_IS_NAN(z.real))
700 errno = EDOM;
701 else
702 errno = 0;
703 return r;
704 }
Christian Heimes53876d92008-04-19 00:31:39 +0000705
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000706 if (fabs(z.real) > CM_LOG_LARGE_DOUBLE) {
707 x_minus_one = z.real - copysign(1., z.real);
708 r.real = cos(z.imag) * sinh(x_minus_one) * Py_MATH_E;
709 r.imag = sin(z.imag) * cosh(x_minus_one) * Py_MATH_E;
710 } else {
711 r.real = cos(z.imag) * sinh(z.real);
712 r.imag = sin(z.imag) * cosh(z.real);
713 }
714 /* detect overflow, and set errno accordingly */
715 if (Py_IS_INFINITY(r.real) || Py_IS_INFINITY(r.imag))
716 errno = ERANGE;
717 else
718 errno = 0;
719 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000720}
721
Guido van Rossumc6e22901998-12-04 19:26:43 +0000722
Christian Heimesa342c012008-04-20 21:01:16 +0000723static Py_complex sqrt_special_values[7][7];
Christian Heimes53876d92008-04-19 00:31:39 +0000724
Brett Cannonb0fc4902014-10-14 17:37:02 -0400725/*[clinic input]
726cmath.sqrt = cmath.acos
727
728Return the square root of z.
729[clinic start generated code]*/
730
Tim Peters14e26402001-02-20 20:15:19 +0000731static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300732cmath_sqrt_impl(PyObject *module, Py_complex z)
733/*[clinic end generated code: output=b6507b3029c339fc input=7088b166fc9a58c7]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000734{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000735 /*
736 Method: use symmetries to reduce to the case when x = z.real and y
737 = z.imag are nonnegative. Then the real part of the result is
738 given by
Christian Heimes53876d92008-04-19 00:31:39 +0000739
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000740 s = sqrt((x + hypot(x, y))/2)
Christian Heimes53876d92008-04-19 00:31:39 +0000741
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000742 and the imaginary part is
Christian Heimes53876d92008-04-19 00:31:39 +0000743
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000744 d = (y/2)/s
Christian Heimes53876d92008-04-19 00:31:39 +0000745
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000746 If either x or y is very large then there's a risk of overflow in
747 computation of the expression x + hypot(x, y). We can avoid this
748 by rewriting the formula for s as:
Christian Heimes53876d92008-04-19 00:31:39 +0000749
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000750 s = 2*sqrt(x/8 + hypot(x/8, y/8))
Christian Heimes53876d92008-04-19 00:31:39 +0000751
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000752 This costs us two extra multiplications/divisions, but avoids the
753 overhead of checking for x and y large.
Christian Heimes53876d92008-04-19 00:31:39 +0000754
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000755 If both x and y are subnormal then hypot(x, y) may also be
756 subnormal, so will lack full precision. We solve this by rescaling
757 x and y by a sufficiently large power of 2 to ensure that x and y
758 are normal.
759 */
Christian Heimes53876d92008-04-19 00:31:39 +0000760
761
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000762 Py_complex r;
763 double s,d;
764 double ax, ay;
Christian Heimes53876d92008-04-19 00:31:39 +0000765
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000766 SPECIAL_VALUE(z, sqrt_special_values);
Christian Heimes53876d92008-04-19 00:31:39 +0000767
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000768 if (z.real == 0. && z.imag == 0.) {
769 r.real = 0.;
770 r.imag = z.imag;
771 return r;
772 }
Christian Heimes53876d92008-04-19 00:31:39 +0000773
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000774 ax = fabs(z.real);
775 ay = fabs(z.imag);
Christian Heimes53876d92008-04-19 00:31:39 +0000776
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000777 if (ax < DBL_MIN && ay < DBL_MIN && (ax > 0. || ay > 0.)) {
778 /* here we catch cases where hypot(ax, ay) is subnormal */
779 ax = ldexp(ax, CM_SCALE_UP);
780 s = ldexp(sqrt(ax + hypot(ax, ldexp(ay, CM_SCALE_UP))),
781 CM_SCALE_DOWN);
782 } else {
783 ax /= 8.;
784 s = 2.*sqrt(ax + hypot(ax, ay/8.));
785 }
786 d = ay/(2.*s);
Christian Heimes53876d92008-04-19 00:31:39 +0000787
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000788 if (z.real >= 0.) {
789 r.real = s;
790 r.imag = copysign(d, z.imag);
791 } else {
792 r.real = d;
793 r.imag = copysign(s, z.imag);
794 }
795 errno = 0;
796 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000797}
798
Guido van Rossumc6e22901998-12-04 19:26:43 +0000799
Brett Cannonb0fc4902014-10-14 17:37:02 -0400800/*[clinic input]
801cmath.tan = cmath.acos
802
803Return the tangent of z.
804[clinic start generated code]*/
Guido van Rossumc6e22901998-12-04 19:26:43 +0000805
Tim Peters14e26402001-02-20 20:15:19 +0000806static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300807cmath_tan_impl(PyObject *module, Py_complex z)
808/*[clinic end generated code: output=7c5f13158a72eb13 input=fc167e528767888e]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000809{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000810 /* tan(z) = -i tanh(iz) */
811 Py_complex s, r;
812 s.real = -z.imag;
813 s.imag = z.real;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400814 s = cmath_tanh_impl(module, s);
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000815 r.real = s.imag;
816 r.imag = -s.real;
817 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000818}
819
Guido van Rossumc6e22901998-12-04 19:26:43 +0000820
Christian Heimes53876d92008-04-19 00:31:39 +0000821/* tanh(infinity + i*y) needs to be dealt with specially */
Christian Heimesa342c012008-04-20 21:01:16 +0000822static Py_complex tanh_special_values[7][7];
Christian Heimes53876d92008-04-19 00:31:39 +0000823
Brett Cannonb0fc4902014-10-14 17:37:02 -0400824/*[clinic input]
825cmath.tanh = cmath.acos
826
827Return the hyperbolic tangent of z.
828[clinic start generated code]*/
829
Tim Peters14e26402001-02-20 20:15:19 +0000830static Py_complex
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300831cmath_tanh_impl(PyObject *module, Py_complex z)
832/*[clinic end generated code: output=36d547ef7aca116c input=22f67f9dc6d29685]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000833{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000834 /* Formula:
Christian Heimes53876d92008-04-19 00:31:39 +0000835
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000836 tanh(x+iy) = (tanh(x)(1+tan(y)^2) + i tan(y)(1-tanh(x))^2) /
837 (1+tan(y)^2 tanh(x)^2)
Christian Heimes53876d92008-04-19 00:31:39 +0000838
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000839 To avoid excessive roundoff error, 1-tanh(x)^2 is better computed
840 as 1/cosh(x)^2. When abs(x) is large, we approximate 1-tanh(x)^2
841 by 4 exp(-2*x) instead, to avoid possible overflow in the
842 computation of cosh(x).
Christian Heimes53876d92008-04-19 00:31:39 +0000843
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000844 */
Christian Heimes53876d92008-04-19 00:31:39 +0000845
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000846 Py_complex r;
847 double tx, ty, cx, txty, denom;
Christian Heimes53876d92008-04-19 00:31:39 +0000848
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000849 /* special treatment for tanh(+/-inf + iy) if y is finite and
850 nonzero */
851 if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) {
852 if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag)
853 && (z.imag != 0.)) {
854 if (z.real > 0) {
855 r.real = 1.0;
856 r.imag = copysign(0.,
857 2.*sin(z.imag)*cos(z.imag));
858 }
859 else {
860 r.real = -1.0;
861 r.imag = copysign(0.,
862 2.*sin(z.imag)*cos(z.imag));
863 }
864 }
865 else {
866 r = tanh_special_values[special_type(z.real)]
867 [special_type(z.imag)];
868 }
869 /* need to set errno = EDOM if z.imag is +/-infinity and
870 z.real is finite */
871 if (Py_IS_INFINITY(z.imag) && Py_IS_FINITE(z.real))
872 errno = EDOM;
873 else
874 errno = 0;
875 return r;
876 }
Christian Heimes53876d92008-04-19 00:31:39 +0000877
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000878 /* danger of overflow in 2.*z.imag !*/
879 if (fabs(z.real) > CM_LOG_LARGE_DOUBLE) {
880 r.real = copysign(1., z.real);
881 r.imag = 4.*sin(z.imag)*cos(z.imag)*exp(-2.*fabs(z.real));
882 } else {
883 tx = tanh(z.real);
884 ty = tan(z.imag);
885 cx = 1./cosh(z.real);
886 txty = tx*ty;
887 denom = 1. + txty*txty;
888 r.real = tx*(1.+ty*ty)/denom;
889 r.imag = ((ty/denom)*cx)*cx;
890 }
891 errno = 0;
892 return r;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000893}
894
Guido van Rossumc6e22901998-12-04 19:26:43 +0000895
Brett Cannonb0fc4902014-10-14 17:37:02 -0400896/*[clinic input]
897cmath.log
898
899 x: Py_complex
900 y_obj: object = NULL
901 /
902
903The logarithm of z to the given base.
904
905If the base not specified, returns the natural logarithm (base e) of z.
906[clinic start generated code]*/
Christian Heimes53876d92008-04-19 00:31:39 +0000907
Raymond Hettingerb67ad7e2004-06-14 07:40:10 +0000908static PyObject *
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300909cmath_log_impl(PyObject *module, Py_complex x, PyObject *y_obj)
910/*[clinic end generated code: output=4effdb7d258e0d94 input=ee0e823a7c6e68ea]*/
Raymond Hettingerb67ad7e2004-06-14 07:40:10 +0000911{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000912 Py_complex y;
Raymond Hettingerb67ad7e2004-06-14 07:40:10 +0000913
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000914 errno = 0;
915 PyFPE_START_PROTECT("complex function", return 0)
916 x = c_log(x);
Brett Cannonb0fc4902014-10-14 17:37:02 -0400917 if (y_obj != NULL) {
918 y = PyComplex_AsCComplex(y_obj);
919 if (PyErr_Occurred()) {
920 return NULL;
921 }
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000922 y = c_log(y);
Antoine Pitroude08cb62014-07-07 19:08:47 -0400923 x = _Py_c_quot(x, y);
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000924 }
925 PyFPE_END_PROTECT(x)
926 if (errno != 0)
927 return math_error();
928 return PyComplex_FromCComplex(x);
Raymond Hettingerb67ad7e2004-06-14 07:40:10 +0000929}
930
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000931
932/* And now the glue to make them available from Python: */
933
Roger E. Masse24070ca1996-12-09 22:59:53 +0000934static PyObject *
Thomas Woutersf3f33dc2000-07-21 06:00:07 +0000935math_error(void)
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000936{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000937 if (errno == EDOM)
938 PyErr_SetString(PyExc_ValueError, "math domain error");
939 else if (errno == ERANGE)
940 PyErr_SetString(PyExc_OverflowError, "math range error");
941 else /* Unexpected math error */
942 PyErr_SetFromErrno(PyExc_ValueError);
943 return NULL;
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000944}
945
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000946
Brett Cannonb0fc4902014-10-14 17:37:02 -0400947/*[clinic input]
948cmath.phase
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000949
Brett Cannonb0fc4902014-10-14 17:37:02 -0400950 z: Py_complex
951 /
952
953Return argument, also known as the phase angle, of a complex.
954[clinic start generated code]*/
Guido van Rossum71aa32f1996-01-12 01:34:57 +0000955
Christian Heimes53876d92008-04-19 00:31:39 +0000956static PyObject *
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300957cmath_phase_impl(PyObject *module, Py_complex z)
958/*[clinic end generated code: output=50725086a7bfd253 input=5cf75228ba94b69d]*/
Christian Heimes53876d92008-04-19 00:31:39 +0000959{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000960 double phi;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400961
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000962 errno = 0;
963 PyFPE_START_PROTECT("arg function", return 0)
964 phi = c_atan2(z);
965 PyFPE_END_PROTECT(phi)
966 if (errno != 0)
967 return math_error();
968 else
969 return PyFloat_FromDouble(phi);
Christian Heimes53876d92008-04-19 00:31:39 +0000970}
971
Brett Cannonb0fc4902014-10-14 17:37:02 -0400972/*[clinic input]
973cmath.polar
974
975 z: Py_complex
976 /
977
978Convert a complex from rectangular coordinates to polar coordinates.
979
980r is the distance from 0 and phi the phase angle.
981[clinic start generated code]*/
Christian Heimes53876d92008-04-19 00:31:39 +0000982
983static PyObject *
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +0300984cmath_polar_impl(PyObject *module, Py_complex z)
985/*[clinic end generated code: output=d0a8147c41dbb654 input=26c353574fd1a861]*/
Christian Heimes53876d92008-04-19 00:31:39 +0000986{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000987 double r, phi;
Brett Cannonb0fc4902014-10-14 17:37:02 -0400988
Antoine Pitrou6bc217d2015-06-23 14:31:11 +0200989 errno = 0;
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000990 PyFPE_START_PROTECT("polar function", return 0)
991 phi = c_atan2(z); /* should not cause any exception */
Antoine Pitroua72f0cd2015-06-23 14:38:13 +0200992 r = _Py_c_abs(z); /* sets errno to ERANGE on overflow */
Antoine Pitrouf95a1b32010-05-09 15:52:27 +0000993 PyFPE_END_PROTECT(r)
994 if (errno != 0)
995 return math_error();
996 else
997 return Py_BuildValue("dd", r, phi);
Christian Heimes53876d92008-04-19 00:31:39 +0000998}
999
Christian Heimes53876d92008-04-19 00:31:39 +00001000/*
1001 rect() isn't covered by the C99 standard, but it's not too hard to
1002 figure out 'spirit of C99' rules for special value handing:
1003
1004 rect(x, t) should behave like exp(log(x) + it) for positive-signed x
1005 rect(x, t) should behave like -exp(log(-x) + it) for negative-signed x
1006 rect(nan, t) should behave like exp(nan + it), except that rect(nan, 0)
1007 gives nan +- i0 with the sign of the imaginary part unspecified.
1008
1009*/
1010
Christian Heimesa342c012008-04-20 21:01:16 +00001011static Py_complex rect_special_values[7][7];
Christian Heimes53876d92008-04-19 00:31:39 +00001012
Brett Cannonb0fc4902014-10-14 17:37:02 -04001013/*[clinic input]
1014cmath.rect
1015
1016 r: double
1017 phi: double
1018 /
1019
1020Convert from polar coordinates to rectangular coordinates.
1021[clinic start generated code]*/
1022
Christian Heimes53876d92008-04-19 00:31:39 +00001023static PyObject *
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +03001024cmath_rect_impl(PyObject *module, double r, double phi)
1025/*[clinic end generated code: output=385a0690925df2d5 input=24c5646d147efd69]*/
Christian Heimes53876d92008-04-19 00:31:39 +00001026{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001027 Py_complex z;
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001028 errno = 0;
1029 PyFPE_START_PROTECT("rect function", return 0)
Christian Heimes53876d92008-04-19 00:31:39 +00001030
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001031 /* deal with special values */
1032 if (!Py_IS_FINITE(r) || !Py_IS_FINITE(phi)) {
1033 /* if r is +/-infinity and phi is finite but nonzero then
1034 result is (+-INF +-INF i), but we need to compute cos(phi)
1035 and sin(phi) to figure out the signs. */
1036 if (Py_IS_INFINITY(r) && (Py_IS_FINITE(phi)
1037 && (phi != 0.))) {
1038 if (r > 0) {
1039 z.real = copysign(INF, cos(phi));
1040 z.imag = copysign(INF, sin(phi));
1041 }
1042 else {
1043 z.real = -copysign(INF, cos(phi));
1044 z.imag = -copysign(INF, sin(phi));
1045 }
1046 }
1047 else {
1048 z = rect_special_values[special_type(r)]
1049 [special_type(phi)];
1050 }
1051 /* need to set errno = EDOM if r is a nonzero number and phi
1052 is infinite */
1053 if (r != 0. && !Py_IS_NAN(r) && Py_IS_INFINITY(phi))
1054 errno = EDOM;
1055 else
1056 errno = 0;
1057 }
Mark Dickinson58ceecf2013-07-20 17:59:13 +01001058 else if (phi == 0.0) {
1059 /* Workaround for buggy results with phi=-0.0 on OS X 10.8. See
1060 bugs.python.org/issue18513. */
1061 z.real = r;
1062 z.imag = r * phi;
1063 errno = 0;
1064 }
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001065 else {
1066 z.real = r * cos(phi);
1067 z.imag = r * sin(phi);
1068 errno = 0;
1069 }
Christian Heimes53876d92008-04-19 00:31:39 +00001070
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001071 PyFPE_END_PROTECT(z)
1072 if (errno != 0)
1073 return math_error();
1074 else
1075 return PyComplex_FromCComplex(z);
Christian Heimes53876d92008-04-19 00:31:39 +00001076}
1077
Brett Cannonb0fc4902014-10-14 17:37:02 -04001078/*[clinic input]
1079cmath.isfinite = cmath.polar
1080
1081Return True if both the real and imaginary parts of z are finite, else False.
1082[clinic start generated code]*/
Christian Heimes53876d92008-04-19 00:31:39 +00001083
1084static PyObject *
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +03001085cmath_isfinite_impl(PyObject *module, Py_complex z)
1086/*[clinic end generated code: output=ac76611e2c774a36 input=848e7ee701895815]*/
Mark Dickinson8e0c9962010-07-11 17:38:24 +00001087{
Mark Dickinson8e0c9962010-07-11 17:38:24 +00001088 return PyBool_FromLong(Py_IS_FINITE(z.real) && Py_IS_FINITE(z.imag));
1089}
1090
Brett Cannonb0fc4902014-10-14 17:37:02 -04001091/*[clinic input]
1092cmath.isnan = cmath.polar
1093
1094Checks if the real or imaginary part of z not a number (NaN).
1095[clinic start generated code]*/
Mark Dickinson8e0c9962010-07-11 17:38:24 +00001096
1097static PyObject *
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +03001098cmath_isnan_impl(PyObject *module, Py_complex z)
1099/*[clinic end generated code: output=e7abf6e0b28beab7 input=71799f5d284c9baf]*/
Christian Heimes53876d92008-04-19 00:31:39 +00001100{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001101 return PyBool_FromLong(Py_IS_NAN(z.real) || Py_IS_NAN(z.imag));
Christian Heimes53876d92008-04-19 00:31:39 +00001102}
1103
Brett Cannonb0fc4902014-10-14 17:37:02 -04001104/*[clinic input]
1105cmath.isinf = cmath.polar
1106
1107Checks if the real or imaginary part of z is infinite.
1108[clinic start generated code]*/
Christian Heimes53876d92008-04-19 00:31:39 +00001109
1110static PyObject *
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +03001111cmath_isinf_impl(PyObject *module, Py_complex z)
1112/*[clinic end generated code: output=502a75a79c773469 input=363df155c7181329]*/
Christian Heimes53876d92008-04-19 00:31:39 +00001113{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001114 return PyBool_FromLong(Py_IS_INFINITY(z.real) ||
1115 Py_IS_INFINITY(z.imag));
Christian Heimes53876d92008-04-19 00:31:39 +00001116}
1117
Tal Einatd5519ed2015-05-31 22:05:00 +03001118/*[clinic input]
1119cmath.isclose -> bool
1120
1121 a: Py_complex
1122 b: Py_complex
1123 *
1124 rel_tol: double = 1e-09
1125 maximum difference for being considered "close", relative to the
1126 magnitude of the input values
1127 abs_tol: double = 0.0
1128 maximum difference for being considered "close", regardless of the
1129 magnitude of the input values
1130
1131Determine whether two complex numbers are close in value.
1132
1133Return True if a is close in value to b, and False otherwise.
1134
1135For the values to be considered close, the difference between them must be
1136smaller than at least one of the tolerances.
1137
1138-inf, inf and NaN behave similarly to the IEEE 754 Standard. That is, NaN is
1139not close to anything, even itself. inf and -inf are only close to themselves.
1140[clinic start generated code]*/
1141
1142static int
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +03001143cmath_isclose_impl(PyObject *module, Py_complex a, Py_complex b,
Tal Einatd5519ed2015-05-31 22:05:00 +03001144 double rel_tol, double abs_tol)
Serhiy Storchaka1a2b24f2016-07-07 17:35:15 +03001145/*[clinic end generated code: output=8a2486cc6e0014d1 input=df9636d7de1d4ac3]*/
Tal Einatd5519ed2015-05-31 22:05:00 +03001146{
1147 double diff;
1148
1149 /* sanity check on the inputs */
1150 if (rel_tol < 0.0 || abs_tol < 0.0 ) {
1151 PyErr_SetString(PyExc_ValueError,
1152 "tolerances must be non-negative");
1153 return -1;
1154 }
1155
1156 if ( (a.real == b.real) && (a.imag == b.imag) ) {
1157 /* short circuit exact equality -- needed to catch two infinities of
1158 the same sign. And perhaps speeds things up a bit sometimes.
1159 */
1160 return 1;
1161 }
1162
1163 /* This catches the case of two infinities of opposite sign, or
1164 one infinity and one finite number. Two infinities of opposite
1165 sign would otherwise have an infinite relative tolerance.
1166 Two infinities of the same sign are caught by the equality check
1167 above.
1168 */
1169
1170 if (Py_IS_INFINITY(a.real) || Py_IS_INFINITY(a.imag) ||
1171 Py_IS_INFINITY(b.real) || Py_IS_INFINITY(b.imag)) {
1172 return 0;
1173 }
1174
1175 /* now do the regular computation
1176 this is essentially the "weak" test from the Boost library
1177 */
1178
1179 diff = _Py_c_abs(_Py_c_diff(a, b));
1180
1181 return (((diff <= rel_tol * _Py_c_abs(b)) ||
1182 (diff <= rel_tol * _Py_c_abs(a))) ||
1183 (diff <= abs_tol));
1184}
Guido van Rossum71aa32f1996-01-12 01:34:57 +00001185
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +00001186PyDoc_STRVAR(module_doc,
Tim Peters14e26402001-02-20 20:15:19 +00001187"This module is always available. It provides access to mathematical\n"
Martin v. Löwis14f8b4c2002-06-13 20:33:02 +00001188"functions for complex numbers.");
Guido van Rossumc6e22901998-12-04 19:26:43 +00001189
Roger E. Masse24070ca1996-12-09 22:59:53 +00001190static PyMethodDef cmath_methods[] = {
Brett Cannonb0fc4902014-10-14 17:37:02 -04001191 CMATH_ACOS_METHODDEF
1192 CMATH_ACOSH_METHODDEF
1193 CMATH_ASIN_METHODDEF
1194 CMATH_ASINH_METHODDEF
1195 CMATH_ATAN_METHODDEF
1196 CMATH_ATANH_METHODDEF
1197 CMATH_COS_METHODDEF
1198 CMATH_COSH_METHODDEF
1199 CMATH_EXP_METHODDEF
Tal Einatd5519ed2015-05-31 22:05:00 +03001200 CMATH_ISCLOSE_METHODDEF
Brett Cannonb0fc4902014-10-14 17:37:02 -04001201 CMATH_ISFINITE_METHODDEF
1202 CMATH_ISINF_METHODDEF
1203 CMATH_ISNAN_METHODDEF
1204 CMATH_LOG_METHODDEF
1205 CMATH_LOG10_METHODDEF
1206 CMATH_PHASE_METHODDEF
1207 CMATH_POLAR_METHODDEF
1208 CMATH_RECT_METHODDEF
1209 CMATH_SIN_METHODDEF
1210 CMATH_SINH_METHODDEF
1211 CMATH_SQRT_METHODDEF
1212 CMATH_TAN_METHODDEF
1213 CMATH_TANH_METHODDEF
1214 {NULL, NULL} /* sentinel */
Guido van Rossum71aa32f1996-01-12 01:34:57 +00001215};
1216
Martin v. Löwis1a214512008-06-11 05:26:20 +00001217
1218static struct PyModuleDef cmathmodule = {
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001219 PyModuleDef_HEAD_INIT,
1220 "cmath",
1221 module_doc,
1222 -1,
1223 cmath_methods,
1224 NULL,
1225 NULL,
1226 NULL,
1227 NULL
Martin v. Löwis1a214512008-06-11 05:26:20 +00001228};
1229
Mark Hammondfe51c6d2002-08-02 02:27:13 +00001230PyMODINIT_FUNC
Martin v. Löwis1a214512008-06-11 05:26:20 +00001231PyInit_cmath(void)
Guido van Rossum71aa32f1996-01-12 01:34:57 +00001232{
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001233 PyObject *m;
Tim Peters14e26402001-02-20 20:15:19 +00001234
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001235 m = PyModule_Create(&cmathmodule);
1236 if (m == NULL)
1237 return NULL;
Fred Drakef4e34842002-04-01 03:45:06 +00001238
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001239 PyModule_AddObject(m, "pi",
1240 PyFloat_FromDouble(Py_MATH_PI));
1241 PyModule_AddObject(m, "e", PyFloat_FromDouble(Py_MATH_E));
Christian Heimesa342c012008-04-20 21:01:16 +00001242
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001243 /* initialize special value tables */
Christian Heimesa342c012008-04-20 21:01:16 +00001244
1245#define INIT_SPECIAL_VALUES(NAME, BODY) { Py_complex* p = (Py_complex*)NAME; BODY }
1246#define C(REAL, IMAG) p->real = REAL; p->imag = IMAG; ++p;
1247
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001248 INIT_SPECIAL_VALUES(acos_special_values, {
1249 C(P34,INF) C(P,INF) C(P,INF) C(P,-INF) C(P,-INF) C(P34,-INF) C(N,INF)
1250 C(P12,INF) C(U,U) C(U,U) C(U,U) C(U,U) C(P12,-INF) C(N,N)
1251 C(P12,INF) C(U,U) C(P12,0.) C(P12,-0.) C(U,U) C(P12,-INF) C(P12,N)
1252 C(P12,INF) C(U,U) C(P12,0.) C(P12,-0.) C(U,U) C(P12,-INF) C(P12,N)
1253 C(P12,INF) C(U,U) C(U,U) C(U,U) C(U,U) C(P12,-INF) C(N,N)
1254 C(P14,INF) C(0.,INF) C(0.,INF) C(0.,-INF) C(0.,-INF) C(P14,-INF) C(N,INF)
1255 C(N,INF) C(N,N) C(N,N) C(N,N) C(N,N) C(N,-INF) C(N,N)
1256 })
Christian Heimesa342c012008-04-20 21:01:16 +00001257
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001258 INIT_SPECIAL_VALUES(acosh_special_values, {
1259 C(INF,-P34) C(INF,-P) C(INF,-P) C(INF,P) C(INF,P) C(INF,P34) C(INF,N)
1260 C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N)
1261 C(INF,-P12) C(U,U) C(0.,-P12) C(0.,P12) C(U,U) C(INF,P12) C(N,N)
1262 C(INF,-P12) C(U,U) C(0.,-P12) C(0.,P12) C(U,U) C(INF,P12) C(N,N)
1263 C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N)
1264 C(INF,-P14) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,P14) C(INF,N)
1265 C(INF,N) C(N,N) C(N,N) C(N,N) C(N,N) C(INF,N) C(N,N)
1266 })
Christian Heimesa342c012008-04-20 21:01:16 +00001267
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001268 INIT_SPECIAL_VALUES(asinh_special_values, {
1269 C(-INF,-P14) C(-INF,-0.) C(-INF,-0.) C(-INF,0.) C(-INF,0.) C(-INF,P14) C(-INF,N)
1270 C(-INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(-INF,P12) C(N,N)
1271 C(-INF,-P12) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(-INF,P12) C(N,N)
1272 C(INF,-P12) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(INF,P12) C(N,N)
1273 C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N)
1274 C(INF,-P14) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,P14) C(INF,N)
1275 C(INF,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(INF,N) C(N,N)
1276 })
Christian Heimesa342c012008-04-20 21:01:16 +00001277
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001278 INIT_SPECIAL_VALUES(atanh_special_values, {
1279 C(-0.,-P12) C(-0.,-P12) C(-0.,-P12) C(-0.,P12) C(-0.,P12) C(-0.,P12) C(-0.,N)
1280 C(-0.,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(-0.,P12) C(N,N)
1281 C(-0.,-P12) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(-0.,P12) C(-0.,N)
1282 C(0.,-P12) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,P12) C(0.,N)
1283 C(0.,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(0.,P12) C(N,N)
1284 C(0.,-P12) C(0.,-P12) C(0.,-P12) C(0.,P12) C(0.,P12) C(0.,P12) C(0.,N)
1285 C(0.,-P12) C(N,N) C(N,N) C(N,N) C(N,N) C(0.,P12) C(N,N)
1286 })
Christian Heimesa342c012008-04-20 21:01:16 +00001287
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001288 INIT_SPECIAL_VALUES(cosh_special_values, {
1289 C(INF,N) C(U,U) C(INF,0.) C(INF,-0.) C(U,U) C(INF,N) C(INF,N)
1290 C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N)
1291 C(N,0.) C(U,U) C(1.,0.) C(1.,-0.) C(U,U) C(N,0.) C(N,0.)
1292 C(N,0.) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(N,0.) C(N,0.)
1293 C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N)
1294 C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N)
1295 C(N,N) C(N,N) C(N,0.) C(N,0.) C(N,N) C(N,N) C(N,N)
1296 })
Christian Heimesa342c012008-04-20 21:01:16 +00001297
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001298 INIT_SPECIAL_VALUES(exp_special_values, {
1299 C(0.,0.) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,0.) C(0.,0.)
1300 C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N)
1301 C(N,N) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(N,N) C(N,N)
1302 C(N,N) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(N,N) C(N,N)
1303 C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N)
1304 C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N)
1305 C(N,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(N,N) C(N,N)
1306 })
Christian Heimesa342c012008-04-20 21:01:16 +00001307
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001308 INIT_SPECIAL_VALUES(log_special_values, {
1309 C(INF,-P34) C(INF,-P) C(INF,-P) C(INF,P) C(INF,P) C(INF,P34) C(INF,N)
1310 C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N)
1311 C(INF,-P12) C(U,U) C(-INF,-P) C(-INF,P) C(U,U) C(INF,P12) C(N,N)
1312 C(INF,-P12) C(U,U) C(-INF,-0.) C(-INF,0.) C(U,U) C(INF,P12) C(N,N)
1313 C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N)
1314 C(INF,-P14) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,P14) C(INF,N)
1315 C(INF,N) C(N,N) C(N,N) C(N,N) C(N,N) C(INF,N) C(N,N)
1316 })
Christian Heimesa342c012008-04-20 21:01:16 +00001317
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001318 INIT_SPECIAL_VALUES(sinh_special_values, {
1319 C(INF,N) C(U,U) C(-INF,-0.) C(-INF,0.) C(U,U) C(INF,N) C(INF,N)
1320 C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N)
1321 C(0.,N) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(0.,N) C(0.,N)
1322 C(0.,N) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,N) C(0.,N)
1323 C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N)
1324 C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N)
1325 C(N,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(N,N) C(N,N)
1326 })
Christian Heimesa342c012008-04-20 21:01:16 +00001327
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001328 INIT_SPECIAL_VALUES(sqrt_special_values, {
1329 C(INF,-INF) C(0.,-INF) C(0.,-INF) C(0.,INF) C(0.,INF) C(INF,INF) C(N,INF)
1330 C(INF,-INF) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,INF) C(N,N)
1331 C(INF,-INF) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(INF,INF) C(N,N)
1332 C(INF,-INF) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(INF,INF) C(N,N)
1333 C(INF,-INF) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,INF) C(N,N)
1334 C(INF,-INF) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,INF) C(INF,N)
1335 C(INF,-INF) C(N,N) C(N,N) C(N,N) C(N,N) C(INF,INF) C(N,N)
1336 })
Christian Heimesa342c012008-04-20 21:01:16 +00001337
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001338 INIT_SPECIAL_VALUES(tanh_special_values, {
1339 C(-1.,0.) C(U,U) C(-1.,-0.) C(-1.,0.) C(U,U) C(-1.,0.) C(-1.,0.)
1340 C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N)
1341 C(N,N) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(N,N) C(N,N)
1342 C(N,N) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(N,N) C(N,N)
1343 C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N)
1344 C(1.,0.) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(1.,0.) C(1.,0.)
1345 C(N,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(N,N) C(N,N)
1346 })
Christian Heimesa342c012008-04-20 21:01:16 +00001347
Antoine Pitrouf95a1b32010-05-09 15:52:27 +00001348 INIT_SPECIAL_VALUES(rect_special_values, {
1349 C(INF,N) C(U,U) C(-INF,0.) C(-INF,-0.) C(U,U) C(INF,N) C(INF,N)
1350 C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N)
1351 C(0.,0.) C(U,U) C(-0.,0.) C(-0.,-0.) C(U,U) C(0.,0.) C(0.,0.)
1352 C(0.,0.) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,0.) C(0.,0.)
1353 C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N)
1354 C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N)
1355 C(N,N) C(N,N) C(N,0.) C(N,0.) C(N,N) C(N,N) C(N,N)
1356 })
1357 return m;
Guido van Rossum71aa32f1996-01-12 01:34:57 +00001358}