| Elliott Hughes | e332f65 | 2018-05-08 15:07:43 -0700 | [diff] [blame] | 1 | /*- |
| 2 | * ==================================================== |
| 3 | * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. |
| 4 | * |
| 5 | * Developed at SunPro, a Sun Microsystems, Inc. business. |
| 6 | * Permission to use, copy, modify, and distribute this |
| 7 | * software is freely granted, provided that this notice |
| 8 | * is preserved. |
| 9 | * ==================================================== |
| 10 | * |
| 11 | * s_sin.c and s_cos.c merged by Steven G. Kargl. Descriptions of the |
| 12 | * algorithms are contained in the original files. |
| 13 | */ |
| 14 | |
| Elliott Hughes | e332f65 | 2018-05-08 15:07:43 -0700 | [diff] [blame] | 15 | #include <float.h> |
| 16 | |
| 17 | #include "math.h" |
| 18 | #define INLINE_REM_PIO2 |
| 19 | #include "math_private.h" |
| 20 | #include "e_rem_pio2.c" |
| 21 | #include "k_sincos.h" |
| 22 | |
| 23 | void |
| 24 | sincos(double x, double *sn, double *cs) |
| 25 | { |
| 26 | double y[2]; |
| 27 | int32_t n, ix; |
| 28 | |
| 29 | /* High word of x. */ |
| 30 | GET_HIGH_WORD(ix, x); |
| 31 | |
| 32 | /* |x| ~< pi/4 */ |
| 33 | ix &= 0x7fffffff; |
| 34 | if (ix <= 0x3fe921fb) { |
| 35 | if (ix < 0x3e400000) { /* |x| < 2**-27 */ |
| 36 | if ((int)x == 0) { /* Generate inexact. */ |
| 37 | *sn = x; |
| 38 | *cs = 1; |
| 39 | return; |
| 40 | } |
| 41 | } |
| 42 | __kernel_sincos(x, 0, 0, sn, cs); |
| 43 | return; |
| 44 | } |
| 45 | |
| 46 | /* If x = Inf or NaN, then sin(x) = NaN and cos(x) = NaN. */ |
| 47 | if (ix >= 0x7ff00000) { |
| 48 | *sn = x - x; |
| 49 | *cs = x - x; |
| 50 | return; |
| 51 | } |
| 52 | |
| 53 | /* Argument reduction. */ |
| 54 | n = __ieee754_rem_pio2(x, y); |
| 55 | |
| 56 | switch(n & 3) { |
| 57 | case 0: |
| 58 | __kernel_sincos(y[0], y[1], 1, sn, cs); |
| 59 | break; |
| 60 | case 1: |
| 61 | __kernel_sincos(y[0], y[1], 1, cs, sn); |
| 62 | *cs = -*cs; |
| 63 | break; |
| 64 | case 2: |
| 65 | __kernel_sincos(y[0], y[1], 1, sn, cs); |
| 66 | *sn = -*sn; |
| 67 | *cs = -*cs; |
| 68 | break; |
| 69 | default: |
| 70 | __kernel_sincos(y[0], y[1], 1, cs, sn); |
| 71 | *sn = -*sn; |
| 72 | } |
| 73 | } |
| 74 | |
| 75 | #if (LDBL_MANT_DIG == 53) |
| 76 | __weak_reference(sincos, sincosl); |
| 77 | #endif |