Convert TEST_ff_i tests from code to data.

[platform/upstream/glibc.git] / math / s_catan.c

diff --git a/math/s_catan.c b/math/s_catan.c
index 46c18bf..4da0a0d 100644 (file)
--- a/math/s_catan.c
+++ b/math/s_catan.c
@@ -20,7 +20,7 @@
 #include <complex.h>
 #include <math.h>
 #include <math_private.h>
-
+#include <float.h>
 
 __complex__ double
 __catan (__complex__ double x)
@@ -61,21 +61,82 @@ __catan (__complex__ double x)
     }
   else
     {
-      double r2, num, den;
+      if (fabs (__real__ x) >= 16.0 / DBL_EPSILON
+         || fabs (__imag__ x) >= 16.0 / DBL_EPSILON)
+       {
+         __real__ res = __copysign (M_PI_2, __real__ x);
+         if (fabs (__real__ x) <= 1.0)
+           __imag__ res = 1.0 / __imag__ x;
+         else if (fabs (__imag__ x) <= 1.0)
+           __imag__ res = __imag__ x / __real__ x / __real__ x;
+         else
+           {
+             double h = __ieee754_hypot (__real__ x / 2.0, __imag__ x / 2.0);
+             __imag__ res = __imag__ x / h / h / 4.0;
+           }
+       }
+      else
+       {
+         double den, absx, absy;
 
-      r2 = __real__ x * __real__ x;
+         absx = fabs (__real__ x);
+         absy = fabs (__imag__ x);
+         if (absx < absy)
+           {
+             double t = absx;
+             absx = absy;
+             absy = t;
+           }
 
-      den = 1 - r2 - __imag__ x * __imag__ x;
+         if (absy < DBL_EPSILON / 2.0)
+           den = (1.0 - absx) * (1.0 + absx);
+         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;
 
-      __real__ res = 0.5 * __ieee754_atan2 (2.0 * __real__ x, den);
+         __real__ res = 0.5 * __ieee754_atan2 (2.0 * __real__ x, den);
 
-      num = __imag__ x + 1.0;
-      num = r2 + num * num;
+         if (fabs (__imag__ x) == 1.0
+             && fabs (__real__ x) < DBL_EPSILON * DBL_EPSILON)
+           __imag__ res = (__copysign (0.5, __imag__ x)
+                           * (M_LN2 - __ieee754_log (fabs (__real__ x))));
+         else
+           {
+             double r2 = 0.0, num, f;
 
-      den = __imag__ x - 1.0;
-      den = r2 + den * den;
+             if (fabs (__real__ x) >= DBL_EPSILON * DBL_EPSILON)
+               r2 = __real__ x * __real__ x;
 
-      __imag__ res = 0.25 * __ieee754_log (num / den);
+             num = __imag__ x + 1.0;
+             num = r2 + num * num;
+
+             den = __imag__ x - 1.0;
+             den = r2 + den * den;
+
+             f = num / den;
+             if (f < 0.5)
+               __imag__ res = 0.25 * __ieee754_log (f);
+             else
+               {
+                 num = 4.0 * __imag__ x;
+                 __imag__ res = 0.25 * __log1p (num / den);
+               }
+           }
+       }
+
+      if (fabs (__real__ res) < DBL_MIN)
+       {
+         volatile double force_underflow = __real__ res * __real__ res;
+         (void) force_underflow;
+       }
+      if (fabs (__imag__ res) < DBL_MIN)
+       {
+         volatile double force_underflow = __imag__ res * __imag__ res;
+         (void) force_underflow;
+       }
     }
 
   return res;