This is the mail archive of the
libc-alpha@sourceware.org
mailing list for the glibc project.
[PING][PATCH 3/5] Consolidate input partitioning into do_cos and do_sin
Ping!
On Tuesday 23 August 2016 11:52 PM, Siddhesh Poyarekar wrote:
> 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(-)
>
> diff --git a/sysdeps/ieee754/dbl-64/s_sin.c b/sysdeps/ieee754/dbl-64/s_sin.c
> index 7f6cd09..e03c75a 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;
> @@ -288,8 +316,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;
> @@ -309,10 +336,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));
> @@ -321,16 +345,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;
> @@ -369,7 +388,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;
>
> @@ -392,10 +410,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));
> @@ -404,16 +419,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;
> @@ -485,11 +495,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) */
>
> @@ -561,10 +567,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) */
>
> @@ -582,10 +585,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));
> @@ -655,12 +655,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;
>
> @@ -679,15 +676,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;
>
> @@ -774,17 +766,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);
> @@ -807,17 +796,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);
> @@ -870,17 +856,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);
> @@ -902,17 +884,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);
> @@ -932,18 +910,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];
>
>