4 * The contents of this file are subject to the terms of the
5 * Common Development and Distribution License (the "License").
6 * You may not use this file except in compliance with the License.
8 * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
9 * or http://www.opensolaris.org/os/licensing.
10 * See the License for the specific language governing permissions
11 * and limitations under the License.
13 * When distributing Covered Code, include this CDDL HEADER in each
14 * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
15 * If applicable, add the following below this CDDL HEADER, with the
16 * fields enclosed by brackets "[]" replaced with your own identifying
17 * information: Portions Copyright [yyyy] [name of copyright owner]
23 * Copyright 2011 Nexenta Systems, Inc. All rights reserved.
26 * Copyright 2006 Sun Microsystems, Inc. All rights reserved.
27 * Use is subject to license terms.
32 * void sincospi(double x, double *s, double *c)
33 * *s = sin(pi*x); *c = cos(pi*x);
35 * Algorithm, 10/17/2002, K.C. Ng
36 * ------------------------------
37 * Let y = |4x|, z = floor(y), and n = (int)(z mod 8.0) (displayed in binary).
38 * 1. If y == z, then x is a multiple of pi/4. Return the following values:
39 * ---------------------------------------------------
40 * n x mod 2 sin(x*pi) cos(x*pi) tan(x*pi)
41 * ---------------------------------------------------
42 * 000 0.00 +0 ___ +1 ___ +0
43 * 001 0.25 +\/0.5 +\/0.5 +1
44 * 010 0.50 +1 ___ +0 ___ +inf
45 * 011 0.75 +\/0.5 -\/0.5 -1
46 * 100 1.00 -0 ___ -1 ___ +0
47 * 101 1.25 -\/0.5 -\/0.5 +1
48 * 110 1.50 -1 ___ -0 ___ +inf
49 * 111 1.75 -\/0.5 +\/0.5 -1
50 * ---------------------------------------------------
52 * ---------------------------------------------------
53 * n t sin(x*pi) cos(x*pi) tan(x*pi)
54 * ---------------------------------------------------
55 * 000 (y-z)/4 sinpi(t) cospi(t) tanpi(t)
56 * 001 (z+1-y)/4 cospi(t) sinpi(t) 1/tanpi(t)
57 * 010 (y-z)/4 cospi(t) -sinpi(t) -1/tanpi(t)
58 * 011 (z+1-y)/4 sinpi(t) -cospi(t) -tanpi(t)
59 * 100 (y-z)/4 -sinpi(t) -cospi(t) tanpi(t)
60 * 101 (z+1-y)/4 -cospi(t) -sinpi(t) 1/tanpi(t)
61 * 110 (y-z)/4 -cospi(t) sinpi(t) -1/tanpi(t)
62 * 111 (z+1-y)/4 -sinpi(t) cospi(t) -tanpi(t)
63 * ---------------------------------------------------
65 * NOTE. This program compute sinpi/cospi(t<0.25) by __k_sin/cos(pi*t, 0.0).
66 * This will return a result with error slightly more than one ulp (but less
67 * than 2 ulp). If one wants accurate result, one may break up pi*t in
68 * high (tpi_h) and low (tpi_l) parts and call __k_sin/cos(tip_h, tip_lo)
73 #include "libm_protos.h"
74 #include "libm_macros.h"
76 #if defined(__SUNPRO_C)
81 pi
= 3.14159265358979323846, /* 400921FB,54442D18 */
82 sqrth_h
= 0.70710678118654757273731092936941422522068023681640625,
83 sqrth_l
= -4.8336466567264565185935844299127932213411660131004e-17;
87 sincospi(double x
, double *s
, double *c
) {
90 int hx
= ((int *) &x
)[HIWORD
];
91 unsigned h
, lx
= ((unsigned *) &x
)[LOWORD
];
93 ix
= hx
& ~0x80000000;
94 n
= (ix
>> 20) - 0x3ff;
95 if (n
>= 51) { /* |x| >= 2**51 */
97 #if defined(FPADD_TRAPS_INCOMPLETE_ON_NAN)
98 *s
= *c
= ix
>= 0x7ff80000 ? x
: x
- x
;
99 /* assumes sparc-like QNaN */
134 else if (n
< -2) /* |x| < 0.25 */
135 *s
= __k_sincos(pi
* fabs(x
), 0.0, c
);
137 /* y = |4x|, z = floor(y), and n = (int)(z mod 8.0) */
138 if (ix
< 0x41C00000) { /* |x| < 2**29 */
140 n
= (int) y
; /* exact */
145 else { /* 2**29 <= |x| < 2**51 */
150 ((unsigned *) &z
)[LOWORD
] = h
;
151 ((int *) &z
)[HIWORD
] = ix
;
155 if (k
) { /* x = N/4 */
157 *s
= *c
= sqrth_h
+ sqrth_l
;
167 y
= (n
& 2) == 0 ? 0.0 : 1.0;
170 if (((n
+ 1) & 4) != 0)
176 if (((n
+ (n
& 1)) & 2) == 0)
177 *s
= __k_sincos(pi
* t
, 0.0, c
);
179 *c
= __k_sincos(pi
* t
, 0.0, s
);
182 if (((n
+ 2) & 4) != 0)