1 /* s_cosl.c -- long double version of s_cos.c.
2  */
3 
4 /*
5  * ====================================================
6  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
7  *
8  * Developed at SunPro, a Sun Microsystems, Inc. business.
9  * Permission to use, copy, modify, and distribute this
10  * software is freely granted, provided that this notice
11  * is preserved.
12  * ====================================================
13  */
14 
15 /* cosl(x)
16  * Return cosine function of x.
17  *
18  * kernel function:
19  *	__kernel_sinl		... sine function on [-pi/4,pi/4]
20  *	__kernel_cosl		... cosine function on [-pi/4,pi/4]
21  *	__ieee754_rem_pio2l	... argument reduction routine
22  *
23  * Method.
24  *      Let S,C and T denote the sin, cos and tan respectively on
25  *	[-PI/4, +PI/4]. Reduce the argument x to y1+y2 = x-k*pi/2
26  *	in [-pi/4 , +pi/4], and let n = k mod 4.
27  *	We have
28  *
29  *          n        sin(x)      cos(x)        tan(x)
30  *     ----------------------------------------------------------
31  *	    0	       S	   C		 T
32  *	    1	       C	  -S		-1/T
33  *	    2	      -S	  -C		 T
34  *	    3	      -C	   S		-1/T
35  *     ----------------------------------------------------------
36  *
37  * Special cases:
38  *      Let trig be any of sin, cos, or tan.
39  *      trig(+-INF)  is NaN, with signals;
40  *      trig(NaN)    is that NaN;
41  *
42  * Accuracy:
43  *	TRIG(x) returns trig(x) nearly rounded
44  */
45 
46 #include <errno.h>
47 #include <math.h>
48 #include <math_private.h>
49 #include <libm-alias-ldouble.h>
50 
__cosl(_Float128 x)51 _Float128 __cosl(_Float128 x)
52 {
53 	_Float128 y[2],z=0;
54 	int64_t n, ix;
55 
56     /* High word of x. */
57 	GET_LDOUBLE_MSW64(ix,x);
58 
59     /* |x| ~< pi/4 */
60 	ix &= 0x7fffffffffffffffLL;
61 	if(ix <= 0x3ffe921fb54442d1LL)
62 	  return __kernel_cosl(x,z);
63 
64     /* cos(Inf or NaN) is NaN */
65 	else if (ix>=0x7fff000000000000LL) {
66 	    if (ix == 0x7fff000000000000LL) {
67 		GET_LDOUBLE_LSW64(n,x);
68 		if (n == 0)
69 		    __set_errno (EDOM);
70 	    }
71 	    return x-x;
72 	}
73 
74     /* argument reduction needed */
75 	else {
76 	    n = __ieee754_rem_pio2l(x,y);
77 	    switch(n&3) {
78 		case 0: return  __kernel_cosl(y[0],y[1]);
79 		case 1: return -__kernel_sinl(y[0],y[1],1);
80 		case 2: return -__kernel_cosl(y[0],y[1]);
81 		default:
82 		        return  __kernel_sinl(y[0],y[1],1);
83 	    }
84 	}
85 }
86 libm_alias_ldouble (__cos, cos)
87