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