|
|
1
|
+/* s_sincosf.c -- float version of s_sincos.c
|
|
|
2
|
+ *
|
|
|
3
|
+ * ====================================================
|
|
|
4
|
+ * This file is derived from fdlibm:
|
|
|
5
|
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
|
|
|
6
|
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
|
|
|
7
|
+ * Permission to use, copy, modify, and distribute this
|
|
|
8
|
+ * software is freely granted, provided that this notice
|
|
|
9
|
+ * is preserved.
|
|
|
10
|
+ *
|
|
|
11
|
+ * ====================================================
|
|
|
12
|
+ * Copyright (C) 2013 Elliot Saba
|
|
|
13
|
+ * Developed at the University of Washington
|
|
|
14
|
+ *
|
|
|
15
|
+ * Permission to use, copy, modify, and distribute this
|
|
|
16
|
+ * software is freely granted, provided that this notice
|
|
|
17
|
+ * is preserved.
|
|
|
18
|
+ * ====================================================
|
|
|
19
|
+*/
|
|
|
20
|
+
|
|
|
21
|
+#include "cdefs-compat.h"
|
|
|
22
|
+
|
|
|
23
|
+#include <float.h>
|
|
|
24
|
+#include <openlibm_math.h>
|
|
|
25
|
+
|
|
|
26
|
+//#define INLINE_KERNEL_COSDF
|
|
|
27
|
+//#define INLINE_KERNEL_SINDF
|
|
|
28
|
+//#define INLINE_REM_PIO2F
|
|
|
29
|
+#include "math_private.h"
|
|
|
30
|
+//#include "e_rem_pio2f.c"
|
|
|
31
|
+//#include "k_cosf.c"
|
|
|
32
|
+//#include "k_sinf.c"
|
|
|
33
|
+
|
|
|
34
|
+
|
|
|
35
|
+/* Constants used in shortcircuits in sincosf */
|
|
|
36
|
+static const double
|
|
|
37
|
+sc1pio2 = 1*M_PI_2, /* 0x3FF921FB, 0x54442D18 */
|
|
|
38
|
+sc2pio2 = 2*M_PI_2, /* 0x400921FB, 0x54442D18 */
|
|
|
39
|
+sc3pio2 = 3*M_PI_2, /* 0x4012D97C, 0x7F3321D2 */
|
|
|
40
|
+sc4pio2 = 4*M_PI_2, /* 0x401921FB, 0x54442D18 */
|
|
|
41
|
+
|
|
|
42
|
+/* Constants used in polynomial approximation of sin/cos */
|
|
|
43
|
+one = 1.0,
|
|
|
44
|
+S1 = -0x15555554cbac77.0p-55, /* -0.166666666416265235595 */
|
|
|
45
|
+S2 = 0x111110896efbb2.0p-59, /* 0.0083333293858894631756 */
|
|
|
46
|
+S3 = -0x1a00f9e2cae774.0p-65, /* -0.000198393348360966317347 */
|
|
|
47
|
+S4 = 0x16cd878c3b46a7.0p-71, /* 0.0000027183114939898219064 */
|
|
|
48
|
+C0 = -0x1ffffffd0c5e81.0p-54, /* -0.499999997251031003120 */
|
|
|
49
|
+C1 = 0x155553e1053a42.0p-57, /* 0.0416666233237390631894 */
|
|
|
50
|
+C2 = -0x16c087e80f1e27.0p-62, /* -0.00138867637746099294692 */
|
|
|
51
|
+C3 = 0x199342e0ee5069.0p-68; /* 0.0000243904487962774090654 */
|
|
|
52
|
+
|
|
|
53
|
+static void
|
|
|
54
|
+__kernel_sincosdf( double x, float * s, float * c )
|
|
|
55
|
+{
|
|
|
56
|
+ double r, w, z, v;
|
|
|
57
|
+ z = x*x;
|
|
|
58
|
+ w = z*z;
|
|
|
59
|
+
|
|
|
60
|
+ /* cos-specific computation; equivalent to calling
|
|
|
61
|
+ __kernel_cos(x,y) and storing in k_c*/
|
|
|
62
|
+ r = C2+z*C3;
|
|
|
63
|
+ double k_c = ((one+z*C0) + w*C1) + (w*z)*r;
|
|
|
64
|
+
|
|
|
65
|
+ /* sin-specific computation; equivalent to calling
|
|
|
66
|
+ __kernel_sin(x,y,1) and storing in k_s*/
|
|
|
67
|
+ r = S3+z*S4;
|
|
|
68
|
+ v = z*x;
|
|
|
69
|
+ double k_s = (x + v*(S1+z*S2)) + v*w*r;
|
|
|
70
|
+
|
|
|
71
|
+ *c = k_c;
|
|
|
72
|
+ *s = k_s;
|
|
|
73
|
+}
|
|
|
74
|
+
|
|
|
75
|
+OLM_DLLEXPORT void
|
|
|
76
|
+sincosf(float x, float * s, float * c) {
|
|
|
77
|
+ // Worst approximation of sin and cos NA
|
|
|
78
|
+ *s = x;
|
|
|
79
|
+ *c = x;
|
|
|
80
|
+
|
|
|
81
|
+ double y;
|
|
|
82
|
+ float k_c, k_s;
|
|
|
83
|
+ int32_t n, hx, ix;
|
|
|
84
|
+
|
|
|
85
|
+ GET_FLOAT_WORD(hx,x);
|
|
|
86
|
+ ix = hx & 0x7fffffff;
|
|
|
87
|
+
|
|
|
88
|
+ if(ix <= 0x3f490fda) { /* |x| ~<= pi/4 */
|
|
|
89
|
+ if(ix<0x39800000) { /* |x| < 2**-12 */
|
|
|
90
|
+ /* Check if x is exactly zero */
|
|
|
91
|
+ if(((int)x)==0) {
|
|
|
92
|
+ *s = x;
|
|
|
93
|
+ *c = 1.0f;
|
|
|
94
|
+ return;
|
|
|
95
|
+ }
|
|
|
96
|
+ }
|
|
|
97
|
+ __kernel_sincosdf(x, s, c);
|
|
|
98
|
+ return;
|
|
|
99
|
+ }
|
|
|
100
|
+ /* |x| ~<= 5*pi/4 */
|
|
|
101
|
+ if (ix<=0x407b53d1) {
|
|
|
102
|
+ /* |x| ~<= 3pi/4 */
|
|
|
103
|
+ if(ix<=0x4016cbe3) {
|
|
|
104
|
+ if(hx>0) {
|
|
|
105
|
+ __kernel_sincosdf( sc1pio2 - x, c, s );
|
|
|
106
|
+ }
|
|
|
107
|
+ else {
|
|
|
108
|
+ __kernel_sincosdf( sc1pio2 + x, c, &k_s );
|
|
|
109
|
+ *s = -k_s;
|
|
|
110
|
+ }
|
|
|
111
|
+ } else {
|
|
|
112
|
+
|
|
|
113
|
+ if(hx>0) {
|
|
|
114
|
+ __kernel_sincosdf( sc2pio2 - x, s, &k_c );
|
|
|
115
|
+ *c = -k_c;
|
|
|
116
|
+ } else {
|
|
|
117
|
+ __kernel_sincosdf( -sc2pio2 - x, s, &k_c );
|
|
|
118
|
+ *c = -k_c;
|
|
|
119
|
+ }
|
|
|
120
|
+ }
|
|
|
121
|
+ return;
|
|
|
122
|
+ }
|
|
|
123
|
+
|
|
|
124
|
+ /* |x| ~<= 9*pi/4 */
|
|
|
125
|
+ if(ix<=0x40e231d5) {
|
|
|
126
|
+ /* |x| ~> 7*pi/4 */
|
|
|
127
|
+ if(ix<=0x40afeddf) {
|
|
|
128
|
+ if(hx>0) {
|
|
|
129
|
+ __kernel_sincosdf( x - sc3pio2, c, &k_s );
|
|
|
130
|
+ *s = -k_s;
|
|
|
131
|
+ } else {
|
|
|
132
|
+ __kernel_sincosdf( x + sc3pio2, &k_c, s );
|
|
|
133
|
+ *c = -k_c;
|
|
|
134
|
+ }
|
|
|
135
|
+ }
|
|
|
136
|
+ else {
|
|
|
137
|
+ if( hx > 0 ) {
|
|
|
138
|
+ __kernel_sincosdf( x - sc4pio2, s, c );
|
|
|
139
|
+ } else {
|
|
|
140
|
+ __kernel_sincosdf( x + sc4pio2, s, c );
|
|
|
141
|
+ }
|
|
|
142
|
+ }
|
|
|
143
|
+ return;
|
|
|
144
|
+ }
|
|
|
145
|
+
|
|
|
146
|
+ /* cos(Inf or NaN) is NaN */
|
|
|
147
|
+ else if(ix>=0x7f800000) {
|
|
|
148
|
+ *c = *s = x-x;
|
|
|
149
|
+ } else {
|
|
|
150
|
+ /* general argument reduction needed */
|
|
|
151
|
+ n = __ieee754_rem_pio2f(x,&y);
|
|
|
152
|
+
|
|
|
153
|
+ switch(n&3) {
|
|
|
154
|
+ case 0:
|
|
|
155
|
+ __kernel_sincosdf( y, s, c );
|
|
|
156
|
+ break;
|
|
|
157
|
+ case 1:
|
|
|
158
|
+ __kernel_sincosdf( -y, c, s );
|
|
|
159
|
+ break;
|
|
|
160
|
+ case 2:
|
|
|
161
|
+ __kernel_sincosdf( -y, s, &k_c);
|
|
|
162
|
+ *c = -k_c;
|
|
|
163
|
+ break;
|
|
|
164
|
+ default:
|
|
|
165
|
+ __kernel_sincosdf( -y, &k_c, &k_s );
|
|
|
166
|
+ *c = -k_c;
|
|
|
167
|
+ *s = -k_s;
|
|
|
168
|
+ break;
|
|
|
169
|
+ }
|
|
|
170
|
+ }
|
|
|
171
|
+
|
|
|
172
|
+} |