diff options
author | Antoine Pitrou <solipsis@pitrou.net> | 2010-05-09 16:14:21 (GMT) |
---|---|---|
committer | Antoine Pitrou <solipsis@pitrou.net> | 2010-05-09 16:14:21 (GMT) |
commit | 7f14f0d8a0228c50d5b5de2acbfe9a64ebc6749a (patch) | |
tree | d25489e9531c01f1e9244012bbfaa929f382883e /Modules/cmathmodule.c | |
parent | b7d943625cf4353f6cb72df16252759f2dbd8e06 (diff) | |
download | cpython-7f14f0d8a0228c50d5b5de2acbfe9a64ebc6749a.zip cpython-7f14f0d8a0228c50d5b5de2acbfe9a64ebc6749a.tar.gz cpython-7f14f0d8a0228c50d5b5de2acbfe9a64ebc6749a.tar.bz2 |
Recorded merge of revisions 81032 via svnmerge from
svn+ssh://pythondev@svn.python.org/python/branches/py3k
................
r81032 | antoine.pitrou | 2010-05-09 17:52:27 +0200 (dim., 09 mai 2010) | 9 lines
Recorded merge of revisions 81029 via svnmerge from
svn+ssh://pythondev@svn.python.org/python/trunk
........
r81029 | antoine.pitrou | 2010-05-09 16:46:46 +0200 (dim., 09 mai 2010) | 3 lines
Untabify C files. Will watch buildbots.
........
................
Diffstat (limited to 'Modules/cmathmodule.c')
-rw-r--r-- | Modules/cmathmodule.c | 1620 |
1 files changed, 810 insertions, 810 deletions
diff --git a/Modules/cmathmodule.c b/Modules/cmathmodule.c index fbf6ece..545e834 100644 --- a/Modules/cmathmodule.c +++ b/Modules/cmathmodule.c @@ -31,7 +31,7 @@ #define CM_LOG_LARGE_DOUBLE (log(CM_LARGE_DOUBLE)) #define CM_SQRT_DBL_MIN (sqrt(DBL_MIN)) -/* +/* CM_SCALE_UP is an odd integer chosen such that multiplication by 2**CM_SCALE_UP is sufficient to turn a subnormal into a normal. CM_SCALE_DOWN is (-(CM_SCALE_UP+1)/2). These scalings are used to compute @@ -62,46 +62,46 @@ static PyObject * math_error(void); */ enum special_types { - ST_NINF, /* 0, negative infinity */ - ST_NEG, /* 1, negative finite number (nonzero) */ - ST_NZERO, /* 2, -0. */ - ST_PZERO, /* 3, +0. */ - ST_POS, /* 4, positive finite number (nonzero) */ - ST_PINF, /* 5, positive infinity */ - ST_NAN /* 6, Not a Number */ + ST_NINF, /* 0, negative infinity */ + ST_NEG, /* 1, negative finite number (nonzero) */ + ST_NZERO, /* 2, -0. */ + ST_PZERO, /* 3, +0. */ + ST_POS, /* 4, positive finite number (nonzero) */ + ST_PINF, /* 5, positive infinity */ + ST_NAN /* 6, Not a Number */ }; static enum special_types special_type(double d) { - if (Py_IS_FINITE(d)) { - if (d != 0) { - if (copysign(1., d) == 1.) - return ST_POS; - else - return ST_NEG; - } - else { - if (copysign(1., d) == 1.) - return ST_PZERO; - else - return ST_NZERO; - } - } - if (Py_IS_NAN(d)) - return ST_NAN; - if (copysign(1., d) == 1.) - return ST_PINF; - else - return ST_NINF; + if (Py_IS_FINITE(d)) { + if (d != 0) { + if (copysign(1., d) == 1.) + return ST_POS; + else + return ST_NEG; + } + else { + if (copysign(1., d) == 1.) + return ST_PZERO; + else + return ST_NZERO; + } + } + if (Py_IS_NAN(d)) + return ST_NAN; + if (copysign(1., d) == 1.) + return ST_PINF; + else + return ST_NINF; } -#define SPECIAL_VALUE(z, table) \ - if (!Py_IS_FINITE((z).real) || !Py_IS_FINITE((z).imag)) { \ - errno = 0; \ - return table[special_type((z).real)] \ - [special_type((z).imag)]; \ - } +#define SPECIAL_VALUE(z, table) \ + if (!Py_IS_FINITE((z).real) || !Py_IS_FINITE((z).imag)) { \ + errno = 0; \ + return table[special_type((z).real)] \ + [special_type((z).imag)]; \ + } #define P Py_MATH_PI #define P14 0.25*Py_MATH_PI @@ -125,34 +125,34 @@ static Py_complex acos_special_values[7][7]; static Py_complex c_acos(Py_complex z) { - Py_complex s1, s2, r; - - SPECIAL_VALUE(z, acos_special_values); - - if (fabs(z.real) > CM_LARGE_DOUBLE || fabs(z.imag) > CM_LARGE_DOUBLE) { - /* avoid unnecessary overflow for large arguments */ - r.real = atan2(fabs(z.imag), z.real); - /* split into cases to make sure that the branch cut has the - correct continuity on systems with unsigned zeros */ - if (z.real < 0.) { - r.imag = -copysign(log(hypot(z.real/2., z.imag/2.)) + - M_LN2*2., z.imag); - } else { - r.imag = copysign(log(hypot(z.real/2., z.imag/2.)) + - M_LN2*2., -z.imag); - } - } else { - s1.real = 1.-z.real; - s1.imag = -z.imag; - s1 = c_sqrt(s1); - s2.real = 1.+z.real; - s2.imag = z.imag; - s2 = c_sqrt(s2); - r.real = 2.*atan2(s1.real, s2.real); - r.imag = asinh(s2.real*s1.imag - s2.imag*s1.real); - } - errno = 0; - return r; + Py_complex s1, s2, r; + + SPECIAL_VALUE(z, acos_special_values); + + if (fabs(z.real) > CM_LARGE_DOUBLE || fabs(z.imag) > CM_LARGE_DOUBLE) { + /* avoid unnecessary overflow for large arguments */ + r.real = atan2(fabs(z.imag), z.real); + /* split into cases to make sure that the branch cut has the + correct continuity on systems with unsigned zeros */ + if (z.real < 0.) { + r.imag = -copysign(log(hypot(z.real/2., z.imag/2.)) + + M_LN2*2., z.imag); + } else { + r.imag = copysign(log(hypot(z.real/2., z.imag/2.)) + + M_LN2*2., -z.imag); + } + } else { + s1.real = 1.-z.real; + s1.imag = -z.imag; + s1 = c_sqrt(s1); + s2.real = 1.+z.real; + s2.imag = z.imag; + s2 = c_sqrt(s2); + r.real = 2.*atan2(s1.real, s2.real); + r.imag = asinh(s2.real*s1.imag - s2.imag*s1.real); + } + errno = 0; + return r; } PyDoc_STRVAR(c_acos_doc, @@ -166,26 +166,26 @@ static Py_complex acosh_special_values[7][7]; static Py_complex c_acosh(Py_complex z) { - Py_complex s1, s2, r; - - SPECIAL_VALUE(z, acosh_special_values); - - if (fabs(z.real) > CM_LARGE_DOUBLE || fabs(z.imag) > CM_LARGE_DOUBLE) { - /* avoid unnecessary overflow for large arguments */ - r.real = log(hypot(z.real/2., z.imag/2.)) + M_LN2*2.; - r.imag = atan2(z.imag, z.real); - } else { - s1.real = z.real - 1.; - s1.imag = z.imag; - s1 = c_sqrt(s1); - s2.real = z.real + 1.; - s2.imag = z.imag; - s2 = c_sqrt(s2); - r.real = asinh(s1.real*s2.real + s1.imag*s2.imag); - r.imag = 2.*atan2(s1.imag, s2.real); - } - errno = 0; - return r; + Py_complex s1, s2, r; + + SPECIAL_VALUE(z, acosh_special_values); + + if (fabs(z.real) > CM_LARGE_DOUBLE || fabs(z.imag) > CM_LARGE_DOUBLE) { + /* avoid unnecessary overflow for large arguments */ + r.real = log(hypot(z.real/2., z.imag/2.)) + M_LN2*2.; + r.imag = atan2(z.imag, z.real); + } else { + s1.real = z.real - 1.; + s1.imag = z.imag; + s1 = c_sqrt(s1); + s2.real = z.real + 1.; + s2.imag = z.imag; + s2 = c_sqrt(s2); + r.real = asinh(s1.real*s2.real + s1.imag*s2.imag); + r.imag = 2.*atan2(s1.imag, s2.real); + } + errno = 0; + return r; } PyDoc_STRVAR(c_acosh_doc, @@ -197,14 +197,14 @@ PyDoc_STRVAR(c_acosh_doc, static Py_complex c_asin(Py_complex z) { - /* asin(z) = -i asinh(iz) */ - Py_complex s, r; - s.real = -z.imag; - s.imag = z.real; - s = c_asinh(s); - r.real = s.imag; - r.imag = -s.real; - return r; + /* asin(z) = -i asinh(iz) */ + Py_complex s, r; + s.real = -z.imag; + s.imag = z.real; + s = c_asinh(s); + r.real = s.imag; + r.imag = -s.real; + return r; } PyDoc_STRVAR(c_asin_doc, @@ -218,31 +218,31 @@ static Py_complex asinh_special_values[7][7]; static Py_complex c_asinh(Py_complex z) { - Py_complex s1, s2, r; - - SPECIAL_VALUE(z, asinh_special_values); - - if (fabs(z.real) > CM_LARGE_DOUBLE || fabs(z.imag) > CM_LARGE_DOUBLE) { - if (z.imag >= 0.) { - r.real = copysign(log(hypot(z.real/2., z.imag/2.)) + - M_LN2*2., z.real); - } else { - r.real = -copysign(log(hypot(z.real/2., z.imag/2.)) + - M_LN2*2., -z.real); - } - r.imag = atan2(z.imag, fabs(z.real)); - } else { - s1.real = 1.+z.imag; - s1.imag = -z.real; - s1 = c_sqrt(s1); - s2.real = 1.-z.imag; - s2.imag = z.real; - s2 = c_sqrt(s2); - r.real = asinh(s1.real*s2.imag-s2.real*s1.imag); - r.imag = atan2(z.imag, s1.real*s2.real-s1.imag*s2.imag); - } - errno = 0; - return r; + Py_complex s1, s2, r; + + SPECIAL_VALUE(z, asinh_special_values); + + if (fabs(z.real) > CM_LARGE_DOUBLE || fabs(z.imag) > CM_LARGE_DOUBLE) { + if (z.imag >= 0.) { + r.real = copysign(log(hypot(z.real/2., z.imag/2.)) + + M_LN2*2., z.real); + } else { + r.real = -copysign(log(hypot(z.real/2., z.imag/2.)) + + M_LN2*2., -z.real); + } + r.imag = atan2(z.imag, fabs(z.real)); + } else { + s1.real = 1.+z.imag; + s1.imag = -z.real; + s1 = c_sqrt(s1); + s2.real = 1.-z.imag; + s2.imag = z.real; + s2 = c_sqrt(s2); + r.real = asinh(s1.real*s2.imag-s2.real*s1.imag); + r.imag = atan2(z.imag, s1.real*s2.real-s1.imag*s2.imag); + } + errno = 0; + return r; } PyDoc_STRVAR(c_asinh_doc, @@ -254,14 +254,14 @@ PyDoc_STRVAR(c_asinh_doc, static Py_complex c_atan(Py_complex z) { - /* atan(z) = -i atanh(iz) */ - Py_complex s, r; - s.real = -z.imag; - s.imag = z.real; - s = c_atanh(s); - r.real = s.imag; - r.imag = -s.real; - return r; + /* atan(z) = -i atanh(iz) */ + Py_complex s, r; + s.real = -z.imag; + s.imag = z.real; + s = c_atanh(s); + r.real = s.imag; + r.imag = -s.real; + return r; } /* Windows screws up atan2 for inf and nan, and alpha Tru64 5.1 doesn't follow @@ -269,29 +269,29 @@ c_atan(Py_complex z) static double c_atan2(Py_complex z) { - if (Py_IS_NAN(z.real) || Py_IS_NAN(z.imag)) - return Py_NAN; - if (Py_IS_INFINITY(z.imag)) { - if (Py_IS_INFINITY(z.real)) { - if (copysign(1., z.real) == 1.) - /* atan2(+-inf, +inf) == +-pi/4 */ - return copysign(0.25*Py_MATH_PI, z.imag); - else - /* atan2(+-inf, -inf) == +-pi*3/4 */ - return copysign(0.75*Py_MATH_PI, z.imag); - } - /* atan2(+-inf, x) == +-pi/2 for finite x */ - return copysign(0.5*Py_MATH_PI, z.imag); - } - if (Py_IS_INFINITY(z.real) || z.imag == 0.) { - if (copysign(1., z.real) == 1.) - /* atan2(+-y, +inf) = atan2(+-0, +x) = +-0. */ - return copysign(0., z.imag); - else - /* atan2(+-y, -inf) = atan2(+-0., -x) = +-pi. */ - return copysign(Py_MATH_PI, z.imag); - } - return atan2(z.imag, z.real); + if (Py_IS_NAN(z.real) || Py_IS_NAN(z.imag)) + return Py_NAN; + if (Py_IS_INFINITY(z.imag)) { + if (Py_IS_INFINITY(z.real)) { + if (copysign(1., z.real) == 1.) + /* atan2(+-inf, +inf) == +-pi/4 */ + return copysign(0.25*Py_MATH_PI, z.imag); + else + /* atan2(+-inf, -inf) == +-pi*3/4 */ + return copysign(0.75*Py_MATH_PI, z.imag); + } + /* atan2(+-inf, x) == +-pi/2 for finite x */ + return copysign(0.5*Py_MATH_PI, z.imag); + } + if (Py_IS_INFINITY(z.real) || z.imag == 0.) { + if (copysign(1., z.real) == 1.) + /* atan2(+-y, +inf) = atan2(+-0, +x) = +-0. */ + return copysign(0., z.imag); + else + /* atan2(+-y, -inf) = atan2(+-0., -x) = +-pi. */ + return copysign(Py_MATH_PI, z.imag); + } + return atan2(z.imag, z.real); } PyDoc_STRVAR(c_atan_doc, @@ -305,48 +305,48 @@ static Py_complex atanh_special_values[7][7]; static Py_complex c_atanh(Py_complex z) { - Py_complex r; - double ay, h; - - SPECIAL_VALUE(z, atanh_special_values); - - /* Reduce to case where z.real >= 0., using atanh(z) = -atanh(-z). */ - if (z.real < 0.) { - return c_neg(c_atanh(c_neg(z))); - } - - ay = fabs(z.imag); - if (z.real > CM_SQRT_LARGE_DOUBLE || ay > CM_SQRT_LARGE_DOUBLE) { - /* - if abs(z) is large then we use the approximation - atanh(z) ~ 1/z +/- i*pi/2 (+/- depending on the sign - of z.imag) - */ - h = hypot(z.real/2., z.imag/2.); /* safe from overflow */ - r.real = z.real/4./h/h; - /* the two negations in the next line cancel each other out - except when working with unsigned zeros: they're there to - ensure that the branch cut has the correct continuity on - systems that don't support signed zeros */ - r.imag = -copysign(Py_MATH_PI/2., -z.imag); - errno = 0; - } else if (z.real == 1. && ay < CM_SQRT_DBL_MIN) { - /* C99 standard says: atanh(1+/-0.) should be inf +/- 0i */ - if (ay == 0.) { - r.real = INF; - r.imag = z.imag; - errno = EDOM; - } else { - r.real = -log(sqrt(ay)/sqrt(hypot(ay, 2.))); - r.imag = copysign(atan2(2., -ay)/2, z.imag); - errno = 0; - } - } else { - r.real = log1p(4.*z.real/((1-z.real)*(1-z.real) + ay*ay))/4.; - r.imag = -atan2(-2.*z.imag, (1-z.real)*(1+z.real) - ay*ay)/2.; - errno = 0; - } - return r; + Py_complex r; + double ay, h; + + SPECIAL_VALUE(z, atanh_special_values); + + /* Reduce to case where z.real >= 0., using atanh(z) = -atanh(-z). */ + if (z.real < 0.) { + return c_neg(c_atanh(c_neg(z))); + } + + ay = fabs(z.imag); + if (z.real > CM_SQRT_LARGE_DOUBLE || ay > CM_SQRT_LARGE_DOUBLE) { + /* + if abs(z) is large then we use the approximation + atanh(z) ~ 1/z +/- i*pi/2 (+/- depending on the sign + of z.imag) + */ + h = hypot(z.real/2., z.imag/2.); /* safe from overflow */ + r.real = z.real/4./h/h; + /* the two negations in the next line cancel each other out + except when working with unsigned zeros: they're there to + ensure that the branch cut has the correct continuity on + systems that don't support signed zeros */ + r.imag = -copysign(Py_MATH_PI/2., -z.imag); + errno = 0; + } else if (z.real == 1. && ay < CM_SQRT_DBL_MIN) { + /* C99 standard says: atanh(1+/-0.) should be inf +/- 0i */ + if (ay == 0.) { + r.real = INF; + r.imag = z.imag; + errno = EDOM; + } else { + r.real = -log(sqrt(ay)/sqrt(hypot(ay, 2.))); + r.imag = copysign(atan2(2., -ay)/2, z.imag); + errno = 0; + } + } else { + r.real = log1p(4.*z.real/((1-z.real)*(1-z.real) + ay*ay))/4.; + r.imag = -atan2(-2.*z.imag, (1-z.real)*(1+z.real) - ay*ay)/2.; + errno = 0; + } + return r; } PyDoc_STRVAR(c_atanh_doc, @@ -358,12 +358,12 @@ PyDoc_STRVAR(c_atanh_doc, static Py_complex c_cos(Py_complex z) { - /* cos(z) = cosh(iz) */ - Py_complex r; - r.real = -z.imag; - r.imag = z.real; - r = c_cosh(r); - return r; + /* cos(z) = cosh(iz) */ + Py_complex r; + r.real = -z.imag; + r.imag = z.real; + r = c_cosh(r); + return r; } PyDoc_STRVAR(c_cos_doc, @@ -378,51 +378,51 @@ static Py_complex cosh_special_values[7][7]; static Py_complex c_cosh(Py_complex z) { - Py_complex r; - double x_minus_one; - - /* special treatment for cosh(+/-inf + iy) if y is not a NaN */ - if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) { - if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) && - (z.imag != 0.)) { - if (z.real > 0) { - r.real = copysign(INF, cos(z.imag)); - r.imag = copysign(INF, sin(z.imag)); - } - else { - r.real = copysign(INF, cos(z.imag)); - r.imag = -copysign(INF, sin(z.imag)); - } - } - else { - r = cosh_special_values[special_type(z.real)] - [special_type(z.imag)]; - } - /* need to set errno = EDOM if y is +/- infinity and x is not - a NaN */ - if (Py_IS_INFINITY(z.imag) && !Py_IS_NAN(z.real)) - errno = EDOM; - else - errno = 0; - return r; - } - - if (fabs(z.real) > CM_LOG_LARGE_DOUBLE) { - /* deal correctly with cases where cosh(z.real) overflows but - cosh(z) does not. */ - x_minus_one = z.real - copysign(1., z.real); - r.real = cos(z.imag) * cosh(x_minus_one) * Py_MATH_E; - r.imag = sin(z.imag) * sinh(x_minus_one) * Py_MATH_E; - } else { - r.real = cos(z.imag) * cosh(z.real); - r.imag = sin(z.imag) * sinh(z.real); - } - /* detect overflow, and set errno accordingly */ - if (Py_IS_INFINITY(r.real) || Py_IS_INFINITY(r.imag)) - errno = ERANGE; - else - errno = 0; - return r; + Py_complex r; + double x_minus_one; + + /* special treatment for cosh(+/-inf + iy) if y is not a NaN */ + if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) { + if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) && + (z.imag != 0.)) { + if (z.real > 0) { + r.real = copysign(INF, cos(z.imag)); + r.imag = copysign(INF, sin(z.imag)); + } + else { + r.real = copysign(INF, cos(z.imag)); + r.imag = -copysign(INF, sin(z.imag)); + } + } + else { + r = cosh_special_values[special_type(z.real)] + [special_type(z.imag)]; + } + /* need to set errno = EDOM if y is +/- infinity and x is not + a NaN */ + if (Py_IS_INFINITY(z.imag) && !Py_IS_NAN(z.real)) + errno = EDOM; + else + errno = 0; + return r; + } + + if (fabs(z.real) > CM_LOG_LARGE_DOUBLE) { + /* deal correctly with cases where cosh(z.real) overflows but + cosh(z) does not. */ + x_minus_one = z.real - copysign(1., z.real); + r.real = cos(z.imag) * cosh(x_minus_one) * Py_MATH_E; + r.imag = sin(z.imag) * sinh(x_minus_one) * Py_MATH_E; + } else { + r.real = cos(z.imag) * cosh(z.real); + r.imag = sin(z.imag) * sinh(z.real); + } + /* detect overflow, and set errno accordingly */ + if (Py_IS_INFINITY(r.real) || Py_IS_INFINITY(r.imag)) + errno = ERANGE; + else + errno = 0; + return r; } PyDoc_STRVAR(c_cosh_doc, @@ -438,51 +438,51 @@ static Py_complex exp_special_values[7][7]; static Py_complex c_exp(Py_complex z) { - Py_complex r; - double l; - - if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) { - if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) - && (z.imag != 0.)) { - if (z.real > 0) { - r.real = copysign(INF, cos(z.imag)); - r.imag = copysign(INF, sin(z.imag)); - } - else { - r.real = copysign(0., cos(z.imag)); - r.imag = copysign(0., sin(z.imag)); - } - } - else { - r = exp_special_values[special_type(z.real)] - [special_type(z.imag)]; - } - /* need to set errno = EDOM if y is +/- infinity and x is not - a NaN and not -infinity */ - if (Py_IS_INFINITY(z.imag) && - (Py_IS_FINITE(z.real) || - (Py_IS_INFINITY(z.real) && z.real > 0))) - errno = EDOM; - else - errno = 0; - return r; - } - - if (z.real > CM_LOG_LARGE_DOUBLE) { - l = exp(z.real-1.); - r.real = l*cos(z.imag)*Py_MATH_E; - r.imag = l*sin(z.imag)*Py_MATH_E; - } else { - l = exp(z.real); - r.real = l*cos(z.imag); - r.imag = l*sin(z.imag); - } - /* detect overflow, and set errno accordingly */ - if (Py_IS_INFINITY(r.real) || Py_IS_INFINITY(r.imag)) - errno = ERANGE; - else - errno = 0; - return r; + Py_complex r; + double l; + + if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) { + if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) + && (z.imag != 0.)) { + if (z.real > 0) { + r.real = copysign(INF, cos(z.imag)); + r.imag = copysign(INF, sin(z.imag)); + } + else { + r.real = copysign(0., cos(z.imag)); + r.imag = copysign(0., sin(z.imag)); + } + } + else { + r = exp_special_values[special_type(z.real)] + [special_type(z.imag)]; + } + /* need to set errno = EDOM if y is +/- infinity and x is not + a NaN and not -infinity */ + if (Py_IS_INFINITY(z.imag) && + (Py_IS_FINITE(z.real) || + (Py_IS_INFINITY(z.real) && z.real > 0))) + errno = EDOM; + else + errno = 0; + return r; + } + + if (z.real > CM_LOG_LARGE_DOUBLE) { + l = exp(z.real-1.); + r.real = l*cos(z.imag)*Py_MATH_E; + r.imag = l*sin(z.imag)*Py_MATH_E; + } else { + l = exp(z.real); + r.real = l*cos(z.imag); + r.imag = l*sin(z.imag); + } + /* detect overflow, and set errno accordingly */ + if (Py_IS_INFINITY(r.real) || Py_IS_INFINITY(r.imag)) + errno = ERANGE; + else + errno = 0; + return r; } PyDoc_STRVAR(c_exp_doc, @@ -496,85 +496,85 @@ static Py_complex log_special_values[7][7]; static Py_complex c_log(Py_complex z) { - /* - The usual formula for the real part is log(hypot(z.real, z.imag)). - There are four situations where this formula is potentially - problematic: - - (1) the absolute value of z is subnormal. Then hypot is subnormal, - so has fewer than the usual number of bits of accuracy, hence may - have large relative error. This then gives a large absolute error - in the log. This can be solved by rescaling z by a suitable power - of 2. - - (2) the absolute value of z is greater than DBL_MAX (e.g. when both - z.real and z.imag are within a factor of 1/sqrt(2) of DBL_MAX) - Again, rescaling solves this. - - (3) the absolute value of z is close to 1. In this case it's - difficult to achieve good accuracy, at least in part because a - change of 1ulp in the real or imaginary part of z can result in a - change of billions of ulps in the correctly rounded answer. - - (4) z = 0. The simplest thing to do here is to call the - floating-point log with an argument of 0, and let its behaviour - (returning -infinity, signaling a floating-point exception, setting - errno, or whatever) determine that of c_log. So the usual formula - is fine here. - - */ - - Py_complex r; - double ax, ay, am, an, h; - - SPECIAL_VALUE(z, log_special_values); - - ax = fabs(z.real); - ay = fabs(z.imag); - - if (ax > CM_LARGE_DOUBLE || ay > CM_LARGE_DOUBLE) { - r.real = log(hypot(ax/2., ay/2.)) + M_LN2; - } else if (ax < DBL_MIN && ay < DBL_MIN) { - if (ax > 0. || ay > 0.) { - /* catch cases where hypot(ax, ay) is subnormal */ - r.real = log(hypot(ldexp(ax, DBL_MANT_DIG), - ldexp(ay, DBL_MANT_DIG))) - DBL_MANT_DIG*M_LN2; - } - else { - /* log(+/-0. +/- 0i) */ - r.real = -INF; - r.imag = atan2(z.imag, z.real); - errno = EDOM; - return r; - } - } else { - h = hypot(ax, ay); - if (0.71 <= h && h <= 1.73) { - am = ax > ay ? ax : ay; /* max(ax, ay) */ - an = ax > ay ? ay : ax; /* min(ax, ay) */ - r.real = log1p((am-1)*(am+1)+an*an)/2.; - } else { - r.real = log(h); - } - } - r.imag = atan2(z.imag, z.real); - errno = 0; - return r; + /* + The usual formula for the real part is log(hypot(z.real, z.imag)). + There are four situations where this formula is potentially + problematic: + + (1) the absolute value of z is subnormal. Then hypot is subnormal, + so has fewer than the usual number of bits of accuracy, hence may + have large relative error. This then gives a large absolute error + in the log. This can be solved by rescaling z by a suitable power + of 2. + + (2) the absolute value of z is greater than DBL_MAX (e.g. when both + z.real and z.imag are within a factor of 1/sqrt(2) of DBL_MAX) + Again, rescaling solves this. + + (3) the absolute value of z is close to 1. In this case it's + difficult to achieve good accuracy, at least in part because a + change of 1ulp in the real or imaginary part of z can result in a + change of billions of ulps in the correctly rounded answer. + + (4) z = 0. The simplest thing to do here is to call the + floating-point log with an argument of 0, and let its behaviour + (returning -infinity, signaling a floating-point exception, setting + errno, or whatever) determine that of c_log. So the usual formula + is fine here. + + */ + + Py_complex r; + double ax, ay, am, an, h; + + SPECIAL_VALUE(z, log_special_values); + + ax = fabs(z.real); + ay = fabs(z.imag); + + if (ax > CM_LARGE_DOUBLE || ay > CM_LARGE_DOUBLE) { + r.real = log(hypot(ax/2., ay/2.)) + M_LN2; + } else if (ax < DBL_MIN && ay < DBL_MIN) { + if (ax > 0. || ay > 0.) { + /* catch cases where hypot(ax, ay) is subnormal */ + r.real = log(hypot(ldexp(ax, DBL_MANT_DIG), + ldexp(ay, DBL_MANT_DIG))) - DBL_MANT_DIG*M_LN2; + } + else { + /* log(+/-0. +/- 0i) */ + r.real = -INF; + r.imag = atan2(z.imag, z.real); + errno = EDOM; + return r; + } + } else { + h = hypot(ax, ay); + if (0.71 <= h && h <= 1.73) { + am = ax > ay ? ax : ay; /* max(ax, ay) */ + an = ax > ay ? ay : ax; /* min(ax, ay) */ + r.real = log1p((am-1)*(am+1)+an*an)/2.; + } else { + r.real = log(h); + } + } + r.imag = atan2(z.imag, z.real); + errno = 0; + return r; } static Py_complex c_log10(Py_complex z) { - Py_complex r; - int errno_save; - - r = c_log(z); - errno_save = errno; /* just in case the divisions affect errno */ - r.real = r.real / M_LN10; - r.imag = r.imag / M_LN10; - errno = errno_save; - return r; + Py_complex r; + int errno_save; + + r = c_log(z); + errno_save = errno; /* just in case the divisions affect errno */ + r.real = r.real / M_LN10; + r.imag = r.imag / M_LN10; + errno = errno_save; + return r; } PyDoc_STRVAR(c_log10_doc, @@ -586,14 +586,14 @@ PyDoc_STRVAR(c_log10_doc, static Py_complex c_sin(Py_complex z) { - /* sin(z) = -i sin(iz) */ - Py_complex s, r; - s.real = -z.imag; - s.imag = z.real; - s = c_sinh(s); - r.real = s.imag; - r.imag = -s.real; - return r; + /* sin(z) = -i sin(iz) */ + Py_complex s, r; + s.real = -z.imag; + s.imag = z.real; + s = c_sinh(s); + r.real = s.imag; + r.imag = -s.real; + return r; } PyDoc_STRVAR(c_sin_doc, @@ -608,50 +608,50 @@ static Py_complex sinh_special_values[7][7]; static Py_complex c_sinh(Py_complex z) { - Py_complex r; - double x_minus_one; - - /* special treatment for sinh(+/-inf + iy) if y is finite and - nonzero */ - if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) { - if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) - && (z.imag != 0.)) { - if (z.real > 0) { - r.real = copysign(INF, cos(z.imag)); - r.imag = copysign(INF, sin(z.imag)); - } - else { - r.real = -copysign(INF, cos(z.imag)); - r.imag = copysign(INF, sin(z.imag)); - } - } - else { - r = sinh_special_values[special_type(z.real)] - [special_type(z.imag)]; - } - /* need to set errno = EDOM if y is +/- infinity and x is not - a NaN */ - if (Py_IS_INFINITY(z.imag) && !Py_IS_NAN(z.real)) - errno = EDOM; - else - errno = 0; - return r; - } - - if (fabs(z.real) > CM_LOG_LARGE_DOUBLE) { - x_minus_one = z.real - copysign(1., z.real); - r.real = cos(z.imag) * sinh(x_minus_one) * Py_MATH_E; - r.imag = sin(z.imag) * cosh(x_minus_one) * Py_MATH_E; - } else { - r.real = cos(z.imag) * sinh(z.real); - r.imag = sin(z.imag) * cosh(z.real); - } - /* detect overflow, and set errno accordingly */ - if (Py_IS_INFINITY(r.real) || Py_IS_INFINITY(r.imag)) - errno = ERANGE; - else - errno = 0; - return r; + Py_complex r; + double x_minus_one; + + /* special treatment for sinh(+/-inf + iy) if y is finite and + nonzero */ + if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) { + if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) + && (z.imag != 0.)) { + if (z.real > 0) { + r.real = copysign(INF, cos(z.imag)); + r.imag = copysign(INF, sin(z.imag)); + } + else { + r.real = -copysign(INF, cos(z.imag)); + r.imag = copysign(INF, sin(z.imag)); + } + } + else { + r = sinh_special_values[special_type(z.real)] + [special_type(z.imag)]; + } + /* need to set errno = EDOM if y is +/- infinity and x is not + a NaN */ + if (Py_IS_INFINITY(z.imag) && !Py_IS_NAN(z.real)) + errno = EDOM; + else + errno = 0; + return r; + } + + if (fabs(z.real) > CM_LOG_LARGE_DOUBLE) { + x_minus_one = z.real - copysign(1., z.real); + r.real = cos(z.imag) * sinh(x_minus_one) * Py_MATH_E; + r.imag = sin(z.imag) * cosh(x_minus_one) * Py_MATH_E; + } else { + r.real = cos(z.imag) * sinh(z.real); + r.imag = sin(z.imag) * cosh(z.real); + } + /* detect overflow, and set errno accordingly */ + if (Py_IS_INFINITY(r.real) || Py_IS_INFINITY(r.imag)) + errno = ERANGE; + else + errno = 0; + return r; } PyDoc_STRVAR(c_sinh_doc, @@ -665,68 +665,68 @@ static Py_complex sqrt_special_values[7][7]; static Py_complex c_sqrt(Py_complex z) { - /* - Method: use symmetries to reduce to the case when x = z.real and y - = z.imag are nonnegative. Then the real part of the result is - given by - - s = sqrt((x + hypot(x, y))/2) - - and the imaginary part is - - d = (y/2)/s - - If either x or y is very large then there's a risk of overflow in - computation of the expression x + hypot(x, y). We can avoid this - by rewriting the formula for s as: - - s = 2*sqrt(x/8 + hypot(x/8, y/8)) - - This costs us two extra multiplications/divisions, but avoids the - overhead of checking for x and y large. - - If both x and y are subnormal then hypot(x, y) may also be - subnormal, so will lack full precision. We solve this by rescaling - x and y by a sufficiently large power of 2 to ensure that x and y - are normal. - */ - - - Py_complex r; - double s,d; - double ax, ay; - - SPECIAL_VALUE(z, sqrt_special_values); - - if (z.real == 0. && z.imag == 0.) { - r.real = 0.; - r.imag = z.imag; - return r; - } - - ax = fabs(z.real); - ay = fabs(z.imag); - - if (ax < DBL_MIN && ay < DBL_MIN && (ax > 0. || ay > 0.)) { - /* here we catch cases where hypot(ax, ay) is subnormal */ - ax = ldexp(ax, CM_SCALE_UP); - s = ldexp(sqrt(ax + hypot(ax, ldexp(ay, CM_SCALE_UP))), - CM_SCALE_DOWN); - } else { - ax /= 8.; - s = 2.*sqrt(ax + hypot(ax, ay/8.)); - } - d = ay/(2.*s); - - if (z.real >= 0.) { - r.real = s; - r.imag = copysign(d, z.imag); - } else { - r.real = d; - r.imag = copysign(s, z.imag); - } - errno = 0; - return r; + /* + Method: use symmetries to reduce to the case when x = z.real and y + = z.imag are nonnegative. Then the real part of the result is + given by + + s = sqrt((x + hypot(x, y))/2) + + and the imaginary part is + + d = (y/2)/s + + If either x or y is very large then there's a risk of overflow in + computation of the expression x + hypot(x, y). We can avoid this + by rewriting the formula for s as: + + s = 2*sqrt(x/8 + hypot(x/8, y/8)) + + This costs us two extra multiplications/divisions, but avoids the + overhead of checking for x and y large. + + If both x and y are subnormal then hypot(x, y) may also be + subnormal, so will lack full precision. We solve this by rescaling + x and y by a sufficiently large power of 2 to ensure that x and y + are normal. + */ + + + Py_complex r; + double s,d; + double ax, ay; + + SPECIAL_VALUE(z, sqrt_special_values); + + if (z.real == 0. && z.imag == 0.) { + r.real = 0.; + r.imag = z.imag; + return r; + } + + ax = fabs(z.real); + ay = fabs(z.imag); + + if (ax < DBL_MIN && ay < DBL_MIN && (ax > 0. || ay > 0.)) { + /* here we catch cases where hypot(ax, ay) is subnormal */ + ax = ldexp(ax, CM_SCALE_UP); + s = ldexp(sqrt(ax + hypot(ax, ldexp(ay, CM_SCALE_UP))), + CM_SCALE_DOWN); + } else { + ax /= 8.; + s = 2.*sqrt(ax + hypot(ax, ay/8.)); + } + d = ay/(2.*s); + + if (z.real >= 0.) { + r.real = s; + r.imag = copysign(d, z.imag); + } else { + r.real = d; + r.imag = copysign(s, z.imag); + } + errno = 0; + return r; } PyDoc_STRVAR(c_sqrt_doc, @@ -738,14 +738,14 @@ PyDoc_STRVAR(c_sqrt_doc, static Py_complex c_tan(Py_complex z) { - /* tan(z) = -i tanh(iz) */ - Py_complex s, r; - s.real = -z.imag; - s.imag = z.real; - s = c_tanh(s); - r.real = s.imag; - r.imag = -s.real; - return r; + /* tan(z) = -i tanh(iz) */ + Py_complex s, r; + s.real = -z.imag; + s.imag = z.real; + s = c_tanh(s); + r.real = s.imag; + r.imag = -s.real; + return r; } PyDoc_STRVAR(c_tan_doc, @@ -760,65 +760,65 @@ static Py_complex tanh_special_values[7][7]; static Py_complex c_tanh(Py_complex z) { - /* Formula: - - tanh(x+iy) = (tanh(x)(1+tan(y)^2) + i tan(y)(1-tanh(x))^2) / - (1+tan(y)^2 tanh(x)^2) - - To avoid excessive roundoff error, 1-tanh(x)^2 is better computed - as 1/cosh(x)^2. When abs(x) is large, we approximate 1-tanh(x)^2 - by 4 exp(-2*x) instead, to avoid possible overflow in the - computation of cosh(x). - - */ - - Py_complex r; - double tx, ty, cx, txty, denom; - - /* special treatment for tanh(+/-inf + iy) if y is finite and - nonzero */ - if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) { - if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) - && (z.imag != 0.)) { - if (z.real > 0) { - r.real = 1.0; - r.imag = copysign(0., - 2.*sin(z.imag)*cos(z.imag)); - } - else { - r.real = -1.0; - r.imag = copysign(0., - 2.*sin(z.imag)*cos(z.imag)); - } - } - else { - r = tanh_special_values[special_type(z.real)] - [special_type(z.imag)]; - } - /* need to set errno = EDOM if z.imag is +/-infinity and - z.real is finite */ - if (Py_IS_INFINITY(z.imag) && Py_IS_FINITE(z.real)) - errno = EDOM; - else - errno = 0; - return r; - } - - /* danger of overflow in 2.*z.imag !*/ - if (fabs(z.real) > CM_LOG_LARGE_DOUBLE) { - r.real = copysign(1., z.real); - r.imag = 4.*sin(z.imag)*cos(z.imag)*exp(-2.*fabs(z.real)); - } else { - tx = tanh(z.real); - ty = tan(z.imag); - cx = 1./cosh(z.real); - txty = tx*ty; - denom = 1. + txty*txty; - r.real = tx*(1.+ty*ty)/denom; - r.imag = ((ty/denom)*cx)*cx; - } - errno = 0; - return r; + /* Formula: + + tanh(x+iy) = (tanh(x)(1+tan(y)^2) + i tan(y)(1-tanh(x))^2) / + (1+tan(y)^2 tanh(x)^2) + + To avoid excessive roundoff error, 1-tanh(x)^2 is better computed + as 1/cosh(x)^2. When abs(x) is large, we approximate 1-tanh(x)^2 + by 4 exp(-2*x) instead, to avoid possible overflow in the + computation of cosh(x). + + */ + + Py_complex r; + double tx, ty, cx, txty, denom; + + /* special treatment for tanh(+/-inf + iy) if y is finite and + nonzero */ + if (!Py_IS_FINITE(z.real) || !Py_IS_FINITE(z.imag)) { + if (Py_IS_INFINITY(z.real) && Py_IS_FINITE(z.imag) + && (z.imag != 0.)) { + if (z.real > 0) { + r.real = 1.0; + r.imag = copysign(0., + 2.*sin(z.imag)*cos(z.imag)); + } + else { + r.real = -1.0; + r.imag = copysign(0., + 2.*sin(z.imag)*cos(z.imag)); + } + } + else { + r = tanh_special_values[special_type(z.real)] + [special_type(z.imag)]; + } + /* need to set errno = EDOM if z.imag is +/-infinity and + z.real is finite */ + if (Py_IS_INFINITY(z.imag) && Py_IS_FINITE(z.real)) + errno = EDOM; + else + errno = 0; + return r; + } + + /* danger of overflow in 2.*z.imag !*/ + if (fabs(z.real) > CM_LOG_LARGE_DOUBLE) { + r.real = copysign(1., z.real); + r.imag = 4.*sin(z.imag)*cos(z.imag)*exp(-2.*fabs(z.real)); + } else { + tx = tanh(z.real); + ty = tan(z.imag); + cx = 1./cosh(z.real); + txty = tx*ty; + denom = 1. + txty*txty; + r.real = tx*(1.+ty*ty)/denom; + r.imag = ((ty/denom)*cx)*cx; + } + errno = 0; + return r; } PyDoc_STRVAR(c_tanh_doc, @@ -830,23 +830,23 @@ PyDoc_STRVAR(c_tanh_doc, static PyObject * cmath_log(PyObject *self, PyObject *args) { - Py_complex x; - Py_complex y; - - if (!PyArg_ParseTuple(args, "D|D", &x, &y)) - return NULL; - - errno = 0; - PyFPE_START_PROTECT("complex function", return 0) - x = c_log(x); - if (PyTuple_GET_SIZE(args) == 2) { - y = c_log(y); - x = c_quot(x, y); - } - PyFPE_END_PROTECT(x) - if (errno != 0) - return math_error(); - return PyComplex_FromCComplex(x); + Py_complex x; + Py_complex y; + + if (!PyArg_ParseTuple(args, "D|D", &x, &y)) + return NULL; + + errno = 0; + PyFPE_START_PROTECT("complex function", return 0) + x = c_log(x); + if (PyTuple_GET_SIZE(args) == 2) { + y = c_log(y); + x = c_quot(x, y); + } + PyFPE_END_PROTECT(x) + if (errno != 0) + return math_error(); + return PyComplex_FromCComplex(x); } PyDoc_STRVAR(cmath_log_doc, @@ -859,42 +859,42 @@ If the base not specified, returns the natural logarithm (base e) of x."); static PyObject * math_error(void) { - if (errno == EDOM) - PyErr_SetString(PyExc_ValueError, "math domain error"); - else if (errno == ERANGE) - PyErr_SetString(PyExc_OverflowError, "math range error"); - else /* Unexpected math error */ - PyErr_SetFromErrno(PyExc_ValueError); - return NULL; + if (errno == EDOM) + PyErr_SetString(PyExc_ValueError, "math domain error"); + else if (errno == ERANGE) + PyErr_SetString(PyExc_OverflowError, "math range error"); + else /* Unexpected math error */ + PyErr_SetFromErrno(PyExc_ValueError); + return NULL; } static PyObject * math_1(PyObject *args, Py_complex (*func)(Py_complex)) { - Py_complex x,r ; - if (!PyArg_ParseTuple(args, "D", &x)) - return NULL; - errno = 0; - PyFPE_START_PROTECT("complex function", return 0); - r = (*func)(x); - PyFPE_END_PROTECT(r); - if (errno == EDOM) { - PyErr_SetString(PyExc_ValueError, "math domain error"); - return NULL; - } - else if (errno == ERANGE) { - PyErr_SetString(PyExc_OverflowError, "math range error"); - return NULL; - } - else { - return PyComplex_FromCComplex(r); - } + Py_complex x,r ; + if (!PyArg_ParseTuple(args, "D", &x)) + return NULL; + errno = 0; + PyFPE_START_PROTECT("complex function", return 0); + r = (*func)(x); + PyFPE_END_PROTECT(r); + if (errno == EDOM) { + PyErr_SetString(PyExc_ValueError, "math domain error"); + return NULL; + } + else if (errno == ERANGE) { + PyErr_SetString(PyExc_OverflowError, "math range error"); + return NULL; + } + else { + return PyComplex_FromCComplex(r); + } } #define FUNC1(stubname, func) \ - static PyObject * stubname(PyObject *self, PyObject *args) { \ - return math_1(args, func); \ - } + static PyObject * stubname(PyObject *self, PyObject *args) { \ + return math_1(args, func); \ + } FUNC1(cmath_acos, c_acos) FUNC1(cmath_acosh, c_acosh) @@ -915,18 +915,18 @@ FUNC1(cmath_tanh, c_tanh) static PyObject * cmath_phase(PyObject *self, PyObject *args) { - Py_complex z; - double phi; - if (!PyArg_ParseTuple(args, "D:phase", &z)) - return NULL; - errno = 0; - PyFPE_START_PROTECT("arg function", return 0) - phi = c_atan2(z); - PyFPE_END_PROTECT(phi) - if (errno != 0) - return math_error(); - else - return PyFloat_FromDouble(phi); + Py_complex z; + double phi; + if (!PyArg_ParseTuple(args, "D:phase", &z)) + return NULL; + errno = 0; + PyFPE_START_PROTECT("arg function", return 0) + phi = c_atan2(z); + PyFPE_END_PROTECT(phi) + if (errno != 0) + return math_error(); + else + return PyFloat_FromDouble(phi); } PyDoc_STRVAR(cmath_phase_doc, @@ -936,18 +936,18 @@ Return argument, also known as the phase angle, of a complex."); static PyObject * cmath_polar(PyObject *self, PyObject *args) { - Py_complex z; - double r, phi; - if (!PyArg_ParseTuple(args, "D:polar", &z)) - return NULL; - PyFPE_START_PROTECT("polar function", return 0) - phi = c_atan2(z); /* should not cause any exception */ - r = c_abs(z); /* sets errno to ERANGE on overflow; otherwise 0 */ - PyFPE_END_PROTECT(r) - if (errno != 0) - return math_error(); - else - return Py_BuildValue("dd", r, phi); + Py_complex z; + double r, phi; + if (!PyArg_ParseTuple(args, "D:polar", &z)) + return NULL; + PyFPE_START_PROTECT("polar function", return 0) + phi = c_atan2(z); /* should not cause any exception */ + r = c_abs(z); /* sets errno to ERANGE on overflow; otherwise 0 */ + PyFPE_END_PROTECT(r) + if (errno != 0) + return math_error(); + else + return Py_BuildValue("dd", r, phi); } PyDoc_STRVAR(cmath_polar_doc, @@ -971,51 +971,51 @@ static Py_complex rect_special_values[7][7]; static PyObject * cmath_rect(PyObject *self, PyObject *args) { - Py_complex z; - double r, phi; - if (!PyArg_ParseTuple(args, "dd:rect", &r, &phi)) - return NULL; - errno = 0; - PyFPE_START_PROTECT("rect function", return 0) - - /* deal with special values */ - if (!Py_IS_FINITE(r) || !Py_IS_FINITE(phi)) { - /* if r is +/-infinity and phi is finite but nonzero then - result is (+-INF +-INF i), but we need to compute cos(phi) - and sin(phi) to figure out the signs. */ - if (Py_IS_INFINITY(r) && (Py_IS_FINITE(phi) - && (phi != 0.))) { - if (r > 0) { - z.real = copysign(INF, cos(phi)); - z.imag = copysign(INF, sin(phi)); - } - else { - z.real = -copysign(INF, cos(phi)); - z.imag = -copysign(INF, sin(phi)); - } - } - else { - z = rect_special_values[special_type(r)] - [special_type(phi)]; - } - /* need to set errno = EDOM if r is a nonzero number and phi - is infinite */ - if (r != 0. && !Py_IS_NAN(r) && Py_IS_INFINITY(phi)) - errno = EDOM; - else - errno = 0; - } - else { - z.real = r * cos(phi); - z.imag = r * sin(phi); - errno = 0; - } - - PyFPE_END_PROTECT(z) - if (errno != 0) - return math_error(); - else - return PyComplex_FromCComplex(z); + Py_complex z; + double r, phi; + if (!PyArg_ParseTuple(args, "dd:rect", &r, &phi)) + return NULL; + errno = 0; + PyFPE_START_PROTECT("rect function", return 0) + + /* deal with special values */ + if (!Py_IS_FINITE(r) || !Py_IS_FINITE(phi)) { + /* if r is +/-infinity and phi is finite but nonzero then + result is (+-INF +-INF i), but we need to compute cos(phi) + and sin(phi) to figure out the signs. */ + if (Py_IS_INFINITY(r) && (Py_IS_FINITE(phi) + && (phi != 0.))) { + if (r > 0) { + z.real = copysign(INF, cos(phi)); + z.imag = copysign(INF, sin(phi)); + } + else { + z.real = -copysign(INF, cos(phi)); + z.imag = -copysign(INF, sin(phi)); + } + } + else { + z = rect_special_values[special_type(r)] + [special_type(phi)]; + } + /* need to set errno = EDOM if r is a nonzero number and phi + is infinite */ + if (r != 0. && !Py_IS_NAN(r) && Py_IS_INFINITY(phi)) + errno = EDOM; + else + errno = 0; + } + else { + z.real = r * cos(phi); + z.imag = r * sin(phi); + errno = 0; + } + + PyFPE_END_PROTECT(z) + if (errno != 0) + return math_error(); + else + return PyComplex_FromCComplex(z); } PyDoc_STRVAR(cmath_rect_doc, @@ -1025,10 +1025,10 @@ Convert from polar coordinates to rectangular coordinates."); static PyObject * cmath_isnan(PyObject *self, PyObject *args) { - Py_complex z; - if (!PyArg_ParseTuple(args, "D:isnan", &z)) - return NULL; - return PyBool_FromLong(Py_IS_NAN(z.real) || Py_IS_NAN(z.imag)); + Py_complex z; + if (!PyArg_ParseTuple(args, "D:isnan", &z)) + return NULL; + return PyBool_FromLong(Py_IS_NAN(z.real) || Py_IS_NAN(z.imag)); } PyDoc_STRVAR(cmath_isnan_doc, @@ -1038,11 +1038,11 @@ Checks if the real or imaginary part of z not a number (NaN)"); static PyObject * cmath_isinf(PyObject *self, PyObject *args) { - Py_complex z; - if (!PyArg_ParseTuple(args, "D:isnan", &z)) - return NULL; - return PyBool_FromLong(Py_IS_INFINITY(z.real) || - Py_IS_INFINITY(z.imag)); + Py_complex z; + if (!PyArg_ParseTuple(args, "D:isnan", &z)) + return NULL; + return PyBool_FromLong(Py_IS_INFINITY(z.real) || + Py_IS_INFINITY(z.imag)); } PyDoc_STRVAR(cmath_isinf_doc, @@ -1055,169 +1055,169 @@ PyDoc_STRVAR(module_doc, "functions for complex numbers."); static PyMethodDef cmath_methods[] = { - {"acos", cmath_acos, METH_VARARGS, c_acos_doc}, - {"acosh", cmath_acosh, METH_VARARGS, c_acosh_doc}, - {"asin", cmath_asin, METH_VARARGS, c_asin_doc}, - {"asinh", cmath_asinh, METH_VARARGS, c_asinh_doc}, - {"atan", cmath_atan, METH_VARARGS, c_atan_doc}, - {"atanh", cmath_atanh, METH_VARARGS, c_atanh_doc}, - {"cos", cmath_cos, METH_VARARGS, c_cos_doc}, - {"cosh", cmath_cosh, METH_VARARGS, c_cosh_doc}, - {"exp", cmath_exp, METH_VARARGS, c_exp_doc}, - {"isinf", cmath_isinf, METH_VARARGS, cmath_isinf_doc}, - {"isnan", cmath_isnan, METH_VARARGS, cmath_isnan_doc}, - {"log", cmath_log, METH_VARARGS, cmath_log_doc}, - {"log10", cmath_log10, METH_VARARGS, c_log10_doc}, - {"phase", cmath_phase, METH_VARARGS, cmath_phase_doc}, - {"polar", cmath_polar, METH_VARARGS, cmath_polar_doc}, - {"rect", cmath_rect, METH_VARARGS, cmath_rect_doc}, - {"sin", cmath_sin, METH_VARARGS, c_sin_doc}, - {"sinh", cmath_sinh, METH_VARARGS, c_sinh_doc}, - {"sqrt", cmath_sqrt, METH_VARARGS, c_sqrt_doc}, - {"tan", cmath_tan, METH_VARARGS, c_tan_doc}, - {"tanh", cmath_tanh, METH_VARARGS, c_tanh_doc}, - {NULL, NULL} /* sentinel */ + {"acos", cmath_acos, METH_VARARGS, c_acos_doc}, + {"acosh", cmath_acosh, METH_VARARGS, c_acosh_doc}, + {"asin", cmath_asin, METH_VARARGS, c_asin_doc}, + {"asinh", cmath_asinh, METH_VARARGS, c_asinh_doc}, + {"atan", cmath_atan, METH_VARARGS, c_atan_doc}, + {"atanh", cmath_atanh, METH_VARARGS, c_atanh_doc}, + {"cos", cmath_cos, METH_VARARGS, c_cos_doc}, + {"cosh", cmath_cosh, METH_VARARGS, c_cosh_doc}, + {"exp", cmath_exp, METH_VARARGS, c_exp_doc}, + {"isinf", cmath_isinf, METH_VARARGS, cmath_isinf_doc}, + {"isnan", cmath_isnan, METH_VARARGS, cmath_isnan_doc}, + {"log", cmath_log, METH_VARARGS, cmath_log_doc}, + {"log10", cmath_log10, METH_VARARGS, c_log10_doc}, + {"phase", cmath_phase, METH_VARARGS, cmath_phase_doc}, + {"polar", cmath_polar, METH_VARARGS, cmath_polar_doc}, + {"rect", cmath_rect, METH_VARARGS, cmath_rect_doc}, + {"sin", cmath_sin, METH_VARARGS, c_sin_doc}, + {"sinh", cmath_sinh, METH_VARARGS, c_sinh_doc}, + {"sqrt", cmath_sqrt, METH_VARARGS, c_sqrt_doc}, + {"tan", cmath_tan, METH_VARARGS, c_tan_doc}, + {"tanh", cmath_tanh, METH_VARARGS, c_tanh_doc}, + {NULL, NULL} /* sentinel */ }; static struct PyModuleDef cmathmodule = { - PyModuleDef_HEAD_INIT, - "cmath", - module_doc, - -1, - cmath_methods, - NULL, - NULL, - NULL, - NULL + PyModuleDef_HEAD_INIT, + "cmath", + module_doc, + -1, + cmath_methods, + NULL, + NULL, + NULL, + NULL }; PyMODINIT_FUNC PyInit_cmath(void) { - PyObject *m; + PyObject *m; - m = PyModule_Create(&cmathmodule); - if (m == NULL) - return NULL; + m = PyModule_Create(&cmathmodule); + if (m == NULL) + return NULL; - PyModule_AddObject(m, "pi", - PyFloat_FromDouble(Py_MATH_PI)); - PyModule_AddObject(m, "e", PyFloat_FromDouble(Py_MATH_E)); + PyModule_AddObject(m, "pi", + PyFloat_FromDouble(Py_MATH_PI)); + PyModule_AddObject(m, "e", PyFloat_FromDouble(Py_MATH_E)); - /* initialize special value tables */ + /* initialize special value tables */ #define INIT_SPECIAL_VALUES(NAME, BODY) { Py_complex* p = (Py_complex*)NAME; BODY } #define C(REAL, IMAG) p->real = REAL; p->imag = IMAG; ++p; - INIT_SPECIAL_VALUES(acos_special_values, { - C(P34,INF) C(P,INF) C(P,INF) C(P,-INF) C(P,-INF) C(P34,-INF) C(N,INF) - C(P12,INF) C(U,U) C(U,U) C(U,U) C(U,U) C(P12,-INF) C(N,N) - C(P12,INF) C(U,U) C(P12,0.) C(P12,-0.) C(U,U) C(P12,-INF) C(P12,N) - C(P12,INF) C(U,U) C(P12,0.) C(P12,-0.) C(U,U) C(P12,-INF) C(P12,N) - C(P12,INF) C(U,U) C(U,U) C(U,U) C(U,U) C(P12,-INF) C(N,N) - C(P14,INF) C(0.,INF) C(0.,INF) C(0.,-INF) C(0.,-INF) C(P14,-INF) C(N,INF) - C(N,INF) C(N,N) C(N,N) C(N,N) C(N,N) C(N,-INF) C(N,N) - }) - - INIT_SPECIAL_VALUES(acosh_special_values, { - C(INF,-P34) C(INF,-P) C(INF,-P) C(INF,P) C(INF,P) C(INF,P34) C(INF,N) - C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) - C(INF,-P12) C(U,U) C(0.,-P12) C(0.,P12) C(U,U) C(INF,P12) C(N,N) - C(INF,-P12) C(U,U) C(0.,-P12) C(0.,P12) C(U,U) C(INF,P12) C(N,N) - C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) - C(INF,-P14) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,P14) C(INF,N) - C(INF,N) C(N,N) C(N,N) C(N,N) C(N,N) C(INF,N) C(N,N) - }) - - INIT_SPECIAL_VALUES(asinh_special_values, { - C(-INF,-P14) C(-INF,-0.) C(-INF,-0.) C(-INF,0.) C(-INF,0.) C(-INF,P14) C(-INF,N) - C(-INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(-INF,P12) C(N,N) - C(-INF,-P12) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(-INF,P12) C(N,N) - C(INF,-P12) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(INF,P12) C(N,N) - C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) - C(INF,-P14) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,P14) C(INF,N) - C(INF,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(INF,N) C(N,N) - }) - - INIT_SPECIAL_VALUES(atanh_special_values, { - C(-0.,-P12) C(-0.,-P12) C(-0.,-P12) C(-0.,P12) C(-0.,P12) C(-0.,P12) C(-0.,N) - C(-0.,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(-0.,P12) C(N,N) - C(-0.,-P12) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(-0.,P12) C(-0.,N) - C(0.,-P12) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,P12) C(0.,N) - C(0.,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(0.,P12) C(N,N) - C(0.,-P12) C(0.,-P12) C(0.,-P12) C(0.,P12) C(0.,P12) C(0.,P12) C(0.,N) - C(0.,-P12) C(N,N) C(N,N) C(N,N) C(N,N) C(0.,P12) C(N,N) - }) - - INIT_SPECIAL_VALUES(cosh_special_values, { - C(INF,N) C(U,U) C(INF,0.) C(INF,-0.) C(U,U) C(INF,N) C(INF,N) - C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) - C(N,0.) C(U,U) C(1.,0.) C(1.,-0.) C(U,U) C(N,0.) C(N,0.) - C(N,0.) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(N,0.) C(N,0.) - C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) - C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N) - C(N,N) C(N,N) C(N,0.) C(N,0.) C(N,N) C(N,N) C(N,N) - }) - - INIT_SPECIAL_VALUES(exp_special_values, { - C(0.,0.) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,0.) C(0.,0.) - C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) - C(N,N) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(N,N) C(N,N) - C(N,N) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(N,N) C(N,N) - C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) - C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N) - C(N,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(N,N) C(N,N) - }) - - INIT_SPECIAL_VALUES(log_special_values, { - C(INF,-P34) C(INF,-P) C(INF,-P) C(INF,P) C(INF,P) C(INF,P34) C(INF,N) - C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) - C(INF,-P12) C(U,U) C(-INF,-P) C(-INF,P) C(U,U) C(INF,P12) C(N,N) - C(INF,-P12) C(U,U) C(-INF,-0.) C(-INF,0.) C(U,U) C(INF,P12) C(N,N) - C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) - C(INF,-P14) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,P14) C(INF,N) - C(INF,N) C(N,N) C(N,N) C(N,N) C(N,N) C(INF,N) C(N,N) - }) - - INIT_SPECIAL_VALUES(sinh_special_values, { - C(INF,N) C(U,U) C(-INF,-0.) C(-INF,0.) C(U,U) C(INF,N) C(INF,N) - C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) - C(0.,N) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(0.,N) C(0.,N) - C(0.,N) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,N) C(0.,N) - C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) - C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N) - C(N,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(N,N) C(N,N) - }) - - INIT_SPECIAL_VALUES(sqrt_special_values, { - C(INF,-INF) C(0.,-INF) C(0.,-INF) C(0.,INF) C(0.,INF) C(INF,INF) C(N,INF) - C(INF,-INF) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,INF) C(N,N) - C(INF,-INF) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(INF,INF) C(N,N) - C(INF,-INF) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(INF,INF) C(N,N) - C(INF,-INF) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,INF) C(N,N) - C(INF,-INF) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,INF) C(INF,N) - C(INF,-INF) C(N,N) C(N,N) C(N,N) C(N,N) C(INF,INF) C(N,N) - }) - - INIT_SPECIAL_VALUES(tanh_special_values, { - C(-1.,0.) C(U,U) C(-1.,-0.) C(-1.,0.) C(U,U) C(-1.,0.) C(-1.,0.) - C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) - C(N,N) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(N,N) C(N,N) - C(N,N) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(N,N) C(N,N) - C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) - C(1.,0.) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(1.,0.) C(1.,0.) - C(N,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(N,N) C(N,N) - }) - - INIT_SPECIAL_VALUES(rect_special_values, { - C(INF,N) C(U,U) C(-INF,0.) C(-INF,-0.) C(U,U) C(INF,N) C(INF,N) - C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) - C(0.,0.) C(U,U) C(-0.,0.) C(-0.,-0.) C(U,U) C(0.,0.) C(0.,0.) - C(0.,0.) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,0.) C(0.,0.) - C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) - C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N) - C(N,N) C(N,N) C(N,0.) C(N,0.) C(N,N) C(N,N) C(N,N) - }) - return m; + INIT_SPECIAL_VALUES(acos_special_values, { + C(P34,INF) C(P,INF) C(P,INF) C(P,-INF) C(P,-INF) C(P34,-INF) C(N,INF) + C(P12,INF) C(U,U) C(U,U) C(U,U) C(U,U) C(P12,-INF) C(N,N) + C(P12,INF) C(U,U) C(P12,0.) C(P12,-0.) C(U,U) C(P12,-INF) C(P12,N) + C(P12,INF) C(U,U) C(P12,0.) C(P12,-0.) C(U,U) C(P12,-INF) C(P12,N) + C(P12,INF) C(U,U) C(U,U) C(U,U) C(U,U) C(P12,-INF) C(N,N) + C(P14,INF) C(0.,INF) C(0.,INF) C(0.,-INF) C(0.,-INF) C(P14,-INF) C(N,INF) + C(N,INF) C(N,N) C(N,N) C(N,N) C(N,N) C(N,-INF) C(N,N) + }) + + INIT_SPECIAL_VALUES(acosh_special_values, { + C(INF,-P34) C(INF,-P) C(INF,-P) C(INF,P) C(INF,P) C(INF,P34) C(INF,N) + C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) + C(INF,-P12) C(U,U) C(0.,-P12) C(0.,P12) C(U,U) C(INF,P12) C(N,N) + C(INF,-P12) C(U,U) C(0.,-P12) C(0.,P12) C(U,U) C(INF,P12) C(N,N) + C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) + C(INF,-P14) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,P14) C(INF,N) + C(INF,N) C(N,N) C(N,N) C(N,N) C(N,N) C(INF,N) C(N,N) + }) + + INIT_SPECIAL_VALUES(asinh_special_values, { + C(-INF,-P14) C(-INF,-0.) C(-INF,-0.) C(-INF,0.) C(-INF,0.) C(-INF,P14) C(-INF,N) + C(-INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(-INF,P12) C(N,N) + C(-INF,-P12) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(-INF,P12) C(N,N) + C(INF,-P12) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(INF,P12) C(N,N) + C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) + C(INF,-P14) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,P14) C(INF,N) + C(INF,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(INF,N) C(N,N) + }) + + INIT_SPECIAL_VALUES(atanh_special_values, { + C(-0.,-P12) C(-0.,-P12) C(-0.,-P12) C(-0.,P12) C(-0.,P12) C(-0.,P12) C(-0.,N) + C(-0.,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(-0.,P12) C(N,N) + C(-0.,-P12) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(-0.,P12) C(-0.,N) + C(0.,-P12) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,P12) C(0.,N) + C(0.,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(0.,P12) C(N,N) + C(0.,-P12) C(0.,-P12) C(0.,-P12) C(0.,P12) C(0.,P12) C(0.,P12) C(0.,N) + C(0.,-P12) C(N,N) C(N,N) C(N,N) C(N,N) C(0.,P12) C(N,N) + }) + + INIT_SPECIAL_VALUES(cosh_special_values, { + C(INF,N) C(U,U) C(INF,0.) C(INF,-0.) C(U,U) C(INF,N) C(INF,N) + C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) + C(N,0.) C(U,U) C(1.,0.) C(1.,-0.) C(U,U) C(N,0.) C(N,0.) + C(N,0.) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(N,0.) C(N,0.) + C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) + C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N) + C(N,N) C(N,N) C(N,0.) C(N,0.) C(N,N) C(N,N) C(N,N) + }) + + INIT_SPECIAL_VALUES(exp_special_values, { + C(0.,0.) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,0.) C(0.,0.) + C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) + C(N,N) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(N,N) C(N,N) + C(N,N) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(N,N) C(N,N) + C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) + C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N) + C(N,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(N,N) C(N,N) + }) + + INIT_SPECIAL_VALUES(log_special_values, { + C(INF,-P34) C(INF,-P) C(INF,-P) C(INF,P) C(INF,P) C(INF,P34) C(INF,N) + C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) + C(INF,-P12) C(U,U) C(-INF,-P) C(-INF,P) C(U,U) C(INF,P12) C(N,N) + C(INF,-P12) C(U,U) C(-INF,-0.) C(-INF,0.) C(U,U) C(INF,P12) C(N,N) + C(INF,-P12) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,P12) C(N,N) + C(INF,-P14) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,P14) C(INF,N) + C(INF,N) C(N,N) C(N,N) C(N,N) C(N,N) C(INF,N) C(N,N) + }) + + INIT_SPECIAL_VALUES(sinh_special_values, { + C(INF,N) C(U,U) C(-INF,-0.) C(-INF,0.) C(U,U) C(INF,N) C(INF,N) + C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) + C(0.,N) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(0.,N) C(0.,N) + C(0.,N) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,N) C(0.,N) + C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) + C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N) + C(N,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(N,N) C(N,N) + }) + + INIT_SPECIAL_VALUES(sqrt_special_values, { + C(INF,-INF) C(0.,-INF) C(0.,-INF) C(0.,INF) C(0.,INF) C(INF,INF) C(N,INF) + C(INF,-INF) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,INF) C(N,N) + C(INF,-INF) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(INF,INF) C(N,N) + C(INF,-INF) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(INF,INF) C(N,N) + C(INF,-INF) C(U,U) C(U,U) C(U,U) C(U,U) C(INF,INF) C(N,N) + C(INF,-INF) C(INF,-0.) C(INF,-0.) C(INF,0.) C(INF,0.) C(INF,INF) C(INF,N) + C(INF,-INF) C(N,N) C(N,N) C(N,N) C(N,N) C(INF,INF) C(N,N) + }) + + INIT_SPECIAL_VALUES(tanh_special_values, { + C(-1.,0.) C(U,U) C(-1.,-0.) C(-1.,0.) C(U,U) C(-1.,0.) C(-1.,0.) + C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) + C(N,N) C(U,U) C(-0.,-0.) C(-0.,0.) C(U,U) C(N,N) C(N,N) + C(N,N) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(N,N) C(N,N) + C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) + C(1.,0.) C(U,U) C(1.,-0.) C(1.,0.) C(U,U) C(1.,0.) C(1.,0.) + C(N,N) C(N,N) C(N,-0.) C(N,0.) C(N,N) C(N,N) C(N,N) + }) + + INIT_SPECIAL_VALUES(rect_special_values, { + C(INF,N) C(U,U) C(-INF,0.) C(-INF,-0.) C(U,U) C(INF,N) C(INF,N) + C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) + C(0.,0.) C(U,U) C(-0.,0.) C(-0.,-0.) C(U,U) C(0.,0.) C(0.,0.) + C(0.,0.) C(U,U) C(0.,-0.) C(0.,0.) C(U,U) C(0.,0.) C(0.,0.) + C(N,N) C(U,U) C(U,U) C(U,U) C(U,U) C(N,N) C(N,N) + C(INF,N) C(U,U) C(INF,-0.) C(INF,0.) C(U,U) C(INF,N) C(INF,N) + C(N,N) C(N,N) C(N,0.) C(N,0.) C(N,N) C(N,N) C(N,N) + }) + return m; } |