public inbox for glibc-cvs@sourceware.org
help / color / mirror / Atom feed
* [glibc/maskray/grte] [PATCH 6/7] sin/cos slow paths: refactor duplicated code into dosin
@ 2021-08-28  0:34 Fangrui Song
  0 siblings, 0 replies; 4+ messages in thread
From: Fangrui Song @ 2021-08-28  0:34 UTC (permalink / raw)
  To: glibc-cvs

https://sourceware.org/git/gitweb.cgi?p=glibc.git;h=c8aaaf67f69d6a9dc92e8fe18dc5c12362be76ea

commit c8aaaf67f69d6a9dc92e8fe18dc5c12362be76ea
Author: Wilco Dijkstra <wdijkstr@arm.com>
Date:   Tue Apr 3 16:46:10 2018 +0100

    [PATCH 6/7] sin/cos slow paths: refactor duplicated code into dosin
    
    Refactor duplicated code into do_sin.  Since all calls to do_sin use copysign to
    set the sign of the result, move it inside do_sin.  Small inputs use a separate
    polynomial, so move this into do_sin as well (the check is based on the more
    conservative case when doing large range reduction, but could be relaxed).
    
            * sysdeps/ieee754/dbl-64/s_sin.c (do_sin): Use TAYLOR_SIN for small
            inputs.  Return correct sign.
            (do_sincos): Remove small input check before do_sin, let do_sin set
            the sign.
            (__sin): Likewise.
            (__cos): Likewise.

Diff:
---
 sysdeps/ieee754/dbl-64/s_sin.c | 40 +++++++++++++---------------------------
 1 file changed, 13 insertions(+), 27 deletions(-)

diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c
index fcb2e6b83d..6e261f24bb 100644
--- a/sysdeps/ieee754/dbl-64/s_sin.c
+++ b/sysdeps/ieee754/dbl-64/s_sin.c
@@ -124,6 +124,11 @@ static inline double
 __always_inline
 do_sin (double x, double dx)
 {
+  double xold = x;
+  /* Max ULP is 0.501 if |x| < 0.126, otherwise ULP is 0.518.  */
+  if (fabs (x) < 0.126)
+    return TAYLOR_SIN (x * x, x, dx);
+
   mynumber u;
 
   if (x <= 0)
@@ -137,7 +142,7 @@ do_sin (double x, double dx)
   c = x * dx + xx * (cs2 + xx * (cs4 + xx * cs6));
   SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs);
   cor = (ssn + s * ccs - sn * c) + cs * s;
-  return sn + cor;
+  return __copysign (sn + cor, xold);
 }
 
 /* Reduce range of x to within PI/2 with abs (x) < 105414350.  The high part
@@ -181,14 +186,8 @@ do_sincos (double a, double da, int4 n)
     /* Max ULP is 0.513.  */
     retval = do_cos (a, da);
   else
-    {
-      double xx = a * a;
-      /* Max ULP is 0.501 if xx < 0.01588, otherwise ULP is 0.518.  */
-      if (xx < 0.01588)
-	retval = TAYLOR_SIN (xx, a, da);
-      else
-	retval = __copysign (do_sin (a, da), a);
-    }
+    /* Max ULP is 0.501 if xx < 0.01588, otherwise ULP is 0.518.  */
+    retval = do_sin (a, da);
 
   return (n & 2) ? -retval : retval;
 }
@@ -207,7 +206,7 @@ SECTION
 __sin (double x)
 {
 #ifndef IN_SINCOS
-  double xx, t, a, da;
+  double t, a, da;
   mynumber u;
   int4 k, m, n;
   double retval = 0;
@@ -228,20 +227,11 @@ __sin (double x)
       math_check_force_underflow (x);
       retval = x;
     }
- /*---------------------------- 2^-26 < |x|< 0.25 ----------------------*/
-  else if (k < 0x3fd00000)
-    {
-      xx = x * x;
-      /* Taylor series.  */
-      t = POLYNOMIAL (xx) * (xx * x);
-      /* Max ULP of x + t is 0.535.  */
-      retval = x + t;
-    }				/*  else  if (k < 0x3fd00000)    */
-/*---------------------------- 0.25<|x|< 0.855469---------------------- */
+/*--------------------------- 2^-26<|x|< 0.855469---------------------- */
   else if (k < 0x3feb6000)
     {
       /* Max ULP is 0.548.  */
-      retval = __copysign (do_sin (x, 0), x);
+      retval = do_sin (x, 0);
     }				/*   else  if (k < 0x3feb6000)    */
 
 /*----------------------- 0.855469  <|x|<2.426265  ----------------------*/
@@ -292,7 +282,7 @@ SECTION
 #endif
 __cos (double x)
 {
-  double y, xx, a, da;
+  double y, a, da;
   mynumber u;
 #ifndef IN_SINCOS
   int4 k, m, n;
@@ -325,13 +315,9 @@ __cos (double x)
       y = hp0 - fabs (x);
       a = y + hp1;
       da = (y - a) + hp1;
-      xx = a * a;
       /* Max ULP is 0.501 if xx < 0.01588 or 0.518 otherwise.
 	 Range reduction uses 106 bits here which is sufficient.  */
-      if (xx < 0.01588)
-	retval = TAYLOR_SIN (xx, a, da);
-      else
-	retval = __copysign (do_sin (a, da), a);
+      retval = do_sin (a, da);
     }				/*   else  if (k < 0x400368fd)    */


^ permalink raw reply	[flat|nested] 4+ messages in thread

* [glibc/maskray/grte] [PATCH 6/7] sin/cos slow paths: refactor duplicated code into dosin
@ 2021-08-28  0:30 Fangrui Song
  0 siblings, 0 replies; 4+ messages in thread
From: Fangrui Song @ 2021-08-28  0:30 UTC (permalink / raw)
  To: glibc-cvs

https://sourceware.org/git/gitweb.cgi?p=glibc.git;h=212cc76402afb2d0494ad902d8956017cc318094

commit 212cc76402afb2d0494ad902d8956017cc318094
Author: Wilco Dijkstra <wdijkstr@arm.com>
Date:   Tue Apr 3 16:46:10 2018 +0100

    [PATCH 6/7] sin/cos slow paths: refactor duplicated code into dosin
    
    Refactor duplicated code into do_sin.  Since all calls to do_sin use copysign to
    set the sign of the result, move it inside do_sin.  Small inputs use a separate
    polynomial, so move this into do_sin as well (the check is based on the more
    conservative case when doing large range reduction, but could be relaxed).
    
            * sysdeps/ieee754/dbl-64/s_sin.c (do_sin): Use TAYLOR_SIN for small
            inputs.  Return correct sign.
            (do_sincos): Remove small input check before do_sin, let do_sin set
            the sign.
            (__sin): Likewise.
            (__cos): Likewise.

Diff:
---
 sysdeps/ieee754/dbl-64/s_sin.c | 40 +++++++++++++---------------------------
 1 file changed, 13 insertions(+), 27 deletions(-)

diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c
index fcb2e6b83d..6e261f24bb 100644
--- a/sysdeps/ieee754/dbl-64/s_sin.c
+++ b/sysdeps/ieee754/dbl-64/s_sin.c
@@ -124,6 +124,11 @@ static inline double
 __always_inline
 do_sin (double x, double dx)
 {
+  double xold = x;
+  /* Max ULP is 0.501 if |x| < 0.126, otherwise ULP is 0.518.  */
+  if (fabs (x) < 0.126)
+    return TAYLOR_SIN (x * x, x, dx);
+
   mynumber u;
 
   if (x <= 0)
@@ -137,7 +142,7 @@ do_sin (double x, double dx)
   c = x * dx + xx * (cs2 + xx * (cs4 + xx * cs6));
   SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs);
   cor = (ssn + s * ccs - sn * c) + cs * s;
-  return sn + cor;
+  return __copysign (sn + cor, xold);
 }
 
 /* Reduce range of x to within PI/2 with abs (x) < 105414350.  The high part
@@ -181,14 +186,8 @@ do_sincos (double a, double da, int4 n)
     /* Max ULP is 0.513.  */
     retval = do_cos (a, da);
   else
-    {
-      double xx = a * a;
-      /* Max ULP is 0.501 if xx < 0.01588, otherwise ULP is 0.518.  */
-      if (xx < 0.01588)
-	retval = TAYLOR_SIN (xx, a, da);
-      else
-	retval = __copysign (do_sin (a, da), a);
-    }
+    /* Max ULP is 0.501 if xx < 0.01588, otherwise ULP is 0.518.  */
+    retval = do_sin (a, da);
 
   return (n & 2) ? -retval : retval;
 }
@@ -207,7 +206,7 @@ SECTION
 __sin (double x)
 {
 #ifndef IN_SINCOS
-  double xx, t, a, da;
+  double t, a, da;
   mynumber u;
   int4 k, m, n;
   double retval = 0;
@@ -228,20 +227,11 @@ __sin (double x)
       math_check_force_underflow (x);
       retval = x;
     }
- /*---------------------------- 2^-26 < |x|< 0.25 ----------------------*/
-  else if (k < 0x3fd00000)
-    {
-      xx = x * x;
-      /* Taylor series.  */
-      t = POLYNOMIAL (xx) * (xx * x);
-      /* Max ULP of x + t is 0.535.  */
-      retval = x + t;
-    }				/*  else  if (k < 0x3fd00000)    */
-/*---------------------------- 0.25<|x|< 0.855469---------------------- */
+/*--------------------------- 2^-26<|x|< 0.855469---------------------- */
   else if (k < 0x3feb6000)
     {
       /* Max ULP is 0.548.  */
-      retval = __copysign (do_sin (x, 0), x);
+      retval = do_sin (x, 0);
     }				/*   else  if (k < 0x3feb6000)    */
 
 /*----------------------- 0.855469  <|x|<2.426265  ----------------------*/
@@ -292,7 +282,7 @@ SECTION
 #endif
 __cos (double x)
 {
-  double y, xx, a, da;
+  double y, a, da;
   mynumber u;
 #ifndef IN_SINCOS
   int4 k, m, n;
@@ -325,13 +315,9 @@ __cos (double x)
       y = hp0 - fabs (x);
       a = y + hp1;
       da = (y - a) + hp1;
-      xx = a * a;
       /* Max ULP is 0.501 if xx < 0.01588 or 0.518 otherwise.
 	 Range reduction uses 106 bits here which is sufficient.  */
-      if (xx < 0.01588)
-	retval = TAYLOR_SIN (xx, a, da);
-      else
-	retval = __copysign (do_sin (a, da), a);
+      retval = do_sin (a, da);
     }				/*   else  if (k < 0x400368fd)    */


^ permalink raw reply	[flat|nested] 4+ messages in thread

* [glibc/maskray/grte] [PATCH 6/7] sin/cos slow paths: refactor duplicated code into dosin
@ 2021-08-27 23:48 Fangrui Song
  0 siblings, 0 replies; 4+ messages in thread
From: Fangrui Song @ 2021-08-27 23:48 UTC (permalink / raw)
  To: glibc-cvs

https://sourceware.org/git/gitweb.cgi?p=glibc.git;h=717749cd2765773a6e13f150bf95271e8fb36e3f

commit 717749cd2765773a6e13f150bf95271e8fb36e3f
Author: Wilco Dijkstra <wdijkstr@arm.com>
Date:   Tue Apr 3 16:46:10 2018 +0100

    [PATCH 6/7] sin/cos slow paths: refactor duplicated code into dosin
    
    Refactor duplicated code into do_sin.  Since all calls to do_sin use copysign to
    set the sign of the result, move it inside do_sin.  Small inputs use a separate
    polynomial, so move this into do_sin as well (the check is based on the more
    conservative case when doing large range reduction, but could be relaxed).
    
            * sysdeps/ieee754/dbl-64/s_sin.c (do_sin): Use TAYLOR_SIN for small
            inputs.  Return correct sign.
            (do_sincos): Remove small input check before do_sin, let do_sin set
            the sign.
            (__sin): Likewise.
            (__cos): Likewise.

Diff:
---
 sysdeps/ieee754/dbl-64/s_sin.c | 40 +++++++++++++---------------------------
 1 file changed, 13 insertions(+), 27 deletions(-)

diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c
index fcb2e6b83d..6e261f24bb 100644
--- a/sysdeps/ieee754/dbl-64/s_sin.c
+++ b/sysdeps/ieee754/dbl-64/s_sin.c
@@ -124,6 +124,11 @@ static inline double
 __always_inline
 do_sin (double x, double dx)
 {
+  double xold = x;
+  /* Max ULP is 0.501 if |x| < 0.126, otherwise ULP is 0.518.  */
+  if (fabs (x) < 0.126)
+    return TAYLOR_SIN (x * x, x, dx);
+
   mynumber u;
 
   if (x <= 0)
@@ -137,7 +142,7 @@ do_sin (double x, double dx)
   c = x * dx + xx * (cs2 + xx * (cs4 + xx * cs6));
   SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs);
   cor = (ssn + s * ccs - sn * c) + cs * s;
-  return sn + cor;
+  return __copysign (sn + cor, xold);
 }
 
 /* Reduce range of x to within PI/2 with abs (x) < 105414350.  The high part
@@ -181,14 +186,8 @@ do_sincos (double a, double da, int4 n)
     /* Max ULP is 0.513.  */
     retval = do_cos (a, da);
   else
-    {
-      double xx = a * a;
-      /* Max ULP is 0.501 if xx < 0.01588, otherwise ULP is 0.518.  */
-      if (xx < 0.01588)
-	retval = TAYLOR_SIN (xx, a, da);
-      else
-	retval = __copysign (do_sin (a, da), a);
-    }
+    /* Max ULP is 0.501 if xx < 0.01588, otherwise ULP is 0.518.  */
+    retval = do_sin (a, da);
 
   return (n & 2) ? -retval : retval;
 }
@@ -207,7 +206,7 @@ SECTION
 __sin (double x)
 {
 #ifndef IN_SINCOS
-  double xx, t, a, da;
+  double t, a, da;
   mynumber u;
   int4 k, m, n;
   double retval = 0;
@@ -228,20 +227,11 @@ __sin (double x)
       math_check_force_underflow (x);
       retval = x;
     }
- /*---------------------------- 2^-26 < |x|< 0.25 ----------------------*/
-  else if (k < 0x3fd00000)
-    {
-      xx = x * x;
-      /* Taylor series.  */
-      t = POLYNOMIAL (xx) * (xx * x);
-      /* Max ULP of x + t is 0.535.  */
-      retval = x + t;
-    }				/*  else  if (k < 0x3fd00000)    */
-/*---------------------------- 0.25<|x|< 0.855469---------------------- */
+/*--------------------------- 2^-26<|x|< 0.855469---------------------- */
   else if (k < 0x3feb6000)
     {
       /* Max ULP is 0.548.  */
-      retval = __copysign (do_sin (x, 0), x);
+      retval = do_sin (x, 0);
     }				/*   else  if (k < 0x3feb6000)    */
 
 /*----------------------- 0.855469  <|x|<2.426265  ----------------------*/
@@ -292,7 +282,7 @@ SECTION
 #endif
 __cos (double x)
 {
-  double y, xx, a, da;
+  double y, a, da;
   mynumber u;
 #ifndef IN_SINCOS
   int4 k, m, n;
@@ -325,13 +315,9 @@ __cos (double x)
       y = hp0 - fabs (x);
       a = y + hp1;
       da = (y - a) + hp1;
-      xx = a * a;
       /* Max ULP is 0.501 if xx < 0.01588 or 0.518 otherwise.
 	 Range reduction uses 106 bits here which is sufficient.  */
-      if (xx < 0.01588)
-	retval = TAYLOR_SIN (xx, a, da);
-      else
-	retval = __copysign (do_sin (a, da), a);
+      retval = do_sin (a, da);
     }				/*   else  if (k < 0x400368fd)    */


^ permalink raw reply	[flat|nested] 4+ messages in thread

* [glibc/maskray/grte] [PATCH 6/7] sin/cos slow paths: refactor duplicated code into dosin
@ 2021-08-27 23:47 Fangrui Song
  0 siblings, 0 replies; 4+ messages in thread
From: Fangrui Song @ 2021-08-27 23:47 UTC (permalink / raw)
  To: glibc-cvs

https://sourceware.org/git/gitweb.cgi?p=glibc.git;h=a94b663178de5b4bc4bf37c391c1b5f8e279c53a

commit a94b663178de5b4bc4bf37c391c1b5f8e279c53a
Author: Wilco Dijkstra <wdijkstr@arm.com>
Date:   Tue Apr 3 16:46:10 2018 +0100

    [PATCH 6/7] sin/cos slow paths: refactor duplicated code into dosin
    
    Refactor duplicated code into do_sin.  Since all calls to do_sin use copysign to
    set the sign of the result, move it inside do_sin.  Small inputs use a separate
    polynomial, so move this into do_sin as well (the check is based on the more
    conservative case when doing large range reduction, but could be relaxed).
    
            * sysdeps/ieee754/dbl-64/s_sin.c (do_sin): Use TAYLOR_SIN for small
            inputs.  Return correct sign.
            (do_sincos): Remove small input check before do_sin, let do_sin set
            the sign.
            (__sin): Likewise.
            (__cos): Likewise.

Diff:
---
 sysdeps/ieee754/dbl-64/s_sin.c | 40 +++++++++++++---------------------------
 1 file changed, 13 insertions(+), 27 deletions(-)

diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c
index fcb2e6b83d..6e261f24bb 100644
--- a/sysdeps/ieee754/dbl-64/s_sin.c
+++ b/sysdeps/ieee754/dbl-64/s_sin.c
@@ -124,6 +124,11 @@ static inline double
 __always_inline
 do_sin (double x, double dx)
 {
+  double xold = x;
+  /* Max ULP is 0.501 if |x| < 0.126, otherwise ULP is 0.518.  */
+  if (fabs (x) < 0.126)
+    return TAYLOR_SIN (x * x, x, dx);
+
   mynumber u;
 
   if (x <= 0)
@@ -137,7 +142,7 @@ do_sin (double x, double dx)
   c = x * dx + xx * (cs2 + xx * (cs4 + xx * cs6));
   SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs);
   cor = (ssn + s * ccs - sn * c) + cs * s;
-  return sn + cor;
+  return __copysign (sn + cor, xold);
 }
 
 /* Reduce range of x to within PI/2 with abs (x) < 105414350.  The high part
@@ -181,14 +186,8 @@ do_sincos (double a, double da, int4 n)
     /* Max ULP is 0.513.  */
     retval = do_cos (a, da);
   else
-    {
-      double xx = a * a;
-      /* Max ULP is 0.501 if xx < 0.01588, otherwise ULP is 0.518.  */
-      if (xx < 0.01588)
-	retval = TAYLOR_SIN (xx, a, da);
-      else
-	retval = __copysign (do_sin (a, da), a);
-    }
+    /* Max ULP is 0.501 if xx < 0.01588, otherwise ULP is 0.518.  */
+    retval = do_sin (a, da);
 
   return (n & 2) ? -retval : retval;
 }
@@ -207,7 +206,7 @@ SECTION
 __sin (double x)
 {
 #ifndef IN_SINCOS
-  double xx, t, a, da;
+  double t, a, da;
   mynumber u;
   int4 k, m, n;
   double retval = 0;
@@ -228,20 +227,11 @@ __sin (double x)
       math_check_force_underflow (x);
       retval = x;
     }
- /*---------------------------- 2^-26 < |x|< 0.25 ----------------------*/
-  else if (k < 0x3fd00000)
-    {
-      xx = x * x;
-      /* Taylor series.  */
-      t = POLYNOMIAL (xx) * (xx * x);
-      /* Max ULP of x + t is 0.535.  */
-      retval = x + t;
-    }				/*  else  if (k < 0x3fd00000)    */
-/*---------------------------- 0.25<|x|< 0.855469---------------------- */
+/*--------------------------- 2^-26<|x|< 0.855469---------------------- */
   else if (k < 0x3feb6000)
     {
       /* Max ULP is 0.548.  */
-      retval = __copysign (do_sin (x, 0), x);
+      retval = do_sin (x, 0);
     }				/*   else  if (k < 0x3feb6000)    */
 
 /*----------------------- 0.855469  <|x|<2.426265  ----------------------*/
@@ -292,7 +282,7 @@ SECTION
 #endif
 __cos (double x)
 {
-  double y, xx, a, da;
+  double y, a, da;
   mynumber u;
 #ifndef IN_SINCOS
   int4 k, m, n;
@@ -325,13 +315,9 @@ __cos (double x)
       y = hp0 - fabs (x);
       a = y + hp1;
       da = (y - a) + hp1;
-      xx = a * a;
       /* Max ULP is 0.501 if xx < 0.01588 or 0.518 otherwise.
 	 Range reduction uses 106 bits here which is sufficient.  */
-      if (xx < 0.01588)
-	retval = TAYLOR_SIN (xx, a, da);
-      else
-	retval = __copysign (do_sin (a, da), a);
+      retval = do_sin (a, da);
     }				/*   else  if (k < 0x400368fd)    */


^ permalink raw reply	[flat|nested] 4+ messages in thread

end of thread, other threads:[~2021-08-28  0:34 UTC | newest]

Thread overview: 4+ messages (download: mbox.gz / follow: Atom feed)
-- links below jump to the message on this page --
2021-08-28  0:34 [glibc/maskray/grte] [PATCH 6/7] sin/cos slow paths: refactor duplicated code into dosin Fangrui Song
  -- strict thread matches above, loose matches on Subject: below --
2021-08-28  0:30 Fangrui Song
2021-08-27 23:48 Fangrui Song
2021-08-27 23:47 Fangrui Song

This is a public inbox, see mirroring instructions
for how to clone and mirror all data and code used for this inbox;
as well as URLs for read-only IMAP folder(s) and NNTP newsgroup(s).