GNU C Library master sources branch master updated. glibc-2.22-620-gb300455

siddhesh@sourceware.org siddhesh@sourceware.org
Mon Dec 21 05:20:00 GMT 2015


This is an automated email from the git hooks/post-receive script. It was
generated because a ref change was pushed to the repository containing
the project "GNU C Library master sources".

The branch, master has been updated
       via  b300455644e2945da05eb49d12d3a037f1408be1 (commit)
       via  f7953c44d545c91442f123300288853284501dab (commit)
       via  a045832debe4191bcffed6e539a5470152e9c9e8 (commit)
       via  760c2eb7da879e05252468a512a582a0e5c180af (commit)
      from  5537f466d620552e7f24d52a9b97ced736d2abb1 (commit)

Those revisions listed above that are new to this repository have
not appeared on any other notification email; so we list those
revisions in full, below.

- Log -----------------------------------------------------------------
http://sourceware.org/git/gitweb.cgi?p=glibc.git;a=commitdiff;h=b300455644e2945da05eb49d12d3a037f1408be1

commit b300455644e2945da05eb49d12d3a037f1408be1
Author: Siddhesh Poyarekar <siddhesh.poyarekar@linaro.org>
Date:   Mon Dec 21 10:43:04 2015 +0530

    Consolidate sincos computation for 2.426265 < |x| < 105414350
    
    Like the previous change, exploit the fact that computation for sin
    and cos is identical except that it is apart by a quadrant.  Also
    remove csloww, csloww1 and csloww2 since they can easily be expressed
    in terms of sloww, sloww1 and sloww2.

diff --git a/ChangeLog b/ChangeLog
index 2228769..c067997 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,5 +1,14 @@
 2015-12-21  Siddhesh Poyarekar  <siddhesh.poyarekar@linaro.org>
 
+	* sysdeps/ieee754/dbl-64/s_sin.c (csloww, csloww1, csloww2):
+	Remove functions.
+	(sloww, sloww1): Accept argument to offset quadrant.
+	(sloww, sloww1, sloww2): Call __mpsin or __mpcos based on
+	quadrant.
+	(__sin, __cos): Consolidate common code into new functions.
+	(reduce_sincos_1, do_sincos_1): New functions.
+	* sysdeps/ieee754/dbl-64/s_sincos.c (__sincos): Use them.
+
 	* sysdeps/ieee754/dbl-64/s_sin.c (__sin, __cos): Move common
 	code to new functions.
 	(reduce_sincos_2, do_sincos_2): New functions.
diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c
index 3b26a61..9c57321 100644
--- a/sysdeps/ieee754/dbl-64/s_sin.c
+++ b/sysdeps/ieee754/dbl-64/s_sin.c
@@ -32,9 +32,6 @@
 /*            bsloww1                                                       */
 /*            bsloww2                                                       */
 /*            cslow2                                                        */
-/*            csloww                                                        */
-/*            csloww1                                                       */
-/*            csloww2                                                       */
 /* FILES NEEDED: dla.h endian.h mpa.h mydefs.h  usncs.h                     */
 /*               branred.c sincos32.c dosincos.c mpa.c                      */
 /*               sincos.tbl                                                 */
@@ -135,17 +132,14 @@ double __mpcos (double x, double dx, bool reduce_range);
 static double slow (double x);
 static double slow1 (double x);
 static double slow2 (double x);
-static double sloww (double x, double dx, double orig);
-static double sloww1 (double x, double dx, double orig, int m);
+static double sloww (double x, double dx, double orig, int n);
+static double sloww1 (double x, double dx, double orig, int m, int n);
 static double sloww2 (double x, double dx, double orig, int n);
 static double bsloww (double x, double dx, double orig, int n);
 static double bsloww1 (double x, double dx, double orig, int n);
 static double bsloww2 (double x, double dx, double orig, int n);
 int __branred (double x, double *a, double *aa);
 static double cslow2 (double x);
-static double csloww (double x, double dx, double orig);
-static double csloww1 (double x, double dx, double orig, int m);
-static double csloww2 (double x, double dx, double orig, int n);
 
 /* Given a number partitioned into U and X such that U is an index into the
    sin/cos table, this macro computes the cosine of the number by combining
@@ -278,6 +272,91 @@ reduce_and_compute (double x, unsigned int k)
 
 static inline int4
 __always_inline
+reduce_sincos_1 (double x, double *a, double *da)
+{
+  mynumber v;
+
+  double t = (x * hpinv + toint);
+  double xn = t - toint;
+  v.x = t;
+  double y = (x - xn * mp1) - xn * mp2;
+  int4 n = v.i[LOW_HALF] & 3;
+  double db = xn * mp3;
+  double b = y - db;
+  db = (y - b) - db;
+
+  *a = b;
+  *da = db;
+
+  return n;
+}
+
+/* Compute sin (A + DA).  cos can be computed by shifting the quadrant N
+   clockwise.  */
+static double
+__always_inline
+do_sincos_1 (double a, double da, double x, int4 n, int4 k)
+{
+  double xx, retval, res, cor, y;
+  mynumber u;
+  int m;
+  double eps = fabs (x) * 1.2e-30;
+
+  int k1 = (n + k) & 3;
+  switch (k1)
+    {			/* quarter of unit circle */
+    case 2:
+      a = -a;
+      da = -da;
+    case 0:
+      xx = a * a;
+      if (xx < 0.01588)
+	{
+	  /* Taylor series.  */
+	  res = TAYLOR_SIN (xx, a, da, cor);
+	  cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps;
+	  retval = (res == res + cor) ? res : sloww (a, da, x, k);
+	}
+      else
+	{
+	  if (a > 0)
+	    m = 1;
+	  else
+	    {
+	      m = 0;
+	      a = -a;
+	      da = -da;
+	    }
+	  u.x = big + a;
+	  y = a - (u.x - big);
+	  res = do_sin (u, y, da, &cor);
+	  cor = (cor > 0) ? 1.035 * cor + eps : 1.035 * cor - eps;
+	  retval = ((res == res + cor) ? ((m) ? res : -res)
+		    : sloww1 (a, da, x, m, k));
+	}
+      break;
+
+    case 1:
+    case 3:
+      if (a < 0)
+	{
+	  a = -a;
+	  da = -da;
+	}
+      u.x = big + a;
+      y = a - (u.x - big) + da;
+      res = do_cos (u, y, &cor);
+      cor = (cor > 0) ? 1.025 * cor + eps : 1.025 * cor - eps;
+      retval = ((res == res + cor) ? ((k1 & 2) ? -res : res)
+		: sloww2 (a, da, x, n));
+      break;
+    }
+
+  return retval;
+}
+
+static inline int4
+__always_inline
 reduce_sincos_2 (double x, double *a, double *da)
 {
   mynumber v;
@@ -386,9 +465,9 @@ SECTION
 #endif
 __sin (double x)
 {
-  double xx, res, t, cor, y, s, c, sn, ssn, cs, ccs, xn, a, da, eps;
-  mynumber u, v;
-  int4 k, m, n;
+  double xx, res, t, cor, y, s, c, sn, ssn, cs, ccs;
+  mynumber u;
+  int4 k, m;
   double retval = 0;
 
 #ifndef IN_SINCOS
@@ -452,73 +531,15 @@ __sin (double x)
       retval = (res == res + 1.020 * cor) ? ((m > 0) ? res : -res) : slow2 (x);
     }				/*   else  if (k < 0x400368fd)    */
 
+#ifndef IN_SINCOS
 /*-------------------------- 2.426265<|x|< 105414350 ----------------------*/
   else if (k < 0x419921FB)
     {
-      t = (x * hpinv + toint);
-      xn = t - toint;
-      v.x = t;
-      y = (x - xn * mp1) - xn * mp2;
-      n = v.i[LOW_HALF] & 3;
-      da = xn * mp3;
-      a = y - da;
-      da = (y - a) - da;
-      eps = fabs (x) * 1.2e-30;
-
-      switch (n)
-	{			/* quarter of unit circle */
-	case 0:
-	case 2:
-	  xx = a * a;
-	  if (n)
-	    {
-	      a = -a;
-	      da = -da;
-	    }
-	  if (xx < 0.01588)
-	    {
-	      /* Taylor series.  */
-	      res = TAYLOR_SIN (xx, a, da, cor);
-	      cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps;
-	      retval = (res == res + cor) ? res : sloww (a, da, x);
-	    }
-	  else
-	    {
-	      if (a > 0)
-		m = 1;
-	      else
-		{
-		  m = 0;
-		  a = -a;
-		  da = -da;
-		}
-	      u.x = big + a;
-	      y = a - (u.x - big);
-	      res = do_sin (u, y, da, &cor);
-	      cor = (cor > 0) ? 1.035 * cor + eps : 1.035 * cor - eps;
-	      retval = ((res == res + cor) ? ((m) ? res : -res)
-			: sloww1 (a, da, x, m));
-	    }
-	  break;
-
-	case 1:
-	case 3:
-	  if (a < 0)
-	    {
-	      a = -a;
-	      da = -da;
-	    }
-	  u.x = big + a;
-	  y = a - (u.x - big) + da;
-	  res = do_cos (u, y, &cor);
-	  cor = (cor > 0) ? 1.025 * cor + eps : 1.025 * cor - eps;
-	  retval = ((res == res + cor) ? ((n & 2) ? -res : res)
-		    : sloww2 (a, da, x, n));
-	  break;
-	}
+      double a, da;
+      int4 n = reduce_sincos_1 (x, &a, &da);
+      retval = do_sincos_1 (a, da, x, n, 0);
     }				/*   else  if (k <  0x419921FB )    */
 
-#ifndef IN_SINCOS
 /*---------------------105414350 <|x|< 281474976710656 --------------------*/
   else if (k < 0x42F00000)
     {
@@ -558,9 +579,9 @@ SECTION
 #endif
 __cos (double x)
 {
-  double y, xx, res, t, cor, xn, a, da, eps;
-  mynumber u, v;
-  int4 k, m, n;
+  double y, xx, res, cor, a, da;
+  mynumber u;
+  int4 k, m;
 
   double retval = 0;
 
@@ -595,7 +616,7 @@ __cos (double x)
 	{
 	  res = TAYLOR_SIN (xx, a, da, cor);
 	  cor = (cor > 0) ? 1.02 * cor + 1.0e-31 : 1.02 * cor - 1.0e-31;
-	  retval = (res == res + cor) ? res : csloww (a, da, x);
+	  retval = (res == res + cor) ? res : sloww (a, da, x, 1);
 	}
       else
 	{
@@ -614,79 +635,20 @@ __cos (double x)
 	  res = do_sin (u, y, da, &cor);
 	  cor = (cor > 0) ? 1.035 * cor + 1.0e-31 : 1.035 * cor - 1.0e-31;
 	  retval = ((res == res + cor) ? ((m) ? res : -res)
-		    : csloww1 (a, da, x, m));
+		    : sloww1 (a, da, x, m, 1));
 	}
 
     }				/*   else  if (k < 0x400368fd)    */
 
 
+#ifndef IN_SINCOS
   else if (k < 0x419921FB)
     {				/* 2.426265<|x|< 105414350 */
-      t = (x * hpinv + toint);
-      xn = t - toint;
-      v.x = t;
-      y = (x - xn * mp1) - xn * mp2;
-      n = v.i[LOW_HALF] & 3;
-      da = xn * mp3;
-      a = y - da;
-      da = (y - a) - da;
-      eps = fabs (x) * 1.2e-30;
-
-      switch (n)
-	{
-	case 1:
-	case 3:
-	  xx = a * a;
-	  if (n == 1)
-	    {
-	      a = -a;
-	      da = -da;
-	    }
-	  if (xx < 0.01588)
-	    {
-	      res = TAYLOR_SIN (xx, a, da, cor);
-	      cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps;
-	      retval = (res == res + cor) ? res : csloww (a, da, x);
-	    }
-	  else
-	    {
-	      if (a > 0)
-		{
-		  m = 1;
-		}
-	      else
-		{
-		  m = 0;
-		  a = -a;
-		  da = -da;
-		}
-	      u.x = big + a;
-	      y = a - (u.x - big);
-	      res = do_sin (u, y, da, &cor);
-	      cor = (cor > 0) ? 1.035 * cor + eps : 1.035 * cor - eps;
-	      retval = ((res == res + cor) ? ((m) ? res : -res)
-			: csloww1 (a, da, x, m));
-	    }
-	  break;
-
-	case 0:
-	case 2:
-	  if (a < 0)
-	    {
-	      a = -a;
-	      da = -da;
-	    }
-	  u.x = big + a;
-	  y = a - (u.x - big) + da;
-	  res = do_cos (u, y, &cor);
-	  cor = (cor > 0) ? 1.025 * cor + eps : 1.025 * cor - eps;
-	  retval = ((res == res + cor) ? ((n) ? -res : res)
-		    : csloww2 (a, da, x, n));
-	  break;
-	}
+      double a, da;
+      int4 n = reduce_sincos_1 (x, &a, &da);
+      retval = do_sincos_1 (a, da, x, n, 1);
     }				/*   else  if (k <  0x419921FB )    */
 
-#ifndef IN_SINCOS
   else if (k < 0x42F00000)
     {
       double a, da;
@@ -805,12 +767,13 @@ slow2 (double x)
 
 static double
 SECTION
-sloww (double x, double dx, double orig)
+sloww (double x, double dx, double orig, int k)
 {
   double y, t, res, cor, w[2], a, da, xn;
   mynumber v;
   int4 n;
   res = TAYLOR_SLOW (x, dx, cor);
+
   if (cor > 0)
     cor = 1.0005 * cor + fabs (orig) * 3.1e-30;
   else
@@ -832,14 +795,15 @@ sloww (double x, double dx, double orig)
   xn = t - toint;
   v.x = t;
   y = (orig - xn * mp1) - xn * mp2;
-  n = v.i[LOW_HALF] & 3;
+  n = (v.i[LOW_HALF] + k) & 3;
   da = xn * pp3;
   t = y - da;
   da = (y - t) - da;
   y = xn * pp4;
   a = t - y;
   da = ((t - a) - y) + da;
-  if (n & 2)
+
+  if (n == 2 || n == 1)
     {
       a = -a;
       da = -da;
@@ -853,7 +817,7 @@ sloww (double x, double dx, double orig)
   if (w[0] == w[0] + cor)
     return (a > 0) ? w[0] : -w[0];
 
-  return __mpsin (orig, 0, true);
+  return (n & 1) ? __mpcos (orig, 0, true) : __mpsin (orig, 0, true);
 }
 
 /***************************************************************************/
@@ -865,7 +829,7 @@ sloww (double x, double dx, double orig)
 
 static double
 SECTION
-sloww1 (double x, double dx, double orig, int m)
+sloww1 (double x, double dx, double orig, int m, int k)
 {
   mynumber u;
   double w[2], y, cor, res;
@@ -887,7 +851,7 @@ sloww1 (double x, double dx, double orig, int m)
   if (w[0] == w[0] + cor)
     return (m > 0) ? w[0] : -w[0];
 
-  return __mpsin (orig, 0, true);
+  return (k == 1) ? __mpcos (orig, 0, true) : __mpsin (orig, 0, true);
 }
 
 /***************************************************************************/
@@ -921,7 +885,7 @@ sloww2 (double x, double dx, double orig, int n)
   if (w[0] == w[0] + cor)
     return (n & 2) ? -w[0] : w[0];
 
-  return __mpsin (orig, 0, true);
+  return (n & 1) ? __mpsin (orig, 0, true) : __mpcos (orig, 0, true);
 }
 
 /***************************************************************************/
@@ -1052,138 +1016,6 @@ cslow2 (double x)
   return __mpcos (x, 0, false);
 }
 
-/***************************************************************************/
-/*  Routine compute cos(x+dx) (Double-Length number) where x is small enough*/
-/* to use Taylor series around zero and   (x+dx) .Routine receive also      */
-/* (right argument) the  original   value of x for computing error of      */
-/* result.And if result not accurate enough routine calls other routines    */
-/***************************************************************************/
-
-
-static double
-SECTION
-csloww (double x, double dx, double orig)
-{
-  double y, t, res, cor, w[2], a, da, xn;
-  mynumber v;
-  int4 n;
-
-  /* Taylor series */
-  res = TAYLOR_SLOW (x, dx, cor);
-
-  if (cor > 0)
-    cor = 1.0005 * cor + fabs (orig) * 3.1e-30;
-  else
-    cor = 1.0005 * cor - fabs (orig) * 3.1e-30;
-
-  if (res == res + cor)
-    return res;
-
-  (x > 0) ? __dubsin (x, dx, w) : __dubsin (-x, -dx, w);
-
-  if (w[1] > 0)
-    cor = 1.000000001 * w[1] + fabs (orig) * 1.1e-30;
-  else
-    cor = 1.000000001 * w[1] - fabs (orig) * 1.1e-30;
-
-  if (w[0] == w[0] + cor)
-    return (x > 0) ? w[0] : -w[0];
-
-  t = (orig * hpinv + toint);
-  xn = t - toint;
-  v.x = t;
-  y = (orig - xn * mp1) - xn * mp2;
-  n = v.i[LOW_HALF] & 3;
-  da = xn * pp3;
-  t = y - da;
-  da = (y - t) - da;
-  y = xn * pp4;
-  a = t - y;
-  da = ((t - a) - y) + da;
-  if (n == 1)
-    {
-      a = -a;
-      da = -da;
-    }
-  (a > 0) ? __dubsin (a, da, w) : __dubsin (-a, -da, w);
-
-  if (w[1] > 0)
-    cor = 1.000000001 * w[1] + fabs (orig) * 1.1e-40;
-  else
-    cor = 1.000000001 * w[1] - fabs (orig) * 1.1e-40;
-
-  if (w[0] == w[0] + cor)
-    return (a > 0) ? w[0] : -w[0];
-
-  return __mpcos (orig, 0, true);
-}
-
-/***************************************************************************/
-/*  Routine compute sin(x+dx)   (Double-Length number) where x in first or  */
-/*  third quarter of unit circle.Routine receive also (right argument) the  */
-/*  original   value of x for computing error of result.And if result not  */
-/* accurate enough routine calls  other routines                            */
-/***************************************************************************/
-
-static double
-SECTION
-csloww1 (double x, double dx, double orig, int m)
-{
-  mynumber u;
-  double w[2], y, cor, res;
-
-  u.x = big + x;
-  y = x - (u.x - big);
-  res = do_sin_slow (u, y, dx, 3.1e-30 * fabs (orig), &cor);
-
-  if (res == res + cor)
-    return (m > 0) ? res : -res;
-
-  __dubsin (x, dx, w);
-  if (w[1] > 0)
-    cor = 1.000000005 * w[1] + 1.1e-30 * fabs (orig);
-  else
-    cor = 1.000000005 * w[1] - 1.1e-30 * fabs (orig);
-
-  if (w[0] == w[0] + cor)
-    return (m > 0) ? w[0] : -w[0];
-
-  return __mpcos (orig, 0, true);
-}
-
-
-/***************************************************************************/
-/*  Routine compute sin(x+dx)   (Double-Length number) where x in second or */
-/*  fourth quarter of unit circle.Routine receive also  the  original value */
-/* and quarter(n= 1or 3)of x for computing error of result.And if result not*/
-/* accurate enough routine calls  other routines                            */
-/***************************************************************************/
-
-static double
-SECTION
-csloww2 (double x, double dx, double orig, int n)
-{
-  mynumber u;
-  double w[2], y, cor, res;
-
-  u.x = big + x;
-  y = x - (u.x - big);
-  res = do_cos_slow (u, y, dx, 3.1e-30 * fabs (orig), &cor);
-
-  if (res == res + cor)
-    return (n) ? -res : res;
-
-  __docos (x, dx, w);
-  if (w[1] > 0)
-    cor = 1.000000005 * w[1] + 1.1e-30 * fabs (orig);
-  else
-    cor = 1.000000005 * w[1] - 1.1e-30 * fabs (orig);
-  if (w[0] == w[0] + cor)
-    return (n) ? -w[0] : w[0];
-
-  return __mpcos (orig, 0, true);
-}
-
 #ifndef __cos
 weak_alias (__cos, cos)
 # ifdef NO_LONG_DOUBLE
diff --git a/sysdeps/ieee754/dbl-64/s_sincos.c b/sysdeps/ieee754/dbl-64/s_sincos.c
index 7f78b41..0290004 100644
--- a/sysdeps/ieee754/dbl-64/s_sincos.c
+++ b/sysdeps/ieee754/dbl-64/s_sincos.c
@@ -69,12 +69,22 @@ __sincos (double x, double *sinx, double *cosx)
   u.x = x;
   k = 0x7fffffff & u.i[HIGH_HALF];
 
-  if (k < 0x419921FB)
+  if (k < 0x400368fd)
     {
       *sinx = __sin_local (x);
       *cosx = __cos_local (x);
       return;
     }
+  if (k < 0x419921FB)
+    {
+      double a, da;
+      int4 n = reduce_sincos_1 (x, &a, &da);
+
+      *sinx = do_sincos_1 (a, da, x, n, 0);
+      *cosx = do_sincos_1 (a, da, x, n, 1);
+
+      return;
+    }
   if (k < 0x42F00000)
     {
       double a, da;

http://sourceware.org/git/gitweb.cgi?p=glibc.git;a=commitdiff;h=f7953c44d545c91442f123300288853284501dab

commit f7953c44d545c91442f123300288853284501dab
Author: Siddhesh Poyarekar <siddhesh.poyarekar@linaro.org>
Date:   Mon Dec 21 10:41:46 2015 +0530

    Consolidate sin and cos code for 105414350 <|x|< 281474976710656
    
    The sin and cos computation for this range of input is identical
    except for a difference in quadrants by 1.  Exploit that fact and the
    common argument reduction to reduce computations for sincos.

diff --git a/ChangeLog b/ChangeLog
index c8281dd..2228769 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,5 +1,10 @@
 2015-12-21  Siddhesh Poyarekar  <siddhesh.poyarekar@linaro.org>
 
+	* sysdeps/ieee754/dbl-64/s_sin.c (__sin, __cos): Move common
+	code to new functions.
+	(reduce_sincos_2, do_sincos_2): New functions.
+	* sysdeps/ieee754/dbl-64/s_sincos.c (__sincos): Use them.
+
 	* sysdeps/ieee754/dbl-64/s_sin.c (__sin) [!IN_SINCOS]: Skip
 	common code for sincos.
 	(__cos) [!IN_SINCOS]: Likewise.
diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c
index b619905..3b26a61 100644
--- a/sysdeps/ieee754/dbl-64/s_sin.c
+++ b/sysdeps/ieee754/dbl-64/s_sin.c
@@ -276,6 +276,104 @@ reduce_and_compute (double x, unsigned int k)
   return retval;
 }
 
+static inline int4
+__always_inline
+reduce_sincos_2 (double x, double *a, double *da)
+{
+  mynumber v;
+
+  double t = (x * hpinv + toint);
+  double xn = t - toint;
+  v.x = t;
+  double xn1 = (xn + 8.0e22) - 8.0e22;
+  double xn2 = xn - xn1;
+  double y = ((((x - xn1 * mp1) - xn1 * mp2) - xn2 * mp1) - xn2 * mp2);
+  int4 n = v.i[LOW_HALF] & 3;
+  double db = xn1 * pp3;
+  t = y - db;
+  db = (y - t) - db;
+  db = (db - xn2 * pp3) - xn * pp4;
+  double b = t + db;
+  db = (t - b) + db;
+
+  *a = b;
+  *da = db;
+
+  return n;
+}
+
+/* Compute sin (A + DA).  cos can be computed by shifting the quadrant N
+   clockwise.  */
+static double
+__always_inline
+do_sincos_2 (double a, double da, double x, int4 n, int4 k)
+{
+  double res, retval, cor, xx;
+  mynumber u;
+
+  double eps = 1.0e-24;
+
+  k = (n + k) & 3;
+
+  switch (k)
+    {
+    case 2:
+      a = -a;
+      da = -da;
+      /* Fall through.  */
+    case 0:
+      xx = a * a;
+      if (xx < 0.01588)
+	{
+	  /* Taylor series.  */
+	  res = TAYLOR_SIN (xx, a, da, cor);
+	  cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps;
+	  retval = (res == res + cor) ? res : bsloww (a, da, x, n);
+	}
+      else
+	{
+	  double t, db, y;
+	  int m;
+	  if (a > 0)
+	    {
+	      m = 1;
+	      t = a;
+	      db = da;
+	    }
+	  else
+	    {
+	      m = 0;
+	      t = -a;
+	      db = -da;
+	    }
+	  u.x = big + t;
+	  y = t - (u.x - big);
+	  res = do_sin (u, y, db, &cor);
+	  cor = (cor > 0) ? 1.035 * cor + eps : 1.035 * cor - eps;
+	  retval = ((res == res + cor) ? ((m) ? res : -res)
+		    : bsloww1 (a, da, x, n));
+	}
+      break;
+
+    case 1:
+    case 3:
+      if (a < 0)
+	{
+	  a = -a;
+	  da = -da;
+	}
+      u.x = big + a;
+      double y = a - (u.x - big) + da;
+      res = do_cos (u, y, &cor);
+      cor = (cor > 0) ? 1.025 * cor + eps : 1.025 * cor - eps;
+      retval = ((res == res + cor) ? ((n & 2) ? -res : res)
+		: bsloww2 (a, da, x, n));
+      break;
+    }
+
+  return retval;
+}
+
 /*******************************************************************/
 /* An ultimate sin routine. Given an IEEE double machine number x   */
 /* it computes the correctly rounded (to nearest) value of sin(x)  */
@@ -288,8 +386,7 @@ SECTION
 #endif
 __sin (double x)
 {
-  double xx, res, t, cor, y, s, c, sn, ssn, cs, ccs, xn, a, da, db, eps, xn1,
-    xn2;
+  double xx, res, t, cor, y, s, c, sn, ssn, cs, ccs, xn, a, da, eps;
   mynumber u, v;
   int4 k, m, n;
   double retval = 0;
@@ -421,83 +518,16 @@ __sin (double x)
 	}
     }				/*   else  if (k <  0x419921FB )    */
 
+#ifndef IN_SINCOS
 /*---------------------105414350 <|x|< 281474976710656 --------------------*/
   else if (k < 0x42F00000)
     {
-      t = (x * hpinv + toint);
-      xn = t - toint;
-      v.x = t;
-      xn1 = (xn + 8.0e22) - 8.0e22;
-      xn2 = xn - xn1;
-      y = ((((x - xn1 * mp1) - xn1 * mp2) - xn2 * mp1) - xn2 * mp2);
-      n = v.i[LOW_HALF] & 3;
-      da = xn1 * pp3;
-      t = y - da;
-      da = (y - t) - da;
-      da = (da - xn2 * pp3) - xn * pp4;
-      a = t + da;
-      da = (t - a) + da;
-      eps = 1.0e-24;
+      double a, da;
 
-      switch (n)
-	{
-	case 0:
-	case 2:
-	  xx = a * a;
-	  if (n)
-	    {
-	      a = -a;
-	      da = -da;
-	    }
-	  if (xx < 0.01588)
-	    {
-	      /* Taylor series.  */
-	      res = TAYLOR_SIN (xx, a, da, cor);
-	      cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps;
-	      retval = (res == res + cor) ? res : bsloww (a, da, x, n);
-	    }
-	  else
-	    {
-	      double t;
-	      if (a > 0)
-		{
-		  m = 1;
-		  t = a;
-		  db = da;
-		}
-	      else
-		{
-		  m = 0;
-		  t = -a;
-		  db = -da;
-		}
-	      u.x = big + t;
-	      y = t - (u.x - big);
-	      res = do_sin (u, y, db, &cor);
-	      cor = (cor > 0) ? 1.035 * cor + eps : 1.035 * cor - eps;
-	      retval = ((res == res + cor) ? ((m) ? res : -res)
-			: bsloww1 (a, da, x, n));
-	    }
-	  break;
-
-	case 1:
-	case 3:
-	  if (a < 0)
-	    {
-	      a = -a;
-	      da = -da;
-	    }
-	  u.x = big + a;
-	  y = a - (u.x - big) + da;
-	  res = do_cos (u, y, &cor);
-	  cor = (cor > 0) ? 1.025 * cor + eps : 1.025 * cor - eps;
-	  retval = ((res == res + cor) ? ((n & 2) ? -res : res)
-		    : bsloww2 (a, da, x, n));
-	  break;
-	}
+      int4 n = reduce_sincos_2 (x, &a, &da);
+      retval = do_sincos_2 (a, da, x, n, 0);
     }				/*   else  if (k <  0x42F00000 )   */
 
-#ifndef IN_SINCOS
 /* -----------------281474976710656 <|x| <2^1024----------------------------*/
   else if (k < 0x7ff00000)
     retval = reduce_and_compute (x, 0);
@@ -528,8 +558,7 @@ SECTION
 #endif
 __cos (double x)
 {
-  double y, xx, res, t, cor, xn, a, da, db, eps, xn1,
-    xn2;
+  double y, xx, res, t, cor, xn, a, da, eps;
   mynumber u, v;
   int4 k, m, n;
 
@@ -657,81 +686,15 @@ __cos (double x)
 	}
     }				/*   else  if (k <  0x419921FB )    */
 
+#ifndef IN_SINCOS
   else if (k < 0x42F00000)
     {
-      t = (x * hpinv + toint);
-      xn = t - toint;
-      v.x = t;
-      xn1 = (xn + 8.0e22) - 8.0e22;
-      xn2 = xn - xn1;
-      y = ((((x - xn1 * mp1) - xn1 * mp2) - xn2 * mp1) - xn2 * mp2);
-      n = v.i[LOW_HALF] & 3;
-      da = xn1 * pp3;
-      t = y - da;
-      da = (y - t) - da;
-      da = (da - xn2 * pp3) - xn * pp4;
-      a = t + da;
-      da = (t - a) + da;
-      eps = 1.0e-24;
-
-      switch (n)
-	{
-	case 1:
-	case 3:
-	  xx = a * a;
-	  if (n == 1)
-	    {
-	      a = -a;
-	      da = -da;
-	    }
-	  if (xx < 0.01588)
-	    {
-	      res = TAYLOR_SIN (xx, a, da, cor);
-	      cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps;
-	      retval = (res == res + cor) ? res : bsloww (a, da, x, n);
-	    }
-	  else
-	    {
-	      double t;
-	      if (a > 0)
-		{
-		  m = 1;
-		  t = a;
-		  db = da;
-		}
-	      else
-		{
-		  m = 0;
-		  t = -a;
-		  db = -da;
-		}
-	      u.x = big + t;
-	      y = t - (u.x - big);
-	      res = do_sin (u, y, db, &cor);
-	      cor = (cor > 0) ? 1.035 * cor + eps : 1.035 * cor - eps;
-	      retval = ((res == res + cor) ? ((m) ? res : -res)
-			: bsloww1 (a, da, x, n));
-	    }
-	  break;
+      double a, da;
 
-	case 0:
-	case 2:
-	  if (a < 0)
-	    {
-	      a = -a;
-	      da = -da;
-	    }
-	  u.x = big + a;
-	  y = a - (u.x - big) + da;
-	  res = do_cos (u, y, &cor);
-	  cor = (cor > 0) ? 1.025 * cor + eps : 1.025 * cor - eps;
-	  retval = ((res == res + cor) ? ((n) ? -res : res)
-		    : bsloww2 (a, da, x, n));
-	  break;
-	}
+      int4 n = reduce_sincos_2 (x, &a, &da);
+      retval = do_sincos_2 (a, da, x, n, 1);
     }				/*   else  if (k <  0x42F00000 )    */
 
-#ifndef IN_SINCOS
   /* 281474976710656 <|x| <2^1024 */
   else if (k < 0x7ff00000)
     retval = reduce_and_compute (x, 1);
diff --git a/sysdeps/ieee754/dbl-64/s_sincos.c b/sysdeps/ieee754/dbl-64/s_sincos.c
index f50ffa6..7f78b41 100644
--- a/sysdeps/ieee754/dbl-64/s_sincos.c
+++ b/sysdeps/ieee754/dbl-64/s_sincos.c
@@ -69,12 +69,22 @@ __sincos (double x, double *sinx, double *cosx)
   u.x = x;
   k = 0x7fffffff & u.i[HIGH_HALF];
 
-  if (k < 0x42F00000)
+  if (k < 0x419921FB)
     {
       *sinx = __sin_local (x);
       *cosx = __cos_local (x);
       return;
     }
+  if (k < 0x42F00000)
+    {
+      double a, da;
+      int4 n = reduce_sincos_2 (x, &a, &da);
+
+      *sinx = do_sincos_2 (a, da, x, n, 0);
+      *cosx = do_sincos_2 (a, da, x, n, 1);
+
+      return;
+    }
   if (k < 0x7ff00000)
     {
       reduce_and_compute_sincos (x, sinx, cosx);

http://sourceware.org/git/gitweb.cgi?p=glibc.git;a=commitdiff;h=a045832debe4191bcffed6e539a5470152e9c9e8

commit a045832debe4191bcffed6e539a5470152e9c9e8
Author: Siddhesh Poyarekar <siddhesh.poyarekar@linaro.org>
Date:   Mon Dec 21 10:40:32 2015 +0530

    Consolidate range reduction in sincos for x > 281474976710656
    
    Range reduction needs to be done only once for sin and cos, so copy
    over all of the relevant functions (__sin, __cos, reduce_and_compute)
    and consolidate common code.

diff --git a/ChangeLog b/ChangeLog
index 0024c1e..c8281dd 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,12 @@
+2015-12-21  Siddhesh Poyarekar  <siddhesh.poyarekar@linaro.org>
+
+	* sysdeps/ieee754/dbl-64/s_sin.c (__sin) [!IN_SINCOS]: Skip
+	common code for sincos.
+	(__cos) [!IN_SINCOS]: Likewise.
+	* sysdeps/ieee754/dbl-64/s_sincos.c (reduce_and_compute_sincos):
+	New function.
+	(__sincos): Use it.
+
 2015-12-20  Aurelien Jarno  <aurelien@aurel32.net>
 
 	* sysdeps/i386/fpu/libm-test-ulps: Move to ....
diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c
index a635a86..b619905 100644
--- a/sysdeps/ieee754/dbl-64/s_sin.c
+++ b/sysdeps/ieee754/dbl-64/s_sin.c
@@ -497,6 +497,7 @@ __sin (double x)
 	}
     }				/*   else  if (k <  0x42F00000 )   */
 
+#ifndef IN_SINCOS
 /* -----------------281474976710656 <|x| <2^1024----------------------------*/
   else if (k < 0x7ff00000)
     retval = reduce_and_compute (x, 0);
@@ -508,6 +509,7 @@ __sin (double x)
 	__set_errno (EDOM);
       retval = x / x;
     }
+#endif
 
   return retval;
 }
@@ -729,6 +731,7 @@ __cos (double x)
 	}
     }				/*   else  if (k <  0x42F00000 )    */
 
+#ifndef IN_SINCOS
   /* 281474976710656 <|x| <2^1024 */
   else if (k < 0x7ff00000)
     retval = reduce_and_compute (x, 1);
@@ -739,6 +742,7 @@ __cos (double x)
 	__set_errno (EDOM);
       retval = x / x;		/* |x| > 2^1024 */
     }
+#endif
 
   return retval;
 }
diff --git a/sysdeps/ieee754/dbl-64/s_sincos.c b/sysdeps/ieee754/dbl-64/s_sincos.c
index 2a3fc06..f50ffa6 100644
--- a/sysdeps/ieee754/dbl-64/s_sincos.c
+++ b/sysdeps/ieee754/dbl-64/s_sincos.c
@@ -27,13 +27,64 @@
 #define IN_SINCOS 1
 #include "s_sin.c"
 
+/* Consolidated version of reduce_and_compute in s_sin.c that does range
+   reduction only once and computes sin and cos together.  */
+static inline void
+__always_inline
+reduce_and_compute_sincos (double x, double *sinx, double *cosx)
+{
+  double a, da;
+  unsigned int n = __branred (x, &a, &da);
+
+  n = n & 3;
+
+  if (n == 1 || n == 2)
+    {
+      a = -a;
+      da = -da;
+    }
+
+  if (n & 1)
+    {
+      double *temp = cosx;
+      cosx = sinx;
+      sinx = temp;
+    }
+
+  if (a * a < 0.01588)
+    *sinx = bsloww (a, da, x, n);
+  else
+    *sinx = bsloww1 (a, da, x, n);
+  *cosx = bsloww2 (a, da, x, n);
+}
+
 void
 __sincos (double x, double *sinx, double *cosx)
 {
+  mynumber u;
+  int k;
+
   SET_RESTORE_ROUND_53BIT (FE_TONEAREST);
 
-  *sinx = __sin (x);
-  *cosx = __cos (x);
+  u.x = x;
+  k = 0x7fffffff & u.i[HIGH_HALF];
+
+  if (k < 0x42F00000)
+    {
+      *sinx = __sin_local (x);
+      *cosx = __cos_local (x);
+      return;
+    }
+  if (k < 0x7ff00000)
+    {
+      reduce_and_compute_sincos (x, sinx, cosx);
+      return;
+    }
+
+  if (isinf (x))
+    __set_errno (EDOM);
+
+  *sinx = *cosx = x / x;
 }
 weak_alias (__sincos, sincos)
 #ifdef NO_LONG_DOUBLE

http://sourceware.org/git/gitweb.cgi?p=glibc.git;a=commitdiff;h=760c2eb7da879e05252468a512a582a0e5c180af

commit 760c2eb7da879e05252468a512a582a0e5c180af
Author: Siddhesh Poyarekar <siddhesh.poyarekar@linaro.org>
Date:   Mon Dec 21 10:38:38 2015 +0530

    Fix up ChangeLog

diff --git a/ChangeLog b/ChangeLog
index f226ce8..0024c1e 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,8 +1,8 @@
 2015-12-20  Aurelien Jarno  <aurelien@aurel32.net>
 
-       * sysdeps/i386/fpu/libm-test-ulps: Move to ....
-       * sysdeps/i386/i686/multiarch/fpu/libm-test-ulps: ...here.
-       * sysdeps/i386/fpu/libm-test-ulps: Regenerate.
+	* sysdeps/i386/fpu/libm-test-ulps: Move to ....
+	* sysdeps/i386/i686/multiarch/fpu/libm-test-ulps: ...here.
+	* sysdeps/i386/fpu/libm-test-ulps: Regenerate.
 
 2015-12-19  Andrew Senkevich  <andrew.senkevich@intel.com>
 

-----------------------------------------------------------------------

Summary of changes:
 ChangeLog                         |   29 ++-
 sysdeps/ieee754/dbl-64/s_sin.c    |  637 +++++++++++++------------------------
 sysdeps/ieee754/dbl-64/s_sincos.c |   75 +++++-
 3 files changed, 317 insertions(+), 424 deletions(-)


hooks/post-receive
-- 
GNU C Library master sources



More information about the Glibc-cvs mailing list