src/libm/s_sin.c
changeset 2756 a98604b691c8
child 3162 dc1eb82ffdaa
     1.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     1.2 +++ b/src/libm/s_sin.c	Mon Sep 15 06:33:23 2008 +0000
     1.3 @@ -0,0 +1,90 @@
     1.4 +/* @(#)s_sin.c 5.1 93/09/24 */
     1.5 +/*
     1.6 + * ====================================================
     1.7 + * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
     1.8 + *
     1.9 + * Developed at SunPro, a Sun Microsystems, Inc. business.
    1.10 + * Permission to use, copy, modify, and distribute this
    1.11 + * software is freely granted, provided that this notice
    1.12 + * is preserved.
    1.13 + * ====================================================
    1.14 + */
    1.15 +
    1.16 +#if defined(LIBM_SCCS) && !defined(lint)
    1.17 +static char rcsid[] = "$NetBSD: s_sin.c,v 1.7 1995/05/10 20:48:15 jtc Exp $";
    1.18 +#endif
    1.19 +
    1.20 +/* sin(x)
    1.21 + * Return sine function of x.
    1.22 + *
    1.23 + * kernel function:
    1.24 + *	__kernel_sin		... sine function on [-pi/4,pi/4]
    1.25 + *	__kernel_cos		... cose function on [-pi/4,pi/4]
    1.26 + *	__ieee754_rem_pio2	... argument reduction routine
    1.27 + *
    1.28 + * Method.
    1.29 + *      Let S,C and T denote the sin, cos and tan respectively on
    1.30 + *	[-PI/4, +PI/4]. Reduce the argument x to y1+y2 = x-k*pi/2
    1.31 + *	in [-pi/4 , +pi/4], and let n = k mod 4.
    1.32 + *	We have
    1.33 + *
    1.34 + *          n        sin(x)      cos(x)        tan(x)
    1.35 + *     ----------------------------------------------------------
    1.36 + *	    0	       S	   C		 T
    1.37 + *	    1	       C	  -S		-1/T
    1.38 + *	    2	      -S	  -C		 T
    1.39 + *	    3	      -C	   S		-1/T
    1.40 + *     ----------------------------------------------------------
    1.41 + *
    1.42 + * Special cases:
    1.43 + *      Let trig be any of sin, cos, or tan.
    1.44 + *      trig(+-INF)  is NaN, with signals;
    1.45 + *      trig(NaN)    is that NaN;
    1.46 + *
    1.47 + * Accuracy:
    1.48 + *	TRIG(x) returns trig(x) nearly rounded
    1.49 + */
    1.50 +
    1.51 +#include "math.h"
    1.52 +#include "math_private.h"
    1.53 +
    1.54 +libm_hidden_proto(sin)
    1.55 +#ifdef __STDC__
    1.56 +     double sin(double x)
    1.57 +#else
    1.58 +     double sin(x)
    1.59 +     double x;
    1.60 +#endif
    1.61 +{
    1.62 +    double y[2], z = 0.0;
    1.63 +    int32_t n, ix;
    1.64 +
    1.65 +    /* High word of x. */
    1.66 +    GET_HIGH_WORD(ix, x);
    1.67 +
    1.68 +    /* |x| ~< pi/4 */
    1.69 +    ix &= 0x7fffffff;
    1.70 +    if (ix <= 0x3fe921fb)
    1.71 +        return __kernel_sin(x, z, 0);
    1.72 +
    1.73 +    /* sin(Inf or NaN) is NaN */
    1.74 +    else if (ix >= 0x7ff00000)
    1.75 +        return x - x;
    1.76 +
    1.77 +    /* argument reduction needed */
    1.78 +    else {
    1.79 +        n = __ieee754_rem_pio2(x, y);
    1.80 +        switch (n & 3) {
    1.81 +        case 0:
    1.82 +            return __kernel_sin(y[0], y[1], 1);
    1.83 +        case 1:
    1.84 +            return __kernel_cos(y[0], y[1]);
    1.85 +        case 2:
    1.86 +            return -__kernel_sin(y[0], y[1], 1);
    1.87 +        default:
    1.88 +            return -__kernel_cos(y[0], y[1]);
    1.89 +        }
    1.90 +    }
    1.91 +}
    1.92 +
    1.93 +libm_hidden_def(sin)