xref: /freebsd/lib/msun/src/s_sin.c (revision 3819e84017471706d911eec3f1b3b6f9eefc139a)
13a8617a8SJordan K. Hubbard /* @(#)s_sin.c 5.1 93/09/24 */
23a8617a8SJordan K. Hubbard /*
33a8617a8SJordan K. Hubbard  * ====================================================
43a8617a8SJordan K. Hubbard  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
53a8617a8SJordan K. Hubbard  *
63a8617a8SJordan K. Hubbard  * Developed at SunPro, a Sun Microsystems, Inc. business.
73a8617a8SJordan K. Hubbard  * Permission to use, copy, modify, and distribute this
83a8617a8SJordan K. Hubbard  * software is freely granted, provided that this notice
93a8617a8SJordan K. Hubbard  * is preserved.
103a8617a8SJordan K. Hubbard  * ====================================================
113a8617a8SJordan K. Hubbard  */
123a8617a8SJordan K. Hubbard 
133a8617a8SJordan K. Hubbard #ifndef lint
147f3dea24SPeter Wemm static char rcsid[] = "$FreeBSD$";
153a8617a8SJordan K. Hubbard #endif
163a8617a8SJordan K. Hubbard 
173a8617a8SJordan K. Hubbard /* sin(x)
183a8617a8SJordan K. Hubbard  * Return sine function of x.
193a8617a8SJordan K. Hubbard  *
203a8617a8SJordan K. Hubbard  * kernel function:
213a8617a8SJordan K. Hubbard  *	__kernel_sin		... sine function on [-pi/4,pi/4]
223a8617a8SJordan K. Hubbard  *	__kernel_cos		... cose function on [-pi/4,pi/4]
233a8617a8SJordan K. Hubbard  *	__ieee754_rem_pio2	... argument reduction routine
243a8617a8SJordan K. Hubbard  *
253a8617a8SJordan K. Hubbard  * Method.
263a8617a8SJordan K. Hubbard  *      Let S,C and T denote the sin, cos and tan respectively on
273a8617a8SJordan K. Hubbard  *	[-PI/4, +PI/4]. Reduce the argument x to y1+y2 = x-k*pi/2
283a8617a8SJordan K. Hubbard  *	in [-pi/4 , +pi/4], and let n = k mod 4.
293a8617a8SJordan K. Hubbard  *	We have
303a8617a8SJordan K. Hubbard  *
313a8617a8SJordan K. Hubbard  *          n        sin(x)      cos(x)        tan(x)
323a8617a8SJordan K. Hubbard  *     ----------------------------------------------------------
333a8617a8SJordan K. Hubbard  *	    0	       S	   C		 T
343a8617a8SJordan K. Hubbard  *	    1	       C	  -S		-1/T
353a8617a8SJordan K. Hubbard  *	    2	      -S	  -C		 T
363a8617a8SJordan K. Hubbard  *	    3	      -C	   S		-1/T
373a8617a8SJordan K. Hubbard  *     ----------------------------------------------------------
383a8617a8SJordan K. Hubbard  *
393a8617a8SJordan K. Hubbard  * Special cases:
403a8617a8SJordan K. Hubbard  *      Let trig be any of sin, cos, or tan.
413a8617a8SJordan K. Hubbard  *      trig(+-INF)  is NaN, with signals;
423a8617a8SJordan K. Hubbard  *      trig(NaN)    is that NaN;
433a8617a8SJordan K. Hubbard  *
443a8617a8SJordan K. Hubbard  * Accuracy:
453a8617a8SJordan K. Hubbard  *	TRIG(x) returns trig(x) nearly rounded
463a8617a8SJordan K. Hubbard  */
473a8617a8SJordan K. Hubbard 
483a8617a8SJordan K. Hubbard #include "math.h"
493a8617a8SJordan K. Hubbard #include "math_private.h"
503a8617a8SJordan K. Hubbard 
5159b19ff1SAlfred Perlstein double
523819e840SPeter Wemm sin(double x)
533a8617a8SJordan K. Hubbard {
543a8617a8SJordan K. Hubbard 	double y[2],z=0.0;
553a8617a8SJordan K. Hubbard 	int32_t n, ix;
563a8617a8SJordan K. Hubbard 
573a8617a8SJordan K. Hubbard     /* High word of x. */
583a8617a8SJordan K. Hubbard 	GET_HIGH_WORD(ix,x);
593a8617a8SJordan K. Hubbard 
603a8617a8SJordan K. Hubbard     /* |x| ~< pi/4 */
613a8617a8SJordan K. Hubbard 	ix &= 0x7fffffff;
623a8617a8SJordan K. Hubbard 	if(ix <= 0x3fe921fb) return __kernel_sin(x,z,0);
633a8617a8SJordan K. Hubbard 
643a8617a8SJordan K. Hubbard     /* sin(Inf or NaN) is NaN */
653a8617a8SJordan K. Hubbard 	else if (ix>=0x7ff00000) return x-x;
663a8617a8SJordan K. Hubbard 
673a8617a8SJordan K. Hubbard     /* argument reduction needed */
683a8617a8SJordan K. Hubbard 	else {
693a8617a8SJordan K. Hubbard 	    n = __ieee754_rem_pio2(x,y);
703a8617a8SJordan K. Hubbard 	    switch(n&3) {
713a8617a8SJordan K. Hubbard 		case 0: return  __kernel_sin(y[0],y[1],1);
723a8617a8SJordan K. Hubbard 		case 1: return  __kernel_cos(y[0],y[1]);
733a8617a8SJordan K. Hubbard 		case 2: return -__kernel_sin(y[0],y[1],1);
743a8617a8SJordan K. Hubbard 		default:
753a8617a8SJordan K. Hubbard 			return -__kernel_cos(y[0],y[1]);
763a8617a8SJordan K. Hubbard 	    }
773a8617a8SJordan K. Hubbard 	}
783a8617a8SJordan K. Hubbard }
79