diff options
Diffstat (limited to 'math/k_casinhl.c')
| -rw-r--r-- | math/k_casinhl.c | 33 |
1 files changed, 33 insertions, 0 deletions
diff --git a/math/k_casinhl.c b/math/k_casinhl.c index eb090ec3ce..8e5bbd4385 100644 --- a/math/k_casinhl.c +++ b/math/k_casinhl.c @@ -84,6 +84,39 @@ __kernel_casinhl (__complex__ long double x, int adj) else __imag__ res = __ieee754_atan2l (s, rx); } + else if (ix > 1.0L && ix < 1.5L && rx < 0.5L) + { + if (rx < LDBL_EPSILON * LDBL_EPSILON) + { + long double ix2m1 = (ix + 1.0L) * (ix - 1.0L); + long double s = __ieee754_sqrtl (ix2m1); + + __real__ res = __log1pl (2.0L * (ix2m1 + ix * s)) / 2.0L; + if (adj) + __imag__ res = __ieee754_atan2l (rx, __copysignl (s, __imag__ x)); + else + __imag__ res = __ieee754_atan2l (s, rx); + } + else + { + long double ix2m1 = (ix + 1.0L) * (ix - 1.0L); + long double rx2 = rx * rx; + long double f = rx2 * (2.0L + rx2 + 2.0L * ix * ix); + long double d = __ieee754_sqrtl (ix2m1 * ix2m1 + f); + long double dp = d + ix2m1; + long double dm = f / dp; + long double r1 = __ieee754_sqrtl ((dm + rx2) / 2.0L); + long double r2 = rx * ix / r1; + + __real__ res + = __log1pl (rx2 + dp + 2.0L * (rx * r1 + ix * r2)) / 2.0L; + if (adj) + __imag__ res = __ieee754_atan2l (rx + r1, __copysignl (ix + r2, + __imag__ x)); + else + __imag__ res = __ieee754_atan2l (ix + r2, rx + r1); + } + } else if (ix == 1.0L && rx < 0.5L) { if (rx < LDBL_EPSILON / 8.0L) |