about summary refs log tree commit diff
path: root/sysdeps
diff options
context:
space:
mode:
Diffstat (limited to 'sysdeps')
-rw-r--r--sysdeps/generic/math_private.h4
-rw-r--r--sysdeps/i386/fpu/libm-test-ulps80
-rw-r--r--sysdeps/ieee754/dbl-64/x2y2m1.c22
-rw-r--r--sysdeps/ieee754/dbl-64/x2y2m1f.c4
-rw-r--r--sysdeps/ieee754/ldbl-128/x2y2m1l.c22
-rw-r--r--sysdeps/ieee754/ldbl-128ibm/x2y2m1l.c24
-rw-r--r--sysdeps/ieee754/ldbl-96/x2y2m1.c4
-rw-r--r--sysdeps/ieee754/ldbl-96/x2y2m1l.c22
-rw-r--r--sysdeps/x86_64/fpu/libm-test-ulps80
9 files changed, 119 insertions, 143 deletions
diff --git a/sysdeps/generic/math_private.h b/sysdeps/generic/math_private.h
index a8f1a8e4e5..cf1865dac8 100644
--- a/sysdeps/generic/math_private.h
+++ b/sysdeps/generic/math_private.h
@@ -365,8 +365,8 @@ extern double __slowpow (double __x, double __y, double __z);
 extern void __docos (double __x, double __dx, double __v[]);
 
 /* Return X^2 + Y^2 - 1, computed without large cancellation error.
-   It is given that 1 > X >= Y >= epsilon / 2, and that either X >=
-   0.75 or Y >= 0.5.  */
+   It is given that 1 > X >= Y >= epsilon / 2, and that X^2 + Y^2 >=
+   0.5.  */
 extern float __x2y2m1f (float x, float y);
 extern double __x2y2m1 (double x, double y);
 extern long double __x2y2m1l (long double x, long double y);
diff --git a/sysdeps/i386/fpu/libm-test-ulps b/sysdeps/i386/fpu/libm-test-ulps
index 32f24d07a2..438a390338 100644
--- a/sysdeps/i386/fpu/libm-test-ulps
+++ b/sysdeps/i386/fpu/libm-test-ulps
@@ -836,12 +836,12 @@ ildouble: 3
 ldouble: 3
 
 Function: Real part of "clog":
-double: 3
-float: 2
-idouble: 3
-ifloat: 2
-ildouble: 4
-ldouble: 4
+double: 2
+float: 1
+idouble: 2
+ifloat: 1
+ildouble: 3
+ldouble: 3
 
 Function: Imaginary part of "clog":
 double: 1
@@ -864,10 +864,10 @@ ildouble: 2
 ldouble: 2
 
 Function: Real part of "clog10_downward":
-double: 5
-float: 4
-idouble: 5
-ifloat: 4
+double: 3
+float: 3
+idouble: 3
+ifloat: 3
 ildouble: 8
 ldouble: 8
 
@@ -876,14 +876,14 @@ double: 1
 float: 1
 idouble: 1
 ifloat: 1
-ildouble: 2
-ldouble: 2
+ildouble: 3
+ldouble: 3
 
 Function: Real part of "clog10_towardzero":
-double: 5
-float: 4
-idouble: 5
-ifloat: 4
+double: 3
+float: 3
+idouble: 3
+ifloat: 3
 ildouble: 8
 ldouble: 8
 
@@ -896,12 +896,12 @@ ildouble: 3
 ldouble: 3
 
 Function: Real part of "clog10_upward":
-double: 5
-float: 5
-idouble: 5
-ifloat: 5
-ildouble: 6
-ldouble: 6
+double: 3
+float: 3
+idouble: 3
+ifloat: 3
+ildouble: 7
+ldouble: 7
 
 Function: Imaginary part of "clog10_upward":
 double: 1
@@ -912,12 +912,12 @@ ildouble: 3
 ldouble: 3
 
 Function: Real part of "clog_downward":
-double: 5
-float: 5
-idouble: 5
-ifloat: 5
-ildouble: 7
-ldouble: 7
+double: 3
+float: 3
+idouble: 3
+ifloat: 3
+ildouble: 5
+ldouble: 5
 
 Function: Imaginary part of "clog_downward":
 double: 1
@@ -928,12 +928,12 @@ ildouble: 1
 ldouble: 1
 
 Function: Real part of "clog_towardzero":
-double: 5
-float: 5
-idouble: 5
-ifloat: 5
-ildouble: 8
-ldouble: 8
+double: 3
+float: 3
+idouble: 3
+ifloat: 3
+ildouble: 5
+ldouble: 5
 
 Function: Imaginary part of "clog_towardzero":
 double: 1
@@ -944,12 +944,12 @@ ildouble: 1
 ldouble: 1
 
 Function: Real part of "clog_upward":
-double: 5
-float: 5
-idouble: 5
-ifloat: 5
-ildouble: 6
-ldouble: 6
+double: 2
+float: 3
+idouble: 2
+ifloat: 3
+ildouble: 4
+ldouble: 4
 
 Function: Imaginary part of "clog_upward":
 double: 1
diff --git a/sysdeps/ieee754/dbl-64/x2y2m1.c b/sysdeps/ieee754/dbl-64/x2y2m1.c
index c96dae59e1..b040097d18 100644
--- a/sysdeps/ieee754/dbl-64/x2y2m1.c
+++ b/sysdeps/ieee754/dbl-64/x2y2m1.c
@@ -80,32 +80,26 @@ compare (const void *p, const void *q)
 }
 
 /* Return X^2 + Y^2 - 1, computed without large cancellation error.
-   It is given that 1 > X >= Y >= epsilon / 2, and that either X >=
-   0.75 or Y >= 0.5.  */
+   It is given that 1 > X >= Y >= epsilon / 2, and that X^2 + Y^2 >=
+   0.5.  */
 
 double
 __x2y2m1 (double x, double y)
 {
-  double vals[4];
+  double vals[5];
   SET_RESTORE_ROUND (FE_TONEAREST);
   mul_split (&vals[1], &vals[0], x, x);
   mul_split (&vals[3], &vals[2], y, y);
-  if (x >= 0.75)
-    vals[1] -= 1.0;
-  else
-    {
-      vals[1] -= 0.5;
-      vals[3] -= 0.5;
-    }
-  qsort (vals, 4, sizeof (double), compare);
+  vals[4] = -1.0;
+  qsort (vals, 5, sizeof (double), compare);
   /* Add up the values so that each element of VALS has absolute value
      at most equal to the last set bit of the next nonzero
      element.  */
-  for (size_t i = 0; i <= 2; i++)
+  for (size_t i = 0; i <= 3; i++)
     {
       add_split (&vals[i + 1], &vals[i], vals[i + 1], vals[i]);
-      qsort (vals + i + 1, 3 - i, sizeof (double), compare);
+      qsort (vals + i + 1, 4 - i, sizeof (double), compare);
     }
   /* Now any error from this addition will be small.  */
-  return vals[3] + vals[2] + vals[1] + vals[0];
+  return vals[4] + vals[3] + vals[2] + vals[1] + vals[0];
 }
diff --git a/sysdeps/ieee754/dbl-64/x2y2m1f.c b/sysdeps/ieee754/dbl-64/x2y2m1f.c
index 43a8acf8a8..835f6a0e45 100644
--- a/sysdeps/ieee754/dbl-64/x2y2m1f.c
+++ b/sysdeps/ieee754/dbl-64/x2y2m1f.c
@@ -21,8 +21,8 @@
 #include <float.h>
 
 /* Return X^2 + Y^2 - 1, computed without large cancellation error.
-   It is given that 1 > X >= Y >= epsilon / 2, and that either X >=
-   0.75 or Y >= 0.5.  */
+   It is given that 1 > X >= Y >= epsilon / 2, and that X^2 + Y^2 >=
+   0.5.  */
 
 float
 __x2y2m1f (float x, float y)
diff --git a/sysdeps/ieee754/ldbl-128/x2y2m1l.c b/sysdeps/ieee754/ldbl-128/x2y2m1l.c
index 11757c6af8..a0498c30bb 100644
--- a/sysdeps/ieee754/ldbl-128/x2y2m1l.c
+++ b/sysdeps/ieee754/ldbl-128/x2y2m1l.c
@@ -80,32 +80,26 @@ compare (const void *p, const void *q)
 }
 
 /* Return X^2 + Y^2 - 1, computed without large cancellation error.
-   It is given that 1 > X >= Y >= epsilon / 2, and that either X >=
-   0.75 or Y >= 0.5.  */
+   It is given that 1 > X >= Y >= epsilon / 2, and that X^2 + Y^2 >=
+   0.5.  */
 
 long double
 __x2y2m1l (long double x, long double y)
 {
-  long double vals[4];
+  long double vals[5];
   SET_RESTORE_ROUNDL (FE_TONEAREST);
   mul_split (&vals[1], &vals[0], x, x);
   mul_split (&vals[3], &vals[2], y, y);
-  if (x >= 0.75L)
-    vals[1] -= 1.0L;
-  else
-    {
-      vals[1] -= 0.5L;
-      vals[3] -= 0.5L;
-    }
-  qsort (vals, 4, sizeof (long double), compare);
+  vals[4] = -1.0L;
+  qsort (vals, 5, sizeof (long double), compare);
   /* Add up the values so that each element of VALS has absolute value
      at most equal to the last set bit of the next nonzero
      element.  */
-  for (size_t i = 0; i <= 2; i++)
+  for (size_t i = 0; i <= 3; i++)
     {
       add_split (&vals[i + 1], &vals[i], vals[i + 1], vals[i]);
-      qsort (vals + i + 1, 3 - i, sizeof (long double), compare);
+      qsort (vals + i + 1, 4 - i, sizeof (long double), compare);
     }
   /* Now any error from this addition will be small.  */
-  return vals[3] + vals[2] + vals[1] + vals[0];
+  return vals[4] + vals[3] + vals[2] + vals[1] + vals[0];
 }
diff --git a/sysdeps/ieee754/ldbl-128ibm/x2y2m1l.c b/sysdeps/ieee754/ldbl-128ibm/x2y2m1l.c
index a001b58ca5..081fb98255 100644
--- a/sysdeps/ieee754/ldbl-128ibm/x2y2m1l.c
+++ b/sysdeps/ieee754/ldbl-128ibm/x2y2m1l.c
@@ -80,13 +80,13 @@ compare (const void *p, const void *q)
 }
 
 /* Return X^2 + Y^2 - 1, computed without large cancellation error.
-   It is given that 1 > X >= Y >= epsilon / 2, and that either X >=
-   0.75 or Y >= 0.5.  */
+   It is given that 1 > X >= Y >= epsilon / 2, and that X^2 + Y^2 >=
+   0.5.  */
 
 long double
 __x2y2m1l (long double x, long double y)
 {
-  double vals[12];
+  double vals[13];
   SET_RESTORE_ROUND (FE_TONEAREST);
   union ibm_extended_long_double xu, yu;
   xu.ld = x;
@@ -105,25 +105,19 @@ __x2y2m1l (long double x, long double y)
   vals[8] *= 2.0;
   vals[9] *= 2.0;
   mul_split (&vals[11], &vals[10], yu.d[1].d, yu.d[1].d);
-  if (xu.d[0].d >= 0.75)
-    vals[1] -= 1.0;
-  else
-    {
-      vals[1] -= 0.5;
-      vals[7] -= 0.5;
-    }
-  qsort (vals, 12, sizeof (double), compare);
+  vals[12] = -1.0;
+  qsort (vals, 13, sizeof (double), compare);
   /* Add up the values so that each element of VALS has absolute value
      at most equal to the last set bit of the next nonzero
      element.  */
-  for (size_t i = 0; i <= 10; i++)
+  for (size_t i = 0; i <= 11; i++)
     {
       add_split (&vals[i + 1], &vals[i], vals[i + 1], vals[i]);
-      qsort (vals + i + 1, 11 - i, sizeof (double), compare);
+      qsort (vals + i + 1, 12 - i, sizeof (double), compare);
     }
   /* Now any error from this addition will be small.  */
-  long double retval = (long double) vals[11];
-  for (size_t i = 10; i != (size_t) -1; i--)
+  long double retval = (long double) vals[12];
+  for (size_t i = 11; i != (size_t) -1; i--)
     retval += (long double) vals[i];
   return retval;
 }
diff --git a/sysdeps/ieee754/ldbl-96/x2y2m1.c b/sysdeps/ieee754/ldbl-96/x2y2m1.c
index a6cc82cb29..2f6b0be8f9 100644
--- a/sysdeps/ieee754/ldbl-96/x2y2m1.c
+++ b/sysdeps/ieee754/ldbl-96/x2y2m1.c
@@ -27,8 +27,8 @@
 #else
 
 /* Return X^2 + Y^2 - 1, computed without large cancellation error.
-   It is given that 1 > X >= Y >= epsilon / 2, and that either X >=
-   0.75 or Y >= 0.5.  */
+   It is given that 1 > X >= Y >= epsilon / 2, and that X^2 + Y^2 >=
+   0.5.  */
 
 double
 __x2y2m1 (double x, double y)
diff --git a/sysdeps/ieee754/ldbl-96/x2y2m1l.c b/sysdeps/ieee754/ldbl-96/x2y2m1l.c
index 11757c6af8..a0498c30bb 100644
--- a/sysdeps/ieee754/ldbl-96/x2y2m1l.c
+++ b/sysdeps/ieee754/ldbl-96/x2y2m1l.c
@@ -80,32 +80,26 @@ compare (const void *p, const void *q)
 }
 
 /* Return X^2 + Y^2 - 1, computed without large cancellation error.
-   It is given that 1 > X >= Y >= epsilon / 2, and that either X >=
-   0.75 or Y >= 0.5.  */
+   It is given that 1 > X >= Y >= epsilon / 2, and that X^2 + Y^2 >=
+   0.5.  */
 
 long double
 __x2y2m1l (long double x, long double y)
 {
-  long double vals[4];
+  long double vals[5];
   SET_RESTORE_ROUNDL (FE_TONEAREST);
   mul_split (&vals[1], &vals[0], x, x);
   mul_split (&vals[3], &vals[2], y, y);
-  if (x >= 0.75L)
-    vals[1] -= 1.0L;
-  else
-    {
-      vals[1] -= 0.5L;
-      vals[3] -= 0.5L;
-    }
-  qsort (vals, 4, sizeof (long double), compare);
+  vals[4] = -1.0L;
+  qsort (vals, 5, sizeof (long double), compare);
   /* Add up the values so that each element of VALS has absolute value
      at most equal to the last set bit of the next nonzero
      element.  */
-  for (size_t i = 0; i <= 2; i++)
+  for (size_t i = 0; i <= 3; i++)
     {
       add_split (&vals[i + 1], &vals[i], vals[i + 1], vals[i]);
-      qsort (vals + i + 1, 3 - i, sizeof (long double), compare);
+      qsort (vals + i + 1, 4 - i, sizeof (long double), compare);
     }
   /* Now any error from this addition will be small.  */
-  return vals[3] + vals[2] + vals[1] + vals[0];
+  return vals[4] + vals[3] + vals[2] + vals[1] + vals[0];
 }
diff --git a/sysdeps/x86_64/fpu/libm-test-ulps b/sysdeps/x86_64/fpu/libm-test-ulps
index 12c3dd1e01..ef3ab70369 100644
--- a/sysdeps/x86_64/fpu/libm-test-ulps
+++ b/sysdeps/x86_64/fpu/libm-test-ulps
@@ -869,11 +869,11 @@ ldouble: 3
 
 Function: Real part of "clog":
 double: 3
-float: 2
+float: 3
 idouble: 3
-ifloat: 2
-ildouble: 4
-ldouble: 4
+ifloat: 3
+ildouble: 3
+ldouble: 3
 
 Function: Imaginary part of "clog":
 float: 1
@@ -883,9 +883,9 @@ ldouble: 1
 
 Function: Real part of "clog10":
 double: 3
-float: 3
+float: 4
 idouble: 3
-ifloat: 3
+ifloat: 4
 ildouble: 4
 ldouble: 4
 
@@ -898,10 +898,10 @@ ildouble: 2
 ldouble: 2
 
 Function: Real part of "clog10_downward":
-double: 6
-float: 6
-idouble: 6
-ifloat: 6
+double: 5
+float: 4
+idouble: 5
+ifloat: 4
 ildouble: 8
 ldouble: 8
 
@@ -910,14 +910,14 @@ double: 2
 float: 4
 idouble: 2
 ifloat: 4
-ildouble: 2
-ldouble: 2
+ildouble: 3
+ldouble: 3
 
 Function: Real part of "clog10_towardzero":
 double: 5
-float: 4
+float: 5
 idouble: 5
-ifloat: 4
+ifloat: 5
 ildouble: 8
 ldouble: 8
 
@@ -930,28 +930,28 @@ ildouble: 3
 ldouble: 3
 
 Function: Real part of "clog10_upward":
-double: 8
+double: 6
 float: 5
-idouble: 8
+idouble: 6
 ifloat: 5
-ildouble: 6
-ldouble: 6
+ildouble: 7
+ldouble: 7
 
 Function: Imaginary part of "clog10_upward":
 double: 2
-float: 3
+float: 4
 idouble: 2
-ifloat: 3
+ifloat: 4
 ildouble: 3
 ldouble: 3
 
 Function: Real part of "clog_downward":
-double: 7
-float: 5
-idouble: 7
-ifloat: 5
-ildouble: 7
-ldouble: 7
+double: 4
+float: 3
+idouble: 4
+ifloat: 3
+ildouble: 5
+ldouble: 5
 
 Function: Imaginary part of "clog_downward":
 double: 1
@@ -962,28 +962,28 @@ ildouble: 1
 ldouble: 1
 
 Function: Real part of "clog_towardzero":
-double: 7
-float: 5
-idouble: 7
-ifloat: 5
-ildouble: 8
-ldouble: 8
+double: 4
+float: 4
+idouble: 4
+ifloat: 4
+ildouble: 5
+ldouble: 5
 
 Function: Imaginary part of "clog_towardzero":
 double: 1
-float: 2
+float: 3
 idouble: 1
-ifloat: 2
+ifloat: 3
 ildouble: 1
 ldouble: 1
 
 Function: Real part of "clog_upward":
-double: 8
-float: 5
-idouble: 8
-ifloat: 5
-ildouble: 6
-ldouble: 6
+double: 4
+float: 3
+idouble: 4
+ifloat: 3
+ildouble: 4
+ldouble: 4
 
 Function: Imaginary part of "clog_upward":
 double: 1