]> git.ipfire.org Git - thirdparty/glibc.git/blobdiff - math/s_catanh.c
Update copyright dates with scripts/update-copyrights.
[thirdparty/glibc.git] / math / s_catanh.c
index 55c56283b3d5c459246dc7687ed1ca6ab9446f8b..281183efddc62565c892519b96969a76ad09c4ca 100644 (file)
@@ -1,5 +1,5 @@
 /* Return arc hyperbole tangent for double value.
-   Copyright (C) 1997-2013 Free Software Foundation, Inc.
+   Copyright (C) 1997-2015 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
 
@@ -29,7 +29,7 @@ __catanh (__complex__ double x)
   int rcls = fpclassify (__real__ x);
   int icls = fpclassify (__imag__ x);
 
-  if (__builtin_expect (rcls <= FP_INFINITE || icls <= FP_INFINITE, 0))
+  if (__glibc_unlikely (rcls <= FP_INFINITE || icls <= FP_INFINITE))
     {
       if (icls == FP_INFINITE)
        {
@@ -50,7 +50,7 @@ __catanh (__complex__ double x)
          __imag__ res = __nan ("");
        }
     }
-  else if (__builtin_expect (rcls == FP_ZERO && icls == FP_ZERO, 0))
+  else if (__glibc_unlikely (rcls == FP_ZERO && icls == FP_ZERO))
     {
       res = x;
     }
@@ -72,24 +72,55 @@ __catanh (__complex__ double x)
        }
       else
        {
-         double i2 = __imag__ x * __imag__ x;
+         if (fabs (__real__ x) == 1.0
+             && fabs (__imag__ x) < DBL_EPSILON * DBL_EPSILON)
+           __real__ res = (__copysign (0.5, __real__ x)
+                           * (M_LN2 - __ieee754_log (fabs (__imag__ x))));
+         else
+           {
+             double i2 = 0.0;
+             if (fabs (__imag__ x) >= DBL_EPSILON * DBL_EPSILON)
+               i2 = __imag__ x * __imag__ x;
 
-         double num = 1.0 + __real__ x;
-         num = i2 + num * num;
+             double num = 1.0 + __real__ x;
+             num = i2 + num * num;
 
-         double den = 1.0 - __real__ x;
-         den = i2 + den * den;
+             double den = 1.0 - __real__ x;
+             den = i2 + den * den;
 
-         double f = num / den;
-         if (f < 0.5)
-           __real__ res = 0.25 * __ieee754_log (f);
-         else
+             double f = num / den;
+             if (f < 0.5)
+               __real__ res = 0.25 * __ieee754_log (f);
+             else
+               {
+                 num = 4.0 * __real__ x;
+                 __real__ res = 0.25 * __log1p (num / den);
+               }
+           }
+
+         double absx, absy, den;
+
+         absx = fabs (__real__ x);
+         absy = fabs (__imag__ x);
+         if (absx < absy)
            {
-             num = 4.0 * __real__ x;
-             __real__ res = 0.25 * __log1p (num / den);
+             double t = absx;
+             absx = absy;
+             absy = t;
            }
 
-         den = 1 - __real__ x * __real__ x - i2;
+         if (absy < DBL_EPSILON / 2.0)
+           {
+             den = (1.0 - absx) * (1.0 + absx);
+             if (den == -0.0)
+               den = 0.0;
+           }
+         else if (absx >= 1.0)
+           den = (1.0 - absx) * (1.0 + absx) - absy * absy;
+         else if (absx >= 0.75 || absy >= 0.5)
+           den = -__x2y2m1 (absx, absy);
+         else
+           den = (1.0 - absx) * (1.0 + absx) - absy * absy;
 
          __imag__ res = 0.5 * __ieee754_atan2 (2.0 * __imag__ x, den);
        }