diff options
Diffstat (limited to 'sysdeps/ieee754/dbl-64/s_sincos.c')
-rw-r--r-- | sysdeps/ieee754/dbl-64/s_sincos.c | 88 |
1 files changed, 74 insertions, 14 deletions
diff --git a/sysdeps/ieee754/dbl-64/s_sincos.c b/sysdeps/ieee754/dbl-64/s_sincos.c index c8a99991cc..c389226b04 100644 --- a/sysdeps/ieee754/dbl-64/s_sincos.c +++ b/sysdeps/ieee754/dbl-64/s_sincos.c @@ -1,5 +1,5 @@ /* Compute sine and cosine of argument. - Copyright (C) 1997-2015 Free Software Foundation, Inc. + Copyright (C) 1997-2016 Free Software Foundation, Inc. This file is part of the GNU C Library. Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. @@ -22,29 +22,89 @@ #include <math_private.h> +#define __sin __sin_local +#define __cos __cos_local +#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) { - int32_t ix; + mynumber u; + int k; - /* High word of x. */ - GET_HIGH_WORD (ix, x); + SET_RESTORE_ROUND_53BIT (FE_TONEAREST); - /* |x| ~< pi/4 */ - ix &= 0x7fffffff; - if (ix >= 0x7ff00000) + u.x = x; + k = 0x7fffffff & u.i[HIGH_HALF]; + + if (k < 0x400368fd) { - /* sin(Inf or NaN) is NaN */ - *sinx = *cosx = x - x; - if (__isinf_ns (x)) - __set_errno (EDOM); + *sinx = __sin_local (x); + *cosx = __cos_local (x); + return; } - else + 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) { - *sinx = __sin (x); - *cosx = __cos (x); + 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); + return; } + + if (isinf (x)) + __set_errno (EDOM); + + *sinx = *cosx = x / x; } weak_alias (__sincos, sincos) #ifdef NO_LONG_DOUBLE |