Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 1 | /* Complex math module */ |
| 2 | |
| 3 | /* much code borrowed from mathmodule.c */ |
| 4 | |
Roger E. Masse | 24070ca | 1996-12-09 22:59:53 +0000 | [diff] [blame] | 5 | #include "Python.h" |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 6 | /* we need DBL_MAX, DBL_MIN, DBL_EPSILON, DBL_MANT_DIG and FLT_RADIX from |
| 7 | float.h. We assume that FLT_RADIX is either 2 or 16. */ |
| 8 | #include <float.h> |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 9 | |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 10 | #if (FLT_RADIX != 2 && FLT_RADIX != 16) |
| 11 | #error "Modules/cmathmodule.c expects FLT_RADIX to be 2 or 16" |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 12 | #endif |
| 13 | |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 14 | #ifndef M_LN2 |
| 15 | #define M_LN2 (0.6931471805599453094) /* natural log of 2 */ |
| 16 | #endif |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 17 | |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 18 | #ifndef M_LN10 |
| 19 | #define M_LN10 (2.302585092994045684) /* natural log of 10 */ |
| 20 | #endif |
| 21 | |
| 22 | /* |
| 23 | CM_LARGE_DOUBLE is used to avoid spurious overflow in the sqrt, log, |
| 24 | inverse trig and inverse hyperbolic trig functions. Its log is used in the |
Ezio Melotti | 1392500 | 2011-03-16 11:05:33 +0200 | [diff] [blame] | 25 | evaluation of exp, cos, cosh, sin, sinh, tan, and tanh to avoid unnecessary |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 26 | overflow. |
| 27 | */ |
| 28 | |
| 29 | #define CM_LARGE_DOUBLE (DBL_MAX/4.) |
| 30 | #define CM_SQRT_LARGE_DOUBLE (sqrt(CM_LARGE_DOUBLE)) |
| 31 | #define CM_LOG_LARGE_DOUBLE (log(CM_LARGE_DOUBLE)) |
| 32 | #define CM_SQRT_DBL_MIN (sqrt(DBL_MIN)) |
| 33 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 34 | /* |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 35 | CM_SCALE_UP is an odd integer chosen such that multiplication by |
| 36 | 2**CM_SCALE_UP is sufficient to turn a subnormal into a normal. |
| 37 | CM_SCALE_DOWN is (-(CM_SCALE_UP+1)/2). These scalings are used to compute |
| 38 | square roots accurately when the real and imaginary parts of the argument |
| 39 | are subnormal. |
| 40 | */ |
| 41 | |
| 42 | #if FLT_RADIX==2 |
| 43 | #define CM_SCALE_UP (2*(DBL_MANT_DIG/2) + 1) |
| 44 | #elif FLT_RADIX==16 |
| 45 | #define CM_SCALE_UP (4*DBL_MANT_DIG+1) |
| 46 | #endif |
| 47 | #define CM_SCALE_DOWN (-(CM_SCALE_UP+1)/2) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 48 | |
| 49 | /* forward declarations */ |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 50 | static Py_complex c_asinh(Py_complex); |
| 51 | static Py_complex c_atanh(Py_complex); |
| 52 | static Py_complex c_cosh(Py_complex); |
| 53 | static Py_complex c_sinh(Py_complex); |
Jeremy Hylton | 938ace6 | 2002-07-17 16:30:39 +0000 | [diff] [blame] | 54 | static Py_complex c_sqrt(Py_complex); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 55 | static Py_complex c_tanh(Py_complex); |
Raymond Hettinger | b67ad7e | 2004-06-14 07:40:10 +0000 | [diff] [blame] | 56 | static PyObject * math_error(void); |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 57 | |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 58 | /* Code to deal with special values (infinities, NaNs, etc.). */ |
| 59 | |
| 60 | /* special_type takes a double and returns an integer code indicating |
| 61 | the type of the double as follows: |
| 62 | */ |
| 63 | |
| 64 | enum special_types { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 65 | ST_NINF, /* 0, negative infinity */ |
| 66 | ST_NEG, /* 1, negative finite number (nonzero) */ |
| 67 | ST_NZERO, /* 2, -0. */ |
| 68 | ST_PZERO, /* 3, +0. */ |
| 69 | ST_POS, /* 4, positive finite number (nonzero) */ |
| 70 | ST_PINF, /* 5, positive infinity */ |
| 71 | ST_NAN /* 6, Not a Number */ |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 72 | }; |
| 73 | |
| 74 | static enum special_types |
| 75 | special_type(double d) |
| 76 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 77 | if (Py_IS_FINITE(d)) { |
| 78 | if (d != 0) { |
| 79 | if (copysign(1., d) == 1.) |
| 80 | return ST_POS; |
| 81 | else |
| 82 | return ST_NEG; |
| 83 | } |
| 84 | else { |
| 85 | if (copysign(1., d) == 1.) |
| 86 | return ST_PZERO; |
| 87 | else |
| 88 | return ST_NZERO; |
| 89 | } |
| 90 | } |
| 91 | if (Py_IS_NAN(d)) |
| 92 | return ST_NAN; |
| 93 | if (copysign(1., d) == 1.) |
| 94 | return ST_PINF; |
| 95 | else |
| 96 | return ST_NINF; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 97 | } |
| 98 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 99 | #define SPECIAL_VALUE(z, table) \ |
| 100 | if (!Py_IS_FINITE((z).real) || !Py_IS_FINITE((z).imag)) { \ |
| 101 | errno = 0; \ |
| 102 | return table[special_type((z).real)] \ |
| 103 | [special_type((z).imag)]; \ |
| 104 | } |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 105 | |
| 106 | #define P Py_MATH_PI |
| 107 | #define P14 0.25*Py_MATH_PI |
| 108 | #define P12 0.5*Py_MATH_PI |
| 109 | #define P34 0.75*Py_MATH_PI |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 110 | #define INF Py_HUGE_VAL |
| 111 | #define N Py_NAN |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 112 | #define U -9.5426319407711027e33 /* unlikely value, used as placeholder */ |
| 113 | |
| 114 | /* First, the C functions that do the real work. Each of the c_* |
| 115 | functions computes and returns the C99 Annex G recommended result |
| 116 | and also sets errno as follows: errno = 0 if no floating-point |
| 117 | exception is associated with the result; errno = EDOM if C99 Annex |
| 118 | G recommends raising divide-by-zero or invalid for this result; and |
| 119 | errno = ERANGE where the overflow floating-point signal should be |
| 120 | raised. |
| 121 | */ |
| 122 | |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 123 | static Py_complex acos_special_values[7][7]; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 124 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 125 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 126 | c_acos(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 127 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 128 | Py_complex s1, s2, r; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 129 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 130 | SPECIAL_VALUE(z, acos_special_values); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 131 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 132 | if (fabs(z.real) > CM_LARGE_DOUBLE || fabs(z.imag) > CM_LARGE_DOUBLE) { |
| 133 | /* avoid unnecessary overflow for large arguments */ |
| 134 | r.real = atan2(fabs(z.imag), z.real); |
| 135 | /* split into cases to make sure that the branch cut has the |
| 136 | correct continuity on systems with unsigned zeros */ |
| 137 | if (z.real < 0.) { |
| 138 | r.imag = -copysign(log(hypot(z.real/2., z.imag/2.)) + |
| 139 | M_LN2*2., z.imag); |
| 140 | } else { |
| 141 | r.imag = copysign(log(hypot(z.real/2., z.imag/2.)) + |
| 142 | M_LN2*2., -z.imag); |
| 143 | } |
| 144 | } else { |
| 145 | s1.real = 1.-z.real; |
| 146 | s1.imag = -z.imag; |
| 147 | s1 = c_sqrt(s1); |
| 148 | s2.real = 1.+z.real; |
| 149 | s2.imag = z.imag; |
| 150 | s2 = c_sqrt(s2); |
| 151 | r.real = 2.*atan2(s1.real, s2.real); |
| 152 | r.imag = asinh(s2.real*s1.imag - s2.imag*s1.real); |
| 153 | } |
| 154 | errno = 0; |
| 155 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 156 | } |
| 157 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 158 | PyDoc_STRVAR(c_acos_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 159 | "acos(x)\n" |
| 160 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 161 | "Return the arc cosine of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 162 | |
| 163 | |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 164 | static Py_complex acosh_special_values[7][7]; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 165 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 166 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 167 | c_acosh(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 168 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 169 | Py_complex s1, s2, r; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 170 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 171 | SPECIAL_VALUE(z, acosh_special_values); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 172 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 173 | if (fabs(z.real) > CM_LARGE_DOUBLE || fabs(z.imag) > CM_LARGE_DOUBLE) { |
| 174 | /* avoid unnecessary overflow for large arguments */ |
| 175 | r.real = log(hypot(z.real/2., z.imag/2.)) + M_LN2*2.; |
| 176 | r.imag = atan2(z.imag, z.real); |
| 177 | } else { |
| 178 | s1.real = z.real - 1.; |
| 179 | s1.imag = z.imag; |
| 180 | s1 = c_sqrt(s1); |
| 181 | s2.real = z.real + 1.; |
| 182 | s2.imag = z.imag; |
| 183 | s2 = c_sqrt(s2); |
| 184 | r.real = asinh(s1.real*s2.real + s1.imag*s2.imag); |
| 185 | r.imag = 2.*atan2(s1.imag, s2.real); |
| 186 | } |
| 187 | errno = 0; |
| 188 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 189 | } |
| 190 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 191 | PyDoc_STRVAR(c_acosh_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 192 | "acosh(x)\n" |
| 193 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 194 | "Return the hyperbolic arccosine of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 195 | |
| 196 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 197 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 198 | c_asin(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 199 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 200 | /* asin(z) = -i asinh(iz) */ |
| 201 | Py_complex s, r; |
| 202 | s.real = -z.imag; |
| 203 | s.imag = z.real; |
| 204 | s = c_asinh(s); |
| 205 | r.real = s.imag; |
| 206 | r.imag = -s.real; |
| 207 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 208 | } |
| 209 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 210 | PyDoc_STRVAR(c_asin_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 211 | "asin(x)\n" |
| 212 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 213 | "Return the arc sine of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 214 | |
| 215 | |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 216 | static Py_complex asinh_special_values[7][7]; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 217 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 218 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 219 | c_asinh(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 220 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 221 | Py_complex s1, s2, r; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 222 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 223 | SPECIAL_VALUE(z, asinh_special_values); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 224 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 225 | if (fabs(z.real) > CM_LARGE_DOUBLE || fabs(z.imag) > CM_LARGE_DOUBLE) { |
| 226 | if (z.imag >= 0.) { |
| 227 | r.real = copysign(log(hypot(z.real/2., z.imag/2.)) + |
| 228 | M_LN2*2., z.real); |
| 229 | } else { |
| 230 | r.real = -copysign(log(hypot(z.real/2., z.imag/2.)) + |
| 231 | M_LN2*2., -z.real); |
| 232 | } |
| 233 | r.imag = atan2(z.imag, fabs(z.real)); |
| 234 | } else { |
| 235 | s1.real = 1.+z.imag; |
| 236 | s1.imag = -z.real; |
| 237 | s1 = c_sqrt(s1); |
| 238 | s2.real = 1.-z.imag; |
| 239 | s2.imag = z.real; |
| 240 | s2 = c_sqrt(s2); |
| 241 | r.real = asinh(s1.real*s2.imag-s2.real*s1.imag); |
| 242 | r.imag = atan2(z.imag, s1.real*s2.real-s1.imag*s2.imag); |
| 243 | } |
| 244 | errno = 0; |
| 245 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 246 | } |
| 247 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 248 | PyDoc_STRVAR(c_asinh_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 249 | "asinh(x)\n" |
| 250 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 251 | "Return the hyperbolic arc sine of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 252 | |
| 253 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 254 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 255 | c_atan(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 256 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 257 | /* atan(z) = -i atanh(iz) */ |
| 258 | Py_complex s, r; |
| 259 | s.real = -z.imag; |
| 260 | s.imag = z.real; |
| 261 | s = c_atanh(s); |
| 262 | r.real = s.imag; |
| 263 | r.imag = -s.real; |
| 264 | return r; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 265 | } |
| 266 | |
Christian Heimes | e57950f | 2008-04-21 13:08:03 +0000 | [diff] [blame] | 267 | /* Windows screws up atan2 for inf and nan, and alpha Tru64 5.1 doesn't follow |
| 268 | C99 for atan2(0., 0.). */ |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 269 | static double |
| 270 | c_atan2(Py_complex z) |
| 271 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 272 | if (Py_IS_NAN(z.real) || Py_IS_NAN(z.imag)) |
| 273 | return Py_NAN; |
| 274 | if (Py_IS_INFINITY(z.imag)) { |
| 275 | if (Py_IS_INFINITY(z.real)) { |
| 276 | if (copysign(1., z.real) == 1.) |
| 277 | /* atan2(+-inf, +inf) == +-pi/4 */ |
| 278 | return copysign(0.25*Py_MATH_PI, z.imag); |
| 279 | else |
| 280 | /* atan2(+-inf, -inf) == +-pi*3/4 */ |
| 281 | return copysign(0.75*Py_MATH_PI, z.imag); |
| 282 | } |
| 283 | /* atan2(+-inf, x) == +-pi/2 for finite x */ |
| 284 | return copysign(0.5*Py_MATH_PI, z.imag); |
| 285 | } |
| 286 | if (Py_IS_INFINITY(z.real) || z.imag == 0.) { |
| 287 | if (copysign(1., z.real) == 1.) |
| 288 | /* atan2(+-y, +inf) = atan2(+-0, +x) = +-0. */ |
| 289 | return copysign(0., z.imag); |
| 290 | else |
| 291 | /* atan2(+-y, -inf) = atan2(+-0., -x) = +-pi. */ |
| 292 | return copysign(Py_MATH_PI, z.imag); |
| 293 | } |
| 294 | return atan2(z.imag, z.real); |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 295 | } |
| 296 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 297 | PyDoc_STRVAR(c_atan_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 298 | "atan(x)\n" |
| 299 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 300 | "Return the arc tangent of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 301 | |
| 302 | |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 303 | static Py_complex atanh_special_values[7][7]; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 304 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 305 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 306 | c_atanh(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 307 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 308 | Py_complex r; |
| 309 | double ay, h; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 310 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 311 | SPECIAL_VALUE(z, atanh_special_values); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 312 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 313 | /* Reduce to case where z.real >= 0., using atanh(z) = -atanh(-z). */ |
| 314 | if (z.real < 0.) { |
| 315 | return c_neg(c_atanh(c_neg(z))); |
| 316 | } |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 317 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 318 | ay = fabs(z.imag); |
| 319 | if (z.real > CM_SQRT_LARGE_DOUBLE || ay > CM_SQRT_LARGE_DOUBLE) { |
| 320 | /* |
| 321 | if abs(z) is large then we use the approximation |
| 322 | atanh(z) ~ 1/z +/- i*pi/2 (+/- depending on the sign |
| 323 | of z.imag) |
| 324 | */ |
| 325 | h = hypot(z.real/2., z.imag/2.); /* safe from overflow */ |
| 326 | r.real = z.real/4./h/h; |
| 327 | /* the two negations in the next line cancel each other out |
| 328 | except when working with unsigned zeros: they're there to |
| 329 | ensure that the branch cut has the correct continuity on |
| 330 | systems that don't support signed zeros */ |
| 331 | r.imag = -copysign(Py_MATH_PI/2., -z.imag); |
| 332 | errno = 0; |
| 333 | } else if (z.real == 1. && ay < CM_SQRT_DBL_MIN) { |
| 334 | /* C99 standard says: atanh(1+/-0.) should be inf +/- 0i */ |
| 335 | if (ay == 0.) { |
| 336 | r.real = INF; |
| 337 | r.imag = z.imag; |
| 338 | errno = EDOM; |
| 339 | } else { |
| 340 | r.real = -log(sqrt(ay)/sqrt(hypot(ay, 2.))); |
| 341 | r.imag = copysign(atan2(2., -ay)/2, z.imag); |
| 342 | errno = 0; |
| 343 | } |
| 344 | } else { |
| 345 | r.real = log1p(4.*z.real/((1-z.real)*(1-z.real) + ay*ay))/4.; |
| 346 | r.imag = -atan2(-2.*z.imag, (1-z.real)*(1+z.real) - ay*ay)/2.; |
| 347 | errno = 0; |
| 348 | } |
| 349 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 350 | } |
| 351 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 352 | PyDoc_STRVAR(c_atanh_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 353 | "atanh(x)\n" |
| 354 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 355 | "Return the hyperbolic arc tangent of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 356 | |
| 357 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 358 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 359 | c_cos(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 360 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 361 | /* cos(z) = cosh(iz) */ |
| 362 | Py_complex r; |
| 363 | r.real = -z.imag; |
| 364 | r.imag = z.real; |
| 365 | r = c_cosh(r); |
| 366 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 367 | } |
| 368 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 369 | PyDoc_STRVAR(c_cos_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 370 | "cos(x)\n" |
Mark Dickinson | 1bd2e29 | 2009-02-28 15:53:24 +0000 | [diff] [blame] | 371 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 372 | "Return the cosine of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 373 | |
| 374 | |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 375 | /* cosh(infinity + i*y) needs to be dealt with specially */ |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 376 | static Py_complex cosh_special_values[7][7]; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 377 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 378 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 379 | c_cosh(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 380 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 381 | Py_complex r; |
| 382 | double x_minus_one; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 383 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 384 | /* special treatment for cosh(+/-inf + iy) if y is not a NaN */ |
| 385 | if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) { |
| 386 | if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) && |
| 387 | (z.imag != 0.)) { |
| 388 | if (z.real > 0) { |
| 389 | r.real = copysign(INF, cos(z.imag)); |
| 390 | r.imag = copysign(INF, sin(z.imag)); |
| 391 | } |
| 392 | else { |
| 393 | r.real = copysign(INF, cos(z.imag)); |
| 394 | r.imag = -copysign(INF, sin(z.imag)); |
| 395 | } |
| 396 | } |
| 397 | else { |
| 398 | r = cosh_special_values[special_type(z.real)] |
| 399 | [special_type(z.imag)]; |
| 400 | } |
| 401 | /* need to set errno = EDOM if y is +/- infinity and x is not |
| 402 | a NaN */ |
| 403 | if (Py_IS_INFINITY(z.imag) && !Py_IS_NAN(z.real)) |
| 404 | errno = EDOM; |
| 405 | else |
| 406 | errno = 0; |
| 407 | return r; |
| 408 | } |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 409 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 410 | if (fabs(z.real) > CM_LOG_LARGE_DOUBLE) { |
| 411 | /* deal correctly with cases where cosh(z.real) overflows but |
| 412 | cosh(z) does not. */ |
| 413 | x_minus_one = z.real - copysign(1., z.real); |
| 414 | r.real = cos(z.imag) * cosh(x_minus_one) * Py_MATH_E; |
| 415 | r.imag = sin(z.imag) * sinh(x_minus_one) * Py_MATH_E; |
| 416 | } else { |
| 417 | r.real = cos(z.imag) * cosh(z.real); |
| 418 | r.imag = sin(z.imag) * sinh(z.real); |
| 419 | } |
| 420 | /* detect overflow, and set errno accordingly */ |
| 421 | if (Py_IS_INFINITY(r.real) || Py_IS_INFINITY(r.imag)) |
| 422 | errno = ERANGE; |
| 423 | else |
| 424 | errno = 0; |
| 425 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 426 | } |
| 427 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 428 | PyDoc_STRVAR(c_cosh_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 429 | "cosh(x)\n" |
Mark Dickinson | 1bd2e29 | 2009-02-28 15:53:24 +0000 | [diff] [blame] | 430 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 431 | "Return the hyperbolic cosine of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 432 | |
| 433 | |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 434 | /* exp(infinity + i*y) and exp(-infinity + i*y) need special treatment for |
| 435 | finite y */ |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 436 | static Py_complex exp_special_values[7][7]; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 437 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 438 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 439 | c_exp(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 440 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 441 | Py_complex r; |
| 442 | double l; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 443 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 444 | if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) { |
| 445 | if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) |
| 446 | && (z.imag != 0.)) { |
| 447 | if (z.real > 0) { |
| 448 | r.real = copysign(INF, cos(z.imag)); |
| 449 | r.imag = copysign(INF, sin(z.imag)); |
| 450 | } |
| 451 | else { |
| 452 | r.real = copysign(0., cos(z.imag)); |
| 453 | r.imag = copysign(0., sin(z.imag)); |
| 454 | } |
| 455 | } |
| 456 | else { |
| 457 | r = exp_special_values[special_type(z.real)] |
| 458 | [special_type(z.imag)]; |
| 459 | } |
| 460 | /* need to set errno = EDOM if y is +/- infinity and x is not |
| 461 | a NaN and not -infinity */ |
| 462 | if (Py_IS_INFINITY(z.imag) && |
| 463 | (Py_IS_FINITE(z.real) || |
| 464 | (Py_IS_INFINITY(z.real) && z.real > 0))) |
| 465 | errno = EDOM; |
| 466 | else |
| 467 | errno = 0; |
| 468 | return r; |
| 469 | } |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 470 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 471 | if (z.real > CM_LOG_LARGE_DOUBLE) { |
| 472 | l = exp(z.real-1.); |
| 473 | r.real = l*cos(z.imag)*Py_MATH_E; |
| 474 | r.imag = l*sin(z.imag)*Py_MATH_E; |
| 475 | } else { |
| 476 | l = exp(z.real); |
| 477 | r.real = l*cos(z.imag); |
| 478 | r.imag = l*sin(z.imag); |
| 479 | } |
| 480 | /* detect overflow, and set errno accordingly */ |
| 481 | if (Py_IS_INFINITY(r.real) || Py_IS_INFINITY(r.imag)) |
| 482 | errno = ERANGE; |
| 483 | else |
| 484 | errno = 0; |
| 485 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 486 | } |
| 487 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 488 | PyDoc_STRVAR(c_exp_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 489 | "exp(x)\n" |
| 490 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 491 | "Return the exponential value e**x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 492 | |
| 493 | |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 494 | static Py_complex log_special_values[7][7]; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 495 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 496 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 497 | c_log(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 498 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 499 | /* |
| 500 | The usual formula for the real part is log(hypot(z.real, z.imag)). |
| 501 | There are four situations where this formula is potentially |
| 502 | problematic: |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 503 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 504 | (1) the absolute value of z is subnormal. Then hypot is subnormal, |
| 505 | so has fewer than the usual number of bits of accuracy, hence may |
| 506 | have large relative error. This then gives a large absolute error |
| 507 | in the log. This can be solved by rescaling z by a suitable power |
| 508 | of 2. |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 509 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 510 | (2) the absolute value of z is greater than DBL_MAX (e.g. when both |
| 511 | z.real and z.imag are within a factor of 1/sqrt(2) of DBL_MAX) |
| 512 | Again, rescaling solves this. |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 513 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 514 | (3) the absolute value of z is close to 1. In this case it's |
| 515 | difficult to achieve good accuracy, at least in part because a |
| 516 | change of 1ulp in the real or imaginary part of z can result in a |
| 517 | change of billions of ulps in the correctly rounded answer. |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 518 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 519 | (4) z = 0. The simplest thing to do here is to call the |
| 520 | floating-point log with an argument of 0, and let its behaviour |
| 521 | (returning -infinity, signaling a floating-point exception, setting |
| 522 | errno, or whatever) determine that of c_log. So the usual formula |
| 523 | is fine here. |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 524 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 525 | */ |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 526 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 527 | Py_complex r; |
| 528 | double ax, ay, am, an, h; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 529 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 530 | SPECIAL_VALUE(z, log_special_values); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 531 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 532 | ax = fabs(z.real); |
| 533 | ay = fabs(z.imag); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 534 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 535 | if (ax > CM_LARGE_DOUBLE || ay > CM_LARGE_DOUBLE) { |
| 536 | r.real = log(hypot(ax/2., ay/2.)) + M_LN2; |
| 537 | } else if (ax < DBL_MIN && ay < DBL_MIN) { |
| 538 | if (ax > 0. || ay > 0.) { |
| 539 | /* catch cases where hypot(ax, ay) is subnormal */ |
| 540 | r.real = log(hypot(ldexp(ax, DBL_MANT_DIG), |
| 541 | ldexp(ay, DBL_MANT_DIG))) - DBL_MANT_DIG*M_LN2; |
| 542 | } |
| 543 | else { |
| 544 | /* log(+/-0. +/- 0i) */ |
| 545 | r.real = -INF; |
| 546 | r.imag = atan2(z.imag, z.real); |
| 547 | errno = EDOM; |
| 548 | return r; |
| 549 | } |
| 550 | } else { |
| 551 | h = hypot(ax, ay); |
| 552 | if (0.71 <= h && h <= 1.73) { |
| 553 | am = ax > ay ? ax : ay; /* max(ax, ay) */ |
| 554 | an = ax > ay ? ay : ax; /* min(ax, ay) */ |
| 555 | r.real = log1p((am-1)*(am+1)+an*an)/2.; |
| 556 | } else { |
| 557 | r.real = log(h); |
| 558 | } |
| 559 | } |
| 560 | r.imag = atan2(z.imag, z.real); |
| 561 | errno = 0; |
| 562 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 563 | } |
| 564 | |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 565 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 566 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 567 | c_log10(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 568 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 569 | Py_complex r; |
| 570 | int errno_save; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 571 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 572 | r = c_log(z); |
| 573 | errno_save = errno; /* just in case the divisions affect errno */ |
| 574 | r.real = r.real / M_LN10; |
| 575 | r.imag = r.imag / M_LN10; |
| 576 | errno = errno_save; |
| 577 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 578 | } |
| 579 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 580 | PyDoc_STRVAR(c_log10_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 581 | "log10(x)\n" |
| 582 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 583 | "Return the base-10 logarithm of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 584 | |
| 585 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 586 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 587 | c_sin(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 588 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 589 | /* sin(z) = -i sin(iz) */ |
| 590 | Py_complex s, r; |
| 591 | s.real = -z.imag; |
| 592 | s.imag = z.real; |
| 593 | s = c_sinh(s); |
| 594 | r.real = s.imag; |
| 595 | r.imag = -s.real; |
| 596 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 597 | } |
| 598 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 599 | PyDoc_STRVAR(c_sin_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 600 | "sin(x)\n" |
| 601 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 602 | "Return the sine of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 603 | |
| 604 | |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 605 | /* sinh(infinity + i*y) needs to be dealt with specially */ |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 606 | static Py_complex sinh_special_values[7][7]; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 607 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 608 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 609 | c_sinh(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 610 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 611 | Py_complex r; |
| 612 | double x_minus_one; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 613 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 614 | /* special treatment for sinh(+/-inf + iy) if y is finite and |
| 615 | nonzero */ |
| 616 | if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) { |
| 617 | if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) |
| 618 | && (z.imag != 0.)) { |
| 619 | if (z.real > 0) { |
| 620 | r.real = copysign(INF, cos(z.imag)); |
| 621 | r.imag = copysign(INF, sin(z.imag)); |
| 622 | } |
| 623 | else { |
| 624 | r.real = -copysign(INF, cos(z.imag)); |
| 625 | r.imag = copysign(INF, sin(z.imag)); |
| 626 | } |
| 627 | } |
| 628 | else { |
| 629 | r = sinh_special_values[special_type(z.real)] |
| 630 | [special_type(z.imag)]; |
| 631 | } |
| 632 | /* need to set errno = EDOM if y is +/- infinity and x is not |
| 633 | a NaN */ |
| 634 | if (Py_IS_INFINITY(z.imag) && !Py_IS_NAN(z.real)) |
| 635 | errno = EDOM; |
| 636 | else |
| 637 | errno = 0; |
| 638 | return r; |
| 639 | } |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 640 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 641 | if (fabs(z.real) > CM_LOG_LARGE_DOUBLE) { |
| 642 | x_minus_one = z.real - copysign(1., z.real); |
| 643 | r.real = cos(z.imag) * sinh(x_minus_one) * Py_MATH_E; |
| 644 | r.imag = sin(z.imag) * cosh(x_minus_one) * Py_MATH_E; |
| 645 | } else { |
| 646 | r.real = cos(z.imag) * sinh(z.real); |
| 647 | r.imag = sin(z.imag) * cosh(z.real); |
| 648 | } |
| 649 | /* detect overflow, and set errno accordingly */ |
| 650 | if (Py_IS_INFINITY(r.real) || Py_IS_INFINITY(r.imag)) |
| 651 | errno = ERANGE; |
| 652 | else |
| 653 | errno = 0; |
| 654 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 655 | } |
| 656 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 657 | PyDoc_STRVAR(c_sinh_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 658 | "sinh(x)\n" |
| 659 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 660 | "Return the hyperbolic sine of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 661 | |
| 662 | |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 663 | static Py_complex sqrt_special_values[7][7]; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 664 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 665 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 666 | c_sqrt(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 667 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 668 | /* |
| 669 | Method: use symmetries to reduce to the case when x = z.real and y |
| 670 | = z.imag are nonnegative. Then the real part of the result is |
| 671 | given by |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 672 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 673 | s = sqrt((x + hypot(x, y))/2) |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 674 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 675 | and the imaginary part is |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 676 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 677 | d = (y/2)/s |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 678 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 679 | If either x or y is very large then there's a risk of overflow in |
| 680 | computation of the expression x + hypot(x, y). We can avoid this |
| 681 | by rewriting the formula for s as: |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 682 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 683 | s = 2*sqrt(x/8 + hypot(x/8, y/8)) |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 684 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 685 | This costs us two extra multiplications/divisions, but avoids the |
| 686 | overhead of checking for x and y large. |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 687 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 688 | If both x and y are subnormal then hypot(x, y) may also be |
| 689 | subnormal, so will lack full precision. We solve this by rescaling |
| 690 | x and y by a sufficiently large power of 2 to ensure that x and y |
| 691 | are normal. |
| 692 | */ |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 693 | |
| 694 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 695 | Py_complex r; |
| 696 | double s,d; |
| 697 | double ax, ay; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 698 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 699 | SPECIAL_VALUE(z, sqrt_special_values); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 700 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 701 | if (z.real == 0. && z.imag == 0.) { |
| 702 | r.real = 0.; |
| 703 | r.imag = z.imag; |
| 704 | return r; |
| 705 | } |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 706 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 707 | ax = fabs(z.real); |
| 708 | ay = fabs(z.imag); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 709 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 710 | if (ax < DBL_MIN && ay < DBL_MIN && (ax > 0. || ay > 0.)) { |
| 711 | /* here we catch cases where hypot(ax, ay) is subnormal */ |
| 712 | ax = ldexp(ax, CM_SCALE_UP); |
| 713 | s = ldexp(sqrt(ax + hypot(ax, ldexp(ay, CM_SCALE_UP))), |
| 714 | CM_SCALE_DOWN); |
| 715 | } else { |
| 716 | ax /= 8.; |
| 717 | s = 2.*sqrt(ax + hypot(ax, ay/8.)); |
| 718 | } |
| 719 | d = ay/(2.*s); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 720 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 721 | if (z.real >= 0.) { |
| 722 | r.real = s; |
| 723 | r.imag = copysign(d, z.imag); |
| 724 | } else { |
| 725 | r.real = d; |
| 726 | r.imag = copysign(s, z.imag); |
| 727 | } |
| 728 | errno = 0; |
| 729 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 730 | } |
| 731 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 732 | PyDoc_STRVAR(c_sqrt_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 733 | "sqrt(x)\n" |
| 734 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 735 | "Return the square root of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 736 | |
| 737 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 738 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 739 | c_tan(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 740 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 741 | /* tan(z) = -i tanh(iz) */ |
| 742 | Py_complex s, r; |
| 743 | s.real = -z.imag; |
| 744 | s.imag = z.real; |
| 745 | s = c_tanh(s); |
| 746 | r.real = s.imag; |
| 747 | r.imag = -s.real; |
| 748 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 749 | } |
| 750 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 751 | PyDoc_STRVAR(c_tan_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 752 | "tan(x)\n" |
| 753 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 754 | "Return the tangent of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 755 | |
| 756 | |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 757 | /* tanh(infinity + i*y) needs to be dealt with specially */ |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 758 | static Py_complex tanh_special_values[7][7]; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 759 | |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 760 | static Py_complex |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 761 | c_tanh(Py_complex z) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 762 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 763 | /* Formula: |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 764 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 765 | tanh(x+iy) = (tanh(x)(1+tan(y)^2) + i tan(y)(1-tanh(x))^2) / |
| 766 | (1+tan(y)^2 tanh(x)^2) |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 767 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 768 | To avoid excessive roundoff error, 1-tanh(x)^2 is better computed |
| 769 | as 1/cosh(x)^2. When abs(x) is large, we approximate 1-tanh(x)^2 |
| 770 | by 4 exp(-2*x) instead, to avoid possible overflow in the |
| 771 | computation of cosh(x). |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 772 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 773 | */ |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 774 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 775 | Py_complex r; |
| 776 | double tx, ty, cx, txty, denom; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 777 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 778 | /* special treatment for tanh(+/-inf + iy) if y is finite and |
| 779 | nonzero */ |
| 780 | if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) { |
| 781 | if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) |
| 782 | && (z.imag != 0.)) { |
| 783 | if (z.real > 0) { |
| 784 | r.real = 1.0; |
| 785 | r.imag = copysign(0., |
| 786 | 2.*sin(z.imag)*cos(z.imag)); |
| 787 | } |
| 788 | else { |
| 789 | r.real = -1.0; |
| 790 | r.imag = copysign(0., |
| 791 | 2.*sin(z.imag)*cos(z.imag)); |
| 792 | } |
| 793 | } |
| 794 | else { |
| 795 | r = tanh_special_values[special_type(z.real)] |
| 796 | [special_type(z.imag)]; |
| 797 | } |
| 798 | /* need to set errno = EDOM if z.imag is +/-infinity and |
| 799 | z.real is finite */ |
| 800 | if (Py_IS_INFINITY(z.imag) && Py_IS_FINITE(z.real)) |
| 801 | errno = EDOM; |
| 802 | else |
| 803 | errno = 0; |
| 804 | return r; |
| 805 | } |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 806 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 807 | /* danger of overflow in 2.*z.imag !*/ |
| 808 | if (fabs(z.real) > CM_LOG_LARGE_DOUBLE) { |
| 809 | r.real = copysign(1., z.real); |
| 810 | r.imag = 4.*sin(z.imag)*cos(z.imag)*exp(-2.*fabs(z.real)); |
| 811 | } else { |
| 812 | tx = tanh(z.real); |
| 813 | ty = tan(z.imag); |
| 814 | cx = 1./cosh(z.real); |
| 815 | txty = tx*ty; |
| 816 | denom = 1. + txty*txty; |
| 817 | r.real = tx*(1.+ty*ty)/denom; |
| 818 | r.imag = ((ty/denom)*cx)*cx; |
| 819 | } |
| 820 | errno = 0; |
| 821 | return r; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 822 | } |
| 823 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 824 | PyDoc_STRVAR(c_tanh_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 825 | "tanh(x)\n" |
| 826 | "\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 827 | "Return the hyperbolic tangent of x."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 828 | |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 829 | |
Raymond Hettinger | b67ad7e | 2004-06-14 07:40:10 +0000 | [diff] [blame] | 830 | static PyObject * |
| 831 | cmath_log(PyObject *self, PyObject *args) |
| 832 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 833 | Py_complex x; |
| 834 | Py_complex y; |
Raymond Hettinger | b67ad7e | 2004-06-14 07:40:10 +0000 | [diff] [blame] | 835 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 836 | if (!PyArg_ParseTuple(args, "D|D", &x, &y)) |
| 837 | return NULL; |
Raymond Hettinger | b67ad7e | 2004-06-14 07:40:10 +0000 | [diff] [blame] | 838 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 839 | errno = 0; |
| 840 | PyFPE_START_PROTECT("complex function", return 0) |
| 841 | x = c_log(x); |
| 842 | if (PyTuple_GET_SIZE(args) == 2) { |
| 843 | y = c_log(y); |
| 844 | x = c_quot(x, y); |
| 845 | } |
| 846 | PyFPE_END_PROTECT(x) |
| 847 | if (errno != 0) |
| 848 | return math_error(); |
| 849 | return PyComplex_FromCComplex(x); |
Raymond Hettinger | b67ad7e | 2004-06-14 07:40:10 +0000 | [diff] [blame] | 850 | } |
| 851 | |
| 852 | PyDoc_STRVAR(cmath_log_doc, |
| 853 | "log(x[, base]) -> the logarithm of x to the given base.\n\ |
| 854 | If the base not specified, returns the natural logarithm (base e) of x."); |
| 855 | |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 856 | |
| 857 | /* And now the glue to make them available from Python: */ |
| 858 | |
Roger E. Masse | 24070ca | 1996-12-09 22:59:53 +0000 | [diff] [blame] | 859 | static PyObject * |
Thomas Wouters | f3f33dc | 2000-07-21 06:00:07 +0000 | [diff] [blame] | 860 | math_error(void) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 861 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 862 | if (errno == EDOM) |
| 863 | PyErr_SetString(PyExc_ValueError, "math domain error"); |
| 864 | else if (errno == ERANGE) |
| 865 | PyErr_SetString(PyExc_OverflowError, "math range error"); |
| 866 | else /* Unexpected math error */ |
| 867 | PyErr_SetFromErrno(PyExc_ValueError); |
| 868 | return NULL; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 869 | } |
| 870 | |
Roger E. Masse | 24070ca | 1996-12-09 22:59:53 +0000 | [diff] [blame] | 871 | static PyObject * |
Peter Schneider-Kamp | f1ca898 | 2000-07-10 09:31:34 +0000 | [diff] [blame] | 872 | math_1(PyObject *args, Py_complex (*func)(Py_complex)) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 873 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 874 | Py_complex x,r ; |
| 875 | if (!PyArg_ParseTuple(args, "D", &x)) |
| 876 | return NULL; |
| 877 | errno = 0; |
| 878 | PyFPE_START_PROTECT("complex function", return 0); |
| 879 | r = (*func)(x); |
| 880 | PyFPE_END_PROTECT(r); |
| 881 | if (errno == EDOM) { |
| 882 | PyErr_SetString(PyExc_ValueError, "math domain error"); |
| 883 | return NULL; |
| 884 | } |
| 885 | else if (errno == ERANGE) { |
| 886 | PyErr_SetString(PyExc_OverflowError, "math range error"); |
| 887 | return NULL; |
| 888 | } |
| 889 | else { |
| 890 | return PyComplex_FromCComplex(r); |
| 891 | } |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 892 | } |
| 893 | |
| 894 | #define FUNC1(stubname, func) \ |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 895 | static PyObject * stubname(PyObject *self, PyObject *args) { \ |
| 896 | return math_1(args, func); \ |
| 897 | } |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 898 | |
| 899 | FUNC1(cmath_acos, c_acos) |
| 900 | FUNC1(cmath_acosh, c_acosh) |
| 901 | FUNC1(cmath_asin, c_asin) |
| 902 | FUNC1(cmath_asinh, c_asinh) |
| 903 | FUNC1(cmath_atan, c_atan) |
| 904 | FUNC1(cmath_atanh, c_atanh) |
| 905 | FUNC1(cmath_cos, c_cos) |
| 906 | FUNC1(cmath_cosh, c_cosh) |
| 907 | FUNC1(cmath_exp, c_exp) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 908 | FUNC1(cmath_log10, c_log10) |
| 909 | FUNC1(cmath_sin, c_sin) |
| 910 | FUNC1(cmath_sinh, c_sinh) |
| 911 | FUNC1(cmath_sqrt, c_sqrt) |
| 912 | FUNC1(cmath_tan, c_tan) |
| 913 | FUNC1(cmath_tanh, c_tanh) |
| 914 | |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 915 | static PyObject * |
| 916 | cmath_phase(PyObject *self, PyObject *args) |
| 917 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 918 | Py_complex z; |
| 919 | double phi; |
| 920 | if (!PyArg_ParseTuple(args, "D:phase", &z)) |
| 921 | return NULL; |
| 922 | errno = 0; |
| 923 | PyFPE_START_PROTECT("arg function", return 0) |
| 924 | phi = c_atan2(z); |
| 925 | PyFPE_END_PROTECT(phi) |
| 926 | if (errno != 0) |
| 927 | return math_error(); |
| 928 | else |
| 929 | return PyFloat_FromDouble(phi); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 930 | } |
| 931 | |
| 932 | PyDoc_STRVAR(cmath_phase_doc, |
| 933 | "phase(z) -> float\n\n\ |
| 934 | Return argument, also known as the phase angle, of a complex."); |
| 935 | |
| 936 | static PyObject * |
| 937 | cmath_polar(PyObject *self, PyObject *args) |
| 938 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 939 | Py_complex z; |
| 940 | double r, phi; |
| 941 | if (!PyArg_ParseTuple(args, "D:polar", &z)) |
| 942 | return NULL; |
| 943 | PyFPE_START_PROTECT("polar function", return 0) |
| 944 | phi = c_atan2(z); /* should not cause any exception */ |
| 945 | r = c_abs(z); /* sets errno to ERANGE on overflow; otherwise 0 */ |
| 946 | PyFPE_END_PROTECT(r) |
| 947 | if (errno != 0) |
| 948 | return math_error(); |
| 949 | else |
| 950 | return Py_BuildValue("dd", r, phi); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 951 | } |
| 952 | |
| 953 | PyDoc_STRVAR(cmath_polar_doc, |
| 954 | "polar(z) -> r: float, phi: float\n\n\ |
| 955 | Convert a complex from rectangular coordinates to polar coordinates. r is\n\ |
| 956 | the distance from 0 and phi the phase angle."); |
| 957 | |
| 958 | /* |
| 959 | rect() isn't covered by the C99 standard, but it's not too hard to |
| 960 | figure out 'spirit of C99' rules for special value handing: |
| 961 | |
| 962 | rect(x, t) should behave like exp(log(x) + it) for positive-signed x |
| 963 | rect(x, t) should behave like -exp(log(-x) + it) for negative-signed x |
| 964 | rect(nan, t) should behave like exp(nan + it), except that rect(nan, 0) |
| 965 | gives nan +- i0 with the sign of the imaginary part unspecified. |
| 966 | |
| 967 | */ |
| 968 | |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 969 | static Py_complex rect_special_values[7][7]; |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 970 | |
| 971 | static PyObject * |
| 972 | cmath_rect(PyObject *self, PyObject *args) |
| 973 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 974 | Py_complex z; |
| 975 | double r, phi; |
| 976 | if (!PyArg_ParseTuple(args, "dd:rect", &r, &phi)) |
| 977 | return NULL; |
| 978 | errno = 0; |
| 979 | PyFPE_START_PROTECT("rect function", return 0) |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 980 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 981 | /* deal with special values */ |
| 982 | if (!Py_IS_FINITE(r) || !Py_IS_FINITE(phi)) { |
| 983 | /* if r is +/-infinity and phi is finite but nonzero then |
| 984 | result is (+-INF +-INF i), but we need to compute cos(phi) |
| 985 | and sin(phi) to figure out the signs. */ |
| 986 | if (Py_IS_INFINITY(r) && (Py_IS_FINITE(phi) |
| 987 | && (phi != 0.))) { |
| 988 | if (r > 0) { |
| 989 | z.real = copysign(INF, cos(phi)); |
| 990 | z.imag = copysign(INF, sin(phi)); |
| 991 | } |
| 992 | else { |
| 993 | z.real = -copysign(INF, cos(phi)); |
| 994 | z.imag = -copysign(INF, sin(phi)); |
| 995 | } |
| 996 | } |
| 997 | else { |
| 998 | z = rect_special_values[special_type(r)] |
| 999 | [special_type(phi)]; |
| 1000 | } |
| 1001 | /* need to set errno = EDOM if r is a nonzero number and phi |
| 1002 | is infinite */ |
| 1003 | if (r != 0. && !Py_IS_NAN(r) && Py_IS_INFINITY(phi)) |
| 1004 | errno = EDOM; |
| 1005 | else |
| 1006 | errno = 0; |
| 1007 | } |
| 1008 | else { |
| 1009 | z.real = r * cos(phi); |
| 1010 | z.imag = r * sin(phi); |
| 1011 | errno = 0; |
| 1012 | } |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 1013 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1014 | PyFPE_END_PROTECT(z) |
| 1015 | if (errno != 0) |
| 1016 | return math_error(); |
| 1017 | else |
| 1018 | return PyComplex_FromCComplex(z); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 1019 | } |
| 1020 | |
| 1021 | PyDoc_STRVAR(cmath_rect_doc, |
| 1022 | "rect(r, phi) -> z: complex\n\n\ |
| 1023 | Convert from polar coordinates to rectangular coordinates."); |
| 1024 | |
| 1025 | static PyObject * |
| 1026 | cmath_isnan(PyObject *self, PyObject *args) |
| 1027 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1028 | Py_complex z; |
| 1029 | if (!PyArg_ParseTuple(args, "D:isnan", &z)) |
| 1030 | return NULL; |
| 1031 | return PyBool_FromLong(Py_IS_NAN(z.real) || Py_IS_NAN(z.imag)); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 1032 | } |
| 1033 | |
| 1034 | PyDoc_STRVAR(cmath_isnan_doc, |
| 1035 | "isnan(z) -> bool\n\ |
| 1036 | Checks if the real or imaginary part of z not a number (NaN)"); |
| 1037 | |
| 1038 | static PyObject * |
| 1039 | cmath_isinf(PyObject *self, PyObject *args) |
| 1040 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1041 | Py_complex z; |
| 1042 | if (!PyArg_ParseTuple(args, "D:isnan", &z)) |
| 1043 | return NULL; |
| 1044 | return PyBool_FromLong(Py_IS_INFINITY(z.real) || |
| 1045 | Py_IS_INFINITY(z.imag)); |
Christian Heimes | 53876d9 | 2008-04-19 00:31:39 +0000 | [diff] [blame] | 1046 | } |
| 1047 | |
| 1048 | PyDoc_STRVAR(cmath_isinf_doc, |
| 1049 | "isinf(z) -> bool\n\ |
| 1050 | Checks if the real or imaginary part of z is infinite."); |
| 1051 | |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 1052 | |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 1053 | PyDoc_STRVAR(module_doc, |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 1054 | "This module is always available. It provides access to mathematical\n" |
Martin v. Löwis | 14f8b4c | 2002-06-13 20:33:02 +0000 | [diff] [blame] | 1055 | "functions for complex numbers."); |
Guido van Rossum | c6e2290 | 1998-12-04 19:26:43 +0000 | [diff] [blame] | 1056 | |
Roger E. Masse | 24070ca | 1996-12-09 22:59:53 +0000 | [diff] [blame] | 1057 | static PyMethodDef cmath_methods[] = { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1058 | {"acos", cmath_acos, METH_VARARGS, c_acos_doc}, |
| 1059 | {"acosh", cmath_acosh, METH_VARARGS, c_acosh_doc}, |
| 1060 | {"asin", cmath_asin, METH_VARARGS, c_asin_doc}, |
| 1061 | {"asinh", cmath_asinh, METH_VARARGS, c_asinh_doc}, |
| 1062 | {"atan", cmath_atan, METH_VARARGS, c_atan_doc}, |
| 1063 | {"atanh", cmath_atanh, METH_VARARGS, c_atanh_doc}, |
| 1064 | {"cos", cmath_cos, METH_VARARGS, c_cos_doc}, |
| 1065 | {"cosh", cmath_cosh, METH_VARARGS, c_cosh_doc}, |
| 1066 | {"exp", cmath_exp, METH_VARARGS, c_exp_doc}, |
| 1067 | {"isinf", cmath_isinf, METH_VARARGS, cmath_isinf_doc}, |
| 1068 | {"isnan", cmath_isnan, METH_VARARGS, cmath_isnan_doc}, |
| 1069 | {"log", cmath_log, METH_VARARGS, cmath_log_doc}, |
| 1070 | {"log10", cmath_log10, METH_VARARGS, c_log10_doc}, |
| 1071 | {"phase", cmath_phase, METH_VARARGS, cmath_phase_doc}, |
| 1072 | {"polar", cmath_polar, METH_VARARGS, cmath_polar_doc}, |
| 1073 | {"rect", cmath_rect, METH_VARARGS, cmath_rect_doc}, |
| 1074 | {"sin", cmath_sin, METH_VARARGS, c_sin_doc}, |
| 1075 | {"sinh", cmath_sinh, METH_VARARGS, c_sinh_doc}, |
| 1076 | {"sqrt", cmath_sqrt, METH_VARARGS, c_sqrt_doc}, |
| 1077 | {"tan", cmath_tan, METH_VARARGS, c_tan_doc}, |
| 1078 | {"tanh", cmath_tanh, METH_VARARGS, c_tanh_doc}, |
| 1079 | {NULL, NULL} /* sentinel */ |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 1080 | }; |
| 1081 | |
Martin v. Löwis | 1a21451 | 2008-06-11 05:26:20 +0000 | [diff] [blame] | 1082 | |
| 1083 | static struct PyModuleDef cmathmodule = { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1084 | PyModuleDef_HEAD_INIT, |
| 1085 | "cmath", |
| 1086 | module_doc, |
| 1087 | -1, |
| 1088 | cmath_methods, |
| 1089 | NULL, |
| 1090 | NULL, |
| 1091 | NULL, |
| 1092 | NULL |
Martin v. Löwis | 1a21451 | 2008-06-11 05:26:20 +0000 | [diff] [blame] | 1093 | }; |
| 1094 | |
Mark Hammond | fe51c6d | 2002-08-02 02:27:13 +0000 | [diff] [blame] | 1095 | PyMODINIT_FUNC |
Martin v. Löwis | 1a21451 | 2008-06-11 05:26:20 +0000 | [diff] [blame] | 1096 | PyInit_cmath(void) |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 1097 | { |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1098 | PyObject *m; |
Tim Peters | 14e2640 | 2001-02-20 20:15:19 +0000 | [diff] [blame] | 1099 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1100 | m = PyModule_Create(&cmathmodule); |
| 1101 | if (m == NULL) |
| 1102 | return NULL; |
Fred Drake | f4e3484 | 2002-04-01 03:45:06 +0000 | [diff] [blame] | 1103 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1104 | PyModule_AddObject(m, "pi", |
| 1105 | PyFloat_FromDouble(Py_MATH_PI)); |
| 1106 | PyModule_AddObject(m, "e", PyFloat_FromDouble(Py_MATH_E)); |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 1107 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1108 | /* initialize special value tables */ |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 1109 | |
| 1110 | #define INIT_SPECIAL_VALUES(NAME, BODY) { Py_complex* p = (Py_complex*)NAME; BODY } |
| 1111 | #define C(REAL, IMAG) p->real = REAL; p->imag = IMAG; ++p; |
| 1112 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1113 | INIT_SPECIAL_VALUES(acos_special_values, { |
| 1114 | C(P34,INF) C(P,INF) C(P,INF) C(P,-INF) C(P,-INF) C(P34,-INF) C(N,INF) |
| 1115 | C(P12,INF) C(U,U) C(U,U) C(U,U) C(U,U) C(P12,-INF) C(N,N) |
| 1116 | C(P12,INF) C(U,U) C(P12,0.) C(P12,-0.) C(U,U) C(P12,-INF) C(P12,N) |
| 1117 | C(P12,INF) C(U,U) C(P12,0.) C(P12,-0.) C(U,U) C(P12,-INF) C(P12,N) |
| 1118 | C(P12,INF) C(U,U) C(U,U) C(U,U) C(U,U) C(P12,-INF) C(N,N) |
| 1119 | C(P14,INF) C(0.,INF) C(0.,INF) C(0.,-INF) C(0.,-INF) C(P14,-INF) C(N,INF) |
| 1120 | C(N,INF) C(N,N) C(N,N) C(N,N) C(N,N) C(N,-INF) C(N,N) |
| 1121 | }) |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 1122 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1123 | INIT_SPECIAL_VALUES(acosh_special_values, { |
| 1124 | C(INF,-P34) C(INF,-P) C(INF,-P) C(INF,P) C(INF,P) C(INF,P34) C(INF,N) |
| 1125 | C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) |
| 1126 | C(INF,-P12) C(U,U) C(0.,-P12) C(0.,P12) C(U,U) C(INF,P12) C(N,N) |
| 1127 | C(INF,-P12) C(U,U) C(0.,-P12) C(0.,P12) C(U,U) C(INF,P12) C(N,N) |
| 1128 | C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) |
| 1129 | C(INF,-P14) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,P14) C(INF,N) |
| 1130 | C(INF,N) C(N,N) C(N,N) C(N,N) C(N,N) C(INF,N) C(N,N) |
| 1131 | }) |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 1132 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1133 | INIT_SPECIAL_VALUES(asinh_special_values, { |
| 1134 | C(-INF,-P14) C(-INF,-0.) C(-INF,-0.) C(-INF,0.) C(-INF,0.) C(-INF,P14) C(-INF,N) |
| 1135 | C(-INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(-INF,P12) C(N,N) |
| 1136 | C(-INF,-P12) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(-INF,P12) C(N,N) |
| 1137 | C(INF,-P12) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(INF,P12) C(N,N) |
| 1138 | C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) |
| 1139 | C(INF,-P14) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,P14) C(INF,N) |
| 1140 | C(INF,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(INF,N) C(N,N) |
| 1141 | }) |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 1142 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1143 | INIT_SPECIAL_VALUES(atanh_special_values, { |
| 1144 | C(-0.,-P12) C(-0.,-P12) C(-0.,-P12) C(-0.,P12) C(-0.,P12) C(-0.,P12) C(-0.,N) |
| 1145 | C(-0.,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(-0.,P12) C(N,N) |
| 1146 | C(-0.,-P12) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(-0.,P12) C(-0.,N) |
| 1147 | C(0.,-P12) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,P12) C(0.,N) |
| 1148 | C(0.,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(0.,P12) C(N,N) |
| 1149 | C(0.,-P12) C(0.,-P12) C(0.,-P12) C(0.,P12) C(0.,P12) C(0.,P12) C(0.,N) |
| 1150 | C(0.,-P12) C(N,N) C(N,N) C(N,N) C(N,N) C(0.,P12) C(N,N) |
| 1151 | }) |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 1152 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1153 | INIT_SPECIAL_VALUES(cosh_special_values, { |
| 1154 | C(INF,N) C(U,U) C(INF,0.) C(INF,-0.) C(U,U) C(INF,N) C(INF,N) |
| 1155 | C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) |
| 1156 | C(N,0.) C(U,U) C(1.,0.) C(1.,-0.) C(U,U) C(N,0.) C(N,0.) |
| 1157 | C(N,0.) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(N,0.) C(N,0.) |
| 1158 | C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) |
| 1159 | C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N) |
| 1160 | C(N,N) C(N,N) C(N,0.) C(N,0.) C(N,N) C(N,N) C(N,N) |
| 1161 | }) |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 1162 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1163 | INIT_SPECIAL_VALUES(exp_special_values, { |
| 1164 | C(0.,0.) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,0.) C(0.,0.) |
| 1165 | C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) |
| 1166 | C(N,N) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(N,N) C(N,N) |
| 1167 | C(N,N) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(N,N) C(N,N) |
| 1168 | C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) |
| 1169 | C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N) |
| 1170 | C(N,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(N,N) C(N,N) |
| 1171 | }) |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 1172 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1173 | INIT_SPECIAL_VALUES(log_special_values, { |
| 1174 | C(INF,-P34) C(INF,-P) C(INF,-P) C(INF,P) C(INF,P) C(INF,P34) C(INF,N) |
| 1175 | C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) |
| 1176 | C(INF,-P12) C(U,U) C(-INF,-P) C(-INF,P) C(U,U) C(INF,P12) C(N,N) |
| 1177 | C(INF,-P12) C(U,U) C(-INF,-0.) C(-INF,0.) C(U,U) C(INF,P12) C(N,N) |
| 1178 | C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) |
| 1179 | C(INF,-P14) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,P14) C(INF,N) |
| 1180 | C(INF,N) C(N,N) C(N,N) C(N,N) C(N,N) C(INF,N) C(N,N) |
| 1181 | }) |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 1182 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1183 | INIT_SPECIAL_VALUES(sinh_special_values, { |
| 1184 | C(INF,N) C(U,U) C(-INF,-0.) C(-INF,0.) C(U,U) C(INF,N) C(INF,N) |
| 1185 | C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) |
| 1186 | C(0.,N) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(0.,N) C(0.,N) |
| 1187 | C(0.,N) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,N) C(0.,N) |
| 1188 | C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) |
| 1189 | C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N) |
| 1190 | C(N,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(N,N) C(N,N) |
| 1191 | }) |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 1192 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1193 | INIT_SPECIAL_VALUES(sqrt_special_values, { |
| 1194 | C(INF,-INF) C(0.,-INF) C(0.,-INF) C(0.,INF) C(0.,INF) C(INF,INF) C(N,INF) |
| 1195 | C(INF,-INF) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,INF) C(N,N) |
| 1196 | C(INF,-INF) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(INF,INF) C(N,N) |
| 1197 | C(INF,-INF) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(INF,INF) C(N,N) |
| 1198 | C(INF,-INF) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,INF) C(N,N) |
| 1199 | C(INF,-INF) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,INF) C(INF,N) |
| 1200 | C(INF,-INF) C(N,N) C(N,N) C(N,N) C(N,N) C(INF,INF) C(N,N) |
| 1201 | }) |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 1202 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1203 | INIT_SPECIAL_VALUES(tanh_special_values, { |
| 1204 | C(-1.,0.) C(U,U) C(-1.,-0.) C(-1.,0.) C(U,U) C(-1.,0.) C(-1.,0.) |
| 1205 | C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) |
| 1206 | C(N,N) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(N,N) C(N,N) |
| 1207 | C(N,N) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(N,N) C(N,N) |
| 1208 | C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) |
| 1209 | C(1.,0.) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(1.,0.) C(1.,0.) |
| 1210 | C(N,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(N,N) C(N,N) |
| 1211 | }) |
Christian Heimes | a342c01 | 2008-04-20 21:01:16 +0000 | [diff] [blame] | 1212 | |
Antoine Pitrou | 7f14f0d | 2010-05-09 16:14:21 +0000 | [diff] [blame] | 1213 | INIT_SPECIAL_VALUES(rect_special_values, { |
| 1214 | C(INF,N) C(U,U) C(-INF,0.) C(-INF,-0.) C(U,U) C(INF,N) C(INF,N) |
| 1215 | C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) |
| 1216 | C(0.,0.) C(U,U) C(-0.,0.) C(-0.,-0.) C(U,U) C(0.,0.) C(0.,0.) |
| 1217 | C(0.,0.) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,0.) C(0.,0.) |
| 1218 | C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) |
| 1219 | C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N) |
| 1220 | C(N,N) C(N,N) C(N,0.) C(N,0.) C(N,N) C(N,N) C(N,N) |
| 1221 | }) |
| 1222 | return m; |
Guido van Rossum | 71aa32f | 1996-01-12 01:34:57 +0000 | [diff] [blame] | 1223 | } |