| OLD | NEW |
| 1 /* origin: FreeBSD /usr/src/lib/msun/src/s_sin.c */ | 1 /* origin: FreeBSD /usr/src/lib/msun/src/s_sin.c */ |
| 2 /* | 2 /* |
| 3 * ==================================================== | 3 * ==================================================== |
| 4 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. | 4 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. |
| 5 * | 5 * |
| 6 * Developed at SunPro, a Sun Microsystems, Inc. business. | 6 * Developed at SunPro, a Sun Microsystems, Inc. business. |
| 7 * Permission to use, copy, modify, and distribute this | 7 * Permission to use, copy, modify, and distribute this |
| 8 * software is freely granted, provided that this notice | 8 * software is freely granted, provided that this notice |
| 9 * is preserved. | 9 * is preserved. |
| 10 * ==================================================== | 10 * ==================================================== |
| 11 */ | 11 */ |
| 12 | 12 |
| 13 #define _GNU_SOURCE | 13 #define _GNU_SOURCE |
| 14 #include "libm.h" | 14 #include "libm.h" |
| 15 | 15 |
| 16 void sincos(double x, double *sin, double *cos) | 16 void sincos(double x, double* sin, double* cos) { |
| 17 { | 17 double y[2], s, c; |
| 18 » double y[2], s, c; | 18 uint32_t ix; |
| 19 » uint32_t ix; | 19 unsigned n; |
| 20 » unsigned n; | |
| 21 | 20 |
| 22 » GET_HIGH_WORD(ix, x); | 21 GET_HIGH_WORD(ix, x); |
| 23 » ix &= 0x7fffffff; | 22 ix &= 0x7fffffff; |
| 24 | 23 |
| 25 » /* |x| ~< pi/4 */ | 24 /* |x| ~< pi/4 */ |
| 26 » if (ix <= 0x3fe921fb) { | 25 if (ix <= 0x3fe921fb) { |
| 27 » » /* if |x| < 2**-27 * sqrt(2) */ | 26 /* if |x| < 2**-27 * sqrt(2) */ |
| 28 » » if (ix < 0x3e46a09e) { | 27 if (ix < 0x3e46a09e) { |
| 29 » » » /* raise inexact if x!=0 and underflow if subnormal */ | 28 /* raise inexact if x!=0 and underflow if subnormal */ |
| 30 » » » FORCE_EVAL(ix < 0x00100000 ? x/0x1p120f : x+0x1p120f); | 29 FORCE_EVAL(ix < 0x00100000 ? x / 0x1p120f : x + 0x1p120f); |
| 31 » » » *sin = x; | 30 *sin = x; |
| 32 » » » *cos = 1.0; | 31 *cos = 1.0; |
| 33 » » » return; | 32 return; |
| 34 » » } | 33 } |
| 35 » » *sin = __sin(x, 0.0, 0); | 34 *sin = __sin(x, 0.0, 0); |
| 36 » » *cos = __cos(x, 0.0); | 35 *cos = __cos(x, 0.0); |
| 37 » » return; | 36 return; |
| 38 » } | 37 } |
| 39 | 38 |
| 40 » /* sincos(Inf or NaN) is NaN */ | 39 /* sincos(Inf or NaN) is NaN */ |
| 41 » if (ix >= 0x7ff00000) { | 40 if (ix >= 0x7ff00000) { |
| 42 » » *sin = *cos = x - x; | 41 *sin = *cos = x - x; |
| 43 » » return; | 42 return; |
| 44 » } | 43 } |
| 45 | 44 |
| 46 » /* argument reduction needed */ | 45 /* argument reduction needed */ |
| 47 » n = __rem_pio2(x, y); | 46 n = __rem_pio2(x, y); |
| 48 » s = __sin(y[0], y[1], 1); | 47 s = __sin(y[0], y[1], 1); |
| 49 » c = __cos(y[0], y[1]); | 48 c = __cos(y[0], y[1]); |
| 50 » switch (n&3) { | 49 switch (n & 3) { |
| 51 » case 0: | 50 case 0: |
| 52 » » *sin = s; | 51 *sin = s; |
| 53 » » *cos = c; | 52 *cos = c; |
| 54 » » break; | 53 break; |
| 55 » case 1: | 54 case 1: |
| 56 » » *sin = c; | 55 *sin = c; |
| 57 » » *cos = -s; | 56 *cos = -s; |
| 58 » » break; | 57 break; |
| 59 » case 2: | 58 case 2: |
| 60 » » *sin = -s; | 59 *sin = -s; |
| 61 » » *cos = -c; | 60 *cos = -c; |
| 62 » » break; | 61 break; |
| 63 » case 3: | 62 case 3: |
| 64 » default: | 63 default: |
| 65 » » *sin = -c; | 64 *sin = -c; |
| 66 » » *cos = s; | 65 *cos = s; |
| 67 » » break; | 66 break; |
| 68 » } | 67 } |
| 69 } | 68 } |
| OLD | NEW |