xref: /freebsd/lib/msun/src/s_sin.c (revision 3a8617a83f16ffc9db4f96e1f0f21af94078e6b1)
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
143a8617a8SJordan K. Hubbard static char rcsid[] = "$Id: s_sin.c,v 1.5 1994/08/18 23:10:17 jtc Exp $";
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 
513a8617a8SJordan K. Hubbard #ifdef __STDC__
523a8617a8SJordan K. Hubbard 	double sin(double x)
533a8617a8SJordan K. Hubbard #else
543a8617a8SJordan K. Hubbard 	double sin(x)
553a8617a8SJordan K. Hubbard 	double x;
563a8617a8SJordan K. Hubbard #endif
573a8617a8SJordan K. Hubbard {
583a8617a8SJordan K. Hubbard 	double y[2],z=0.0;
593a8617a8SJordan K. Hubbard 	int32_t n, ix;
603a8617a8SJordan K. Hubbard 
613a8617a8SJordan K. Hubbard     /* High word of x. */
623a8617a8SJordan K. Hubbard 	GET_HIGH_WORD(ix,x);
633a8617a8SJordan K. Hubbard 
643a8617a8SJordan K. Hubbard     /* |x| ~< pi/4 */
653a8617a8SJordan K. Hubbard 	ix &= 0x7fffffff;
663a8617a8SJordan K. Hubbard 	if(ix <= 0x3fe921fb) return __kernel_sin(x,z,0);
673a8617a8SJordan K. Hubbard 
683a8617a8SJordan K. Hubbard     /* sin(Inf or NaN) is NaN */
693a8617a8SJordan K. Hubbard 	else if (ix>=0x7ff00000) return x-x;
703a8617a8SJordan K. Hubbard 
713a8617a8SJordan K. Hubbard     /* argument reduction needed */
723a8617a8SJordan K. Hubbard 	else {
733a8617a8SJordan K. Hubbard 	    n = __ieee754_rem_pio2(x,y);
743a8617a8SJordan K. Hubbard 	    switch(n&3) {
753a8617a8SJordan K. Hubbard 		case 0: return  __kernel_sin(y[0],y[1],1);
763a8617a8SJordan K. Hubbard 		case 1: return  __kernel_cos(y[0],y[1]);
773a8617a8SJordan K. Hubbard 		case 2: return -__kernel_sin(y[0],y[1],1);
783a8617a8SJordan K. Hubbard 		default:
793a8617a8SJordan K. Hubbard 			return -__kernel_cos(y[0],y[1]);
803a8617a8SJordan K. Hubbard 	    }
813a8617a8SJordan K. Hubbard 	}
823a8617a8SJordan K. Hubbard }
83