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