| Message ID | DB6PR0801MB2053407B5FE9C9DF2174313683DE0@DB6PR0801MB2053.eurprd08.prod.outlook.com |
|---|---|
| State | New |
| Headers | show |
| Series | [1/6] Remove slow paths from sin/cos | expand |
On Fri, 2018-03-09 at 15:50 +0000, Wilco Dijkstra wrote: > For huge inputs use the improved do_sincos function as well. Now no > cases use > the correction factor returned by do_sin, do_cos and TAYLOR_SIN, so > remove it. > > ChangeLog: > 2018-03-09 Wilco Dijkstra <wdijkstr@arm.com> > > * sysdeps/ieee754/dbl-64/s_sin.c (TAYLOR_SIN): Remove cor > parameter. > (do_cos): Remove corp parameter and calculations. > (do_sin): Likewise. > (do_sincos): Remove cor variable. > (__sin): Use do_sincos for huge inputs. > (__cos): Likewise. > * sysdeps/ieee754/dbl-64/s_sincos.c (__sincos): Likewise. Did this patch get mangled in transit or something? I could apply patches 1 through 3 with no problem (and no fuzzing) to ToT but when I got to 4 I get errors. I may have mangled something on this end but I treated this patch just like 1 through 3 and those worked fine. 5 and 6 have issues too but I think that is because 4 did not apply cleanly. Steve Ellcey sellcey@cavium.com % patch -p1 < ~/diff.4 patching file sysdeps/ieee754/dbl-64/s_sin.c Hunk #6 succeeded at 323 (offset -10 lines). Hunk #7 succeeded at 354 (offset -10 lines). Hunk #8 succeeded at 384 (offset -10 lines). Hunk #9 succeeded at 392 (offset -10 lines). Hunk #10 succeeded at 405 (offset -10 lines). Hunk #11 succeeded at 440 (offset -10 lines). Hunk #12 FAILED at 471. Hunk #13 succeeded at 473 (offset -10 lines). Hunk #14 succeeded at 488 (offset -10 lines). 1 out of 14 hunks FAILED -- saving rejects to file sysdeps/ieee754/dbl- 64/s_sin.c.rej patching file sysdeps/ieee754/dbl-64/s_sincos.c patch: **** malformed patch at line 240:
Steve Ellcey wrote: > Did this patch get mangled in transit or something? I could apply > patches 1 through 3 with no problem (and no fuzzing) to ToT but when I > got to 4 I get errors. I may have mangled something on this end but I > treated this patch just like 1 through 3 and those worked fine. 5 and > 6 have issues too but I think that is because 4 did not apply cleanly. Sorry, it looks like a few lines were missing due to cut&paste, here is the full version: diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c index 3b748821f6e5f817dc234ec7f96d951910299e21..5966282db60224528fea2bf55a05dd4120ab12a9 100644 --- a/sysdeps/ieee754/dbl-64/s_sin.c +++ b/sysdeps/ieee754/dbl-64/s_sin.c @@ -67,11 +67,10 @@ The constants s1, s2, s3, etc. are pre-computed values of 1/3!, 1/5! and so on. The result is returned to LHS and correction in COR. */ -#define TAYLOR_SIN(xx, a, da, cor) \ +#define TAYLOR_SIN(xx, a, da) \ ({ \ double t = ((POLYNOMIAL (xx) * (a) - 0.5 * (da)) * (xx) + (da)); \ double res = (a) + t; \ - (cor) = ((a) - res) + t; \ res; \ }) @@ -145,10 +144,10 @@ static double cslow2 (double x); /* Given a number partitioned into X and DX, this function computes the cosine of the number by combining the sin and cos of X (as computed by a variation of the Taylor series) with the values looked up from the sin/cos table to - get the result in RES and a correction value in COR. */ + get the result. */ static inline double __always_inline -do_cos (double x, double dx, double *corp) +do_cos (double x, double dx) { mynumber u; @@ -158,16 +157,13 @@ do_cos (double x, double dx, double *corp) u.x = big + fabs (x); x = fabs (x) - (u.x - big) + dx; - double xx, s, sn, ssn, c, cs, ccs, res, cor; + double xx, s, sn, ssn, c, cs, ccs, cor; xx = x * x; s = x + x * xx * (sn3 + xx * sn5); c = xx * (cs2 + xx * (cs4 + xx * cs6)); SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); cor = (ccs - s * ssn - cs * c) - sn * s; - res = cs + cor; - cor = (cs - res) + cor; - *corp = cor; - return res; + return cs + cor; } /* A more precise variant of DO_COS. EPS is the adjustment to the correction @@ -207,10 +203,10 @@ do_cos_slow (double x, double dx, double eps, double *corp) /* Given a number partitioned into X and DX, this function computes the sine of the number by combining the sin and cos of X (as computed by a variation of the Taylor series) with the values looked up from the sin/cos table to get - the result in RES and a correction value in COR. */ + the result. */ static inline double __always_inline -do_sin (double x, double dx, double *corp) +do_sin (double x, double dx) { mynumber u; @@ -219,16 +215,13 @@ do_sin (double x, double dx, double *corp) u.x = big + fabs (x); x = fabs (x) - (u.x - big); - double xx, s, sn, ssn, c, cs, ccs, cor, res; + double xx, s, sn, ssn, c, cs, ccs, cor; xx = x * x; s = x + (dx + x * xx * (sn3 + xx * sn5)); 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; - res = sn + cor; - cor = (sn - res) + cor; - *corp = cor; - return res; + return sn + cor; } /* A more precise variant of DO_SIN. EPS is the adjustment to the correction @@ -340,19 +333,19 @@ static double __always_inline do_sincos (double a, double da, int4 n) { - double retval, cor; + double retval; if (n & 1) /* Max ULP is 0.513. */ - retval = do_cos (a, da, &cor); + 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, cor); + retval = TAYLOR_SIN (xx, a, da); else - retval = __copysign (do_sin (a, da, &cor), a); + retval = __copysign (do_sin (a, da), a); } return (n & 2) ? -retval : retval; @@ -371,7 +364,7 @@ SECTION #endif __sin (double x) { - double xx, t, a, da, cor; + double xx, t, a, da; mynumber u; int4 k, m, n; double retval = 0; @@ -401,7 +394,7 @@ __sin (double x) else if (k < 0x3feb6000) { /* Max ULP is 0.548. */ - retval = __copysign (do_sin (x, 0, &cor), x); + retval = __copysign (do_sin (x, 0), x); } /* else if (k < 0x3feb6000) */ /*----------------------- 0.855469 <|x|<2.426265 ----------------------*/ @@ -409,7 +402,7 @@ __sin (double x) { t = hp0 - fabs (x); /* Max ULP is 0.51. */ - retval = __copysign (do_cos (t, hp1, &cor), x); + retval = __copysign (do_cos (t, hp1), x); } /* else if (k < 0x400368fd) */ #ifndef IN_SINCOS @@ -422,8 +415,10 @@ __sin (double x) /* --------------------105414350 <|x| <2^1024------------------------------*/ else if (k < 0x7ff00000) - retval = reduce_and_compute (x, false); - + { + n = __branred (x, &a, &da); + retval = do_sincos (a, da, n); + } /*--------------------- |x| > 2^1024 ----------------------------------*/ else { @@ -455,7 +450,7 @@ SECTION #endif __cos (double x) { - double y, xx, cor, a, da; + double y, xx, a, da; mynumber u; int4 k, m, n; @@ -476,7 +471,7 @@ __cos (double x) else if (k < 0x3feb6000) { /* 2^-27 < |x| < 0.855469 */ /* Max ULP is 0.51. */ - retval = do_cos (x, 0, &cor); + retval = do_cos (x, 0); } /* else if (k < 0x3feb6000) */ else if (k < 0x400368fd) @@ -488,9 +483,9 @@ __cos (double x) /* 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, cor); + retval = TAYLOR_SIN (xx, a, da); else - retval = __copysign (do_sin (a, da, &cor), a); + retval = __copysign (do_sin (a, da), a); } /* else if (k < 0x400368fd) */ @@ -503,7 +498,10 @@ __cos (double x) /* 105414350 <|x| <2^1024 */ else if (k < 0x7ff00000) - retval = reduce_and_compute (x, true); + { + n = __branred (x, &a, &da); + retval = do_sincos (a, da, n + 1); + } else { diff --git a/sysdeps/ieee754/dbl-64/s_sincos.c b/sysdeps/ieee754/dbl-64/s_sincos.c index 4f032d2e42593ccde22169b374728386dd8fca8e..4335ecbba3c9894e61c087ac970b392fa73abfab 100644 --- a/sysdeps/ieee754/dbl-64/s_sincos.c +++ b/sysdeps/ieee754/dbl-64/s_sincos.c @@ -28,37 +28,6 @@ #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) { @@ -88,8 +57,11 @@ __sincos (double x, double *sinx, double *cosx) } if (k < 0x7ff00000) { - reduce_and_compute_sincos (x, sinx, cosx); - return; + double a, da; + int4 n = __branred (x, &a, &da); + + *sinx = do_sincos (a, da, n); + *cosx = do_sincos (a, da, n + 1); } if (isinf (x))
diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c index 3b748821f6e5f817dc234ec7f96d951910299e21..5966282db60224528fea2bf55a05dd4120ab12a9 100644 --- a/sysdeps/ieee754/dbl-64/s_sin.c +++ b/sysdeps/ieee754/dbl-64/s_sin.c @@ -67,11 +67,10 @@ The constants s1, s2, s3, etc. are pre-computed values of 1/3!, 1/5! and so on. The result is returned to LHS and correction in COR. */ -#define TAYLOR_SIN(xx, a, da, cor) \ +#define TAYLOR_SIN(xx, a, da) \ ({ \ double t = ((POLYNOMIAL (xx) * (a) - 0.5 * (da)) * (xx) + (da)); \ double res = (a) + t; \ - (cor) = ((a) - res) + t; \ res; \ }) @@ -145,10 +144,10 @@ static double cslow2 (double x); /* Given a number partitioned into X and DX, this function computes the cosine of the number by combining the sin and cos of X (as computed by a variation of the Taylor series) with the values looked up from the sin/cos table to - get the result in RES and a correction value in COR. */ + get the result. */ static inline double __always_inline -do_cos (double x, double dx, double *corp) +do_cos (double x, double dx) { mynumber u; @@ -158,16 +157,13 @@ do_cos (double x, double dx, double *corp) u.x = big + fabs (x); x = fabs (x) - (u.x - big) + dx; - double xx, s, sn, ssn, c, cs, ccs, res, cor; + double xx, s, sn, ssn, c, cs, ccs, cor; xx = x * x; s = x + x * xx * (sn3 + xx * sn5); c = xx * (cs2 + xx * (cs4 + xx * cs6)); SINCOS_TABLE_LOOKUP (u, sn, ssn, cs, ccs); cor = (ccs - s * ssn - cs * c) - sn * s; - res = cs + cor; - cor = (cs - res) + cor; - *corp = cor; - return res; + return cs + cor; } /* A more precise variant of DO_COS. EPS is the adjustment to the correction @@ -207,10 +203,10 @@ do_cos_slow (double x, double dx, double eps, double *corp) /* Given a number partitioned into X and DX, this function computes the sine of the number by combining the sin and cos of X (as computed by a variation of the Taylor series) with the values looked up from the sin/cos table to get - the result in RES and a correction value in COR. */ + the result. */ static inline double __always_inline -do_sin (double x, double dx, double *corp) +do_sin (double x, double dx) { mynumber u; @@ -219,16 +215,13 @@ do_sin (double x, double dx, double *corp) u.x = big + fabs (x); x = fabs (x) - (u.x - big); - double xx, s, sn, ssn, c, cs, ccs, cor, res; + double xx, s, sn, ssn, c, cs, ccs, cor; xx = x * x; s = x + (dx + x * xx * (sn3 + xx * sn5)); 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; - res = sn + cor; - cor = (sn - res) + cor; - *corp = cor; - return res; + return sn + cor; } /* A more precise variant of DO_SIN. EPS is the adjustment to the correction @@ -340,19 +333,19 @@ static double __always_inline do_sincos (double a, double da, int4 n) { - double retval, cor; + double retval; if (n & 1) /* Max ULP is 0.513. */ - retval = do_cos (a, da, &cor); + 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, cor); + retval = TAYLOR_SIN (xx, a, da); else - retval = __copysign (do_sin (a, da, &cor), a); + retval = __copysign (do_sin (a, da), a); } return (n & 2) ? -retval : retval; @@ -371,7 +364,7 @@ SECTION #endif __sin (double x) { - double xx, t, a, da, cor; + double xx, t, a, da; mynumber u; int4 k, m, n; double retval = 0; @@ -401,7 +394,7 @@ __sin (double x) else if (k < 0x3feb6000) { /* Max ULP is 0.548. */ - retval = __copysign (do_sin (x, 0, &cor), x); + retval = __copysign (do_sin (x, 0), x); } /* else if (k < 0x3feb6000) */ /*----------------------- 0.855469 <|x|<2.426265 ----------------------*/ @@ -409,7 +402,7 @@ __sin (double x) { t = hp0 - fabs (x); /* Max ULP is 0.51. */ - retval = __copysign (do_cos (t, hp1, &cor), x); + retval = __copysign (do_cos (t, hp1), x); } /* else if (k < 0x400368fd) */ #ifndef IN_SINCOS @@ -422,8 +415,10 @@ __sin (double x) /* --------------------105414350 <|x| <2^1024------------------------------*/ else if (k < 0x7ff00000) - retval = reduce_and_compute (x, false); - + { + n = __branred (x, &a, &da); + retval = do_sincos (a, da, n); + } /*--------------------- |x| > 2^1024 ----------------------------------*/ else { @@ -455,7 +450,7 @@ SECTION #endif __cos (double x) { - double y, xx, cor, a, da; + double y, xx, a, da; mynumber u; int4 k, m, n; @@ -476,7 +471,7 @@ __cos (double x) else if (k < 0x3feb6000) { /* 2^-27 < |x| < 0.855469 */ /* Max ULP is 0.51. */ - retval = do_cos (x, 0, &cor); + retval = do_cos (x, 0); } /* else if (k < 0x3feb6000) */ else if (k < 0x400368fd) @@ -488,9 +483,9 @@ __cos (double x) /* 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, cor); + retval = TAYLOR_SIN (xx, a, da); else - retval = __copysign (do_sin (a, da, &cor), a); + retval = __copysign (do_sin (a, da), a); } /* else if (k < 0x400368fd) */ @@ -503,7 +498,10 @@ __cos (double x) /* 105414350 <|x| <2^1024 */ else if (k < 0x7ff00000) - retval = reduce_and_compute (x, true); + { + n = __branred (x, &a, &da); + retval = do_sincos (a, da, n + 1); + } else { diff --git a/sysdeps/ieee754/dbl-64/s_sincos.c b/sysdeps/ieee754/dbl-64/s_sincos.c index 4f032d2e42593ccde22169b374728386dd8fca8e..4335ecbba3c9894e61c087ac970b392fa73abfab 100644 --- a/sysdeps/ieee754/dbl-64/s_sincos.c +++ b/sysdeps/ieee754/dbl-64/s_sincos.c @@ -28,37 +28,6 @@ #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) { @@ -88,8 +57,11 @@ __sincos (double x, double *sinx, double *cosx) } if (k < 0x7ff00000) { - reduce_and_compute_sincos (x, sinx, cosx); - return; + double a, da; + int4 n = __branred (x, &a, &da); + + *sinx = do_sincos (a, da, n);