From 758e79ec89c4e285bafcef601702cd27790fd78e Mon Sep 17 00:00:00 2001 From: Siddhesh Poyarekar Date: Thu, 1 Sep 2016 22:08:51 +0530 Subject: Consolidate input partitioning into do_cos and do_sin All calls to do_cos are preceded by code that partitions x into a larger double that gives an offset into the sincos table and a smaller double that is used in a polynomial computation. Consolidate all of them into do_cos and do_sin to reduce code duplication. * sysdeps/ieee754/dbl-64/s_sin.c (do_cos): Accept X and DX as input arguments. Consolidate input partitioning from callers here. (do_cos_slow): Likewise. (do_sin): Likewise. (do_sin_slow): Likewise. (do_sincos_1): Remove the no longer necessary input partitioning. (do_sincos_2): Likewise. (__sin): Likewise. (__cos): Likewise. (slow1): Likewise. (slow2): Likewise. (sloww1): Likewise. (sloww2): Likewise. (bsloww1): Likewise. (bsloww2): Likewise. (cslow2): Likewise. --- sysdeps/ieee754/dbl-64/s_sin.c | 191 ++++++++++++++++++----------------------- 1 file changed, 82 insertions(+), 109 deletions(-) (limited to 'sysdeps/ieee754') diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c index 9ecba05db2..888270f91d 100644 --- a/sysdeps/ieee754/dbl-64/s_sin.c +++ b/sysdeps/ieee754/dbl-64/s_sin.c @@ -141,14 +141,21 @@ static double bsloww2 (double x, double dx, double orig, int n); int __branred (double x, double *a, double *aa); static double cslow2 (double x); -/* 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 - 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. */ +/* 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. */ static double -do_cos (mynumber u, double x, double *corp) +do_cos (double x, double dx, double *corp) { + mynumber u; + + if (x < 0) + dx = -dx; + + u.x = big + fabs (x); + x = fabs (x) - (u.x - big) + dx; + double xx, s, sn, ssn, c, cs, ccs, res, cor; xx = x * x; s = x + x * xx * (sn3 + xx * sn5); @@ -161,11 +168,19 @@ do_cos (mynumber u, double x, double *corp) return res; } -/* A more precise variant of DO_COS where the number is partitioned into U, X - and DX. EPS is the adjustment to the correction COR. */ +/* A more precise variant of DO_COS. EPS is the adjustment to the correction + COR. */ static double -do_cos_slow (mynumber u, double x, double dx, double eps, double *corp) +do_cos_slow (double x, double dx, double eps, double *corp) { + mynumber u; + + if (x <= 0) + dx = -dx; + + u.x = big + fabs (x); + x = fabs (x) - (u.x - big); + double xx, y, x1, x2, e1, e2, res, cor; double s, sn, ssn, c, cs, ccs; xx = x * x; @@ -186,14 +201,20 @@ do_cos_slow (mynumber u, double x, double dx, double eps, double *corp) return res; } -/* Given a number partitioned into U and X and DX such that U is an index into - the sin/cos table, this macro 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. */ +/* 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. */ static double -do_sin (mynumber u, double x, double dx, double *corp) +do_sin (double x, double dx, double *corp) { + mynumber u; + + if (x <= 0) + dx = -dx; + u.x = big + fabs (x); + x = fabs (x) - (u.x - big); + double xx, s, sn, ssn, c, cs, ccs, cor, res; xx = x * x; s = x + (dx + x * xx * (sn3 + xx * sn5)); @@ -206,11 +227,18 @@ do_sin (mynumber u, double x, double dx, double *corp) return res; } -/* A more precise variant of res = do_sin where the number is partitioned into U, X - and DX. EPS is the adjustment to the correction COR. */ +/* A more precise variant of DO_SIN. EPS is the adjustment to the correction + COR. */ static double -do_sin_slow (mynumber u, double x, double dx, double eps, double *corp) +do_sin_slow (double x, double dx, double eps, double *corp) { + mynumber u; + + if (x <= 0) + dx = -dx; + u.x = big + fabs (x); + x = fabs (x) - (u.x - big); + double xx, y, x1, x2, c1, c2, res, cor; double s, sn, ssn, c, cs, ccs; xx = x * x; @@ -289,8 +317,7 @@ 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; + double xx, retval, res, cor; double eps = fabs (x) * 1.2e-30; int k1 = (n + k) & 3; @@ -311,10 +338,7 @@ do_sincos_1 (double a, double da, double x, int4 n, int4 k) } else { - double db = (a > 0 ? da : -da); - u.x = big + fabs (a); - y = fabs (a) - (u.x - big); - res = do_sin (u, y, db, &cor); + res = do_sin (a, da, &cor); cor = (cor > 0) ? 1.035 * cor + eps : 1.035 * cor - eps; retval = ((res == res + cor) ? ((a > 0) ? res : -res) : sloww1 (a, da, x, k)); @@ -323,16 +347,11 @@ do_sincos_1 (double a, double da, double x, int4 n, int4 k) case 1: case 3: - { - double db = (a > 0 ? da : -da); - u.x = big + fabs (a); - y = fabs (a) - (u.x - big) + db; - 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; - } + res = do_cos (a, da, &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; @@ -371,7 +390,6 @@ __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; @@ -394,10 +412,7 @@ do_sincos_2 (double a, double da, double x, int4 n, int4 k) } else { - double db = (a > 0 ? da : -da); - u.x = big + fabs (a); - double y = fabs (a) - (u.x - big); - res = do_sin (u, y, db, &cor); + res = do_sin (a, da, &cor); cor = (cor > 0) ? 1.035 * cor + eps : 1.035 * cor - eps; retval = ((res == res + cor) ? ((a > 0) ? res : -res) : bsloww1 (a, da, x, n)); @@ -406,16 +421,11 @@ do_sincos_2 (double a, double da, double x, int4 n, int4 k) case 1: case 3: - { - double db = (a > 0 ? da : -da); - u.x = big + fabs (a); - double y = fabs (a) - (u.x - big) + db; - 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; - } + res = do_cos (a, da, &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; @@ -487,11 +497,7 @@ __sin (double x) { t = hp0 - fabs (x); - u.x = big + fabs (t); - y = fabs (t) - (u.x - big); - y = ((t >= 0) ? hp1 : -hp1) + y; - - res = do_cos (u, y, &cor); + res = do_cos (t, hp1, &cor); retval = (res == res + 1.020 * cor) ? ((m > 0) ? res : -res) : slow2 (x); } /* else if (k < 0x400368fd) */ @@ -563,10 +569,7 @@ __cos (double x) else if (k < 0x3feb6000) { /* 2^-27 < |x| < 0.855469 */ - y = fabs (x); - u.x = big + y; - y = y - (u.x - big); - res = do_cos (u, y, &cor); + res = do_cos (x, 0, &cor); retval = (res == res + 1.020 * cor) ? res : cslow2 (x); } /* else if (k < 0x3feb6000) */ @@ -584,10 +587,7 @@ __cos (double x) } else { - double db = (a > 0 ? da : -da); - u.x = big + fabs (a); - y = fabs (a) - (u.x - big); - res = do_sin (u, y, db, &cor); + res = do_sin (a, da, &cor); cor = (cor > 0) ? 1.035 * cor + 1.0e-31 : 1.035 * cor - 1.0e-31; retval = ((res == res + cor) ? ((a > 0) ? res : -res) : sloww1 (a, da, x, 1)); @@ -657,12 +657,9 @@ static double SECTION slow1 (double x) { - mynumber u; - double w[2], y, cor, res; - y = fabs (x); - u.x = big + y; - y = y - (u.x - big); - res = do_sin_slow (u, y, 0, 0, &cor); + double w[2], cor, res; + + res = do_sin_slow (x, 0, 0, &cor); if (res == res + cor) return (x > 0) ? res : -res; @@ -681,15 +678,10 @@ static double SECTION slow2 (double x) { - mynumber u; - double w[2], y, y1, y2, cor, res, del; + double w[2], y, y1, y2, cor, res; double t = hp0 - fabs (x); - u.x = big + fabs (t); - y = fabs (t) - (u.x - big); - del = (t >= 0) ? hp1 : -hp1; - - res = do_cos_slow (u, y, del, 0, &cor); + res = do_cos_slow (t, hp1, 0, &cor); if (res == res + cor) return (x > 0) ? res : -res; @@ -776,17 +768,14 @@ static double SECTION sloww1 (double x, double dx, double orig, int k) { - mynumber u; - double w[2], y, cor, res; + double w[2], cor, res; - u.x = big + fabs (x); - y = fabs (x) - (u.x - big); - dx = (x > 0 ? dx : -dx); - res = do_sin_slow (u, y, dx, 3.1e-30 * fabs (orig), &cor); + res = do_sin_slow (x, dx, 3.1e-30 * fabs (orig), &cor); if (res == res + cor) return (x > 0) ? res : -res; + dx = (x > 0 ? dx : -dx); __dubsin (fabs (x), dx, w); double eps = 1.1e-30 * fabs (orig); @@ -809,17 +798,14 @@ static double SECTION sloww2 (double x, double dx, double orig, int n) { - mynumber u; - double w[2], y, cor, res; + double w[2], cor, res; - u.x = big + fabs (x); - y = fabs (x) - (u.x - big); - dx = (x > 0 ? dx : -dx); - res = do_cos_slow (u, y, dx, 3.1e-30 * fabs (orig), &cor); + res = do_cos_slow (x, dx, 3.1e-30 * fabs (orig), &cor); if (res == res + cor) return (n & 2) ? -res : res; + dx = x > 0 ? dx : -dx; __docos (fabs (x), dx, w); double eps = 1.1e-30 * fabs (orig); @@ -872,17 +858,13 @@ static double SECTION bsloww1 (double x, double dx, double orig, int n) { - mynumber u; - double w[2], y, cor, res; + double w[2], cor, res; - y = fabs (x); - u.x = big + y; - y = y - (u.x - big); - dx = (x > 0) ? dx : -dx; - res = do_sin_slow (u, y, dx, 1.1e-24, &cor); + res = do_sin_slow (x, dx, 1.1e-24, &cor); if (res == res + cor) return (x > 0) ? res : -res; + dx = (x > 0) ? dx : -dx; __dubsin (fabs (x), dx, w); cor = 1.000000005 * w[1] + ((w[1] > 0) ? 1.1e-24 : -1.1e-24); @@ -904,17 +886,13 @@ static double SECTION bsloww2 (double x, double dx, double orig, int n) { - mynumber u; - double w[2], y, cor, res; + double w[2], cor, res; - y = fabs (x); - u.x = big + y; - y = y - (u.x - big); - dx = (x > 0) ? dx : -dx; - res = do_cos_slow (u, y, dx, 1.1e-24, &cor); + res = do_cos_slow (x, dx, 1.1e-24, &cor); if (res == res + cor) return (n & 2) ? -res : res; + dx = (x > 0) ? dx : -dx; __docos (fabs (x), dx, w); cor = 1.000000005 * w[1] + ((w[1] > 0) ? 1.1e-24 : -1.1e-24); @@ -934,18 +912,13 @@ static double SECTION cslow2 (double x) { - mynumber u; - double w[2], y, cor, res; + double w[2], cor, res; - y = fabs (x); - u.x = big + y; - y = y - (u.x - big); - res = do_cos_slow (u, y, 0, 0, &cor); + res = do_cos_slow (x, 0, 0, &cor); if (res == res + cor) return res; - y = fabs (x); - __docos (y, 0, w); + __docos (fabs (x), 0, w); if (w[0] == w[0] + 1.000000005 * w[1]) return w[0]; -- cgit 1.4.1