5 * Common Development and Distribution License (the "License").
6 * You may not use this file except in compliance with the License.
7 *
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.
12 *
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]
18 *
19 * CDDL HEADER END
20 */
21
22 /*
23 * Copyright 2011 Nexenta Systems, Inc. All rights reserved.
24 */
25 /*
26 * Copyright 2006 Sun Microsystems, Inc. All rights reserved.
27 * Use is subject to license terms.
28 */
29
30 #pragma weak nearbyintf = __nearbyintf
31
32 #include "libm.h"
33 #include <fenv.h>
34
35 float
36 __nearbyintf(float x) {
37 union {
38 unsigned i;
39 float f;
40 } xx;
41 unsigned hx, sx, i, frac;
42 int rm;
43
44 xx.f = x;
45 sx = xx.i & 0x80000000;
46 hx = xx.i & ~0x80000000;
47
48 /* handle trivial cases */
49 if (hx >= 0x4b000000) { /* x is nan, inf, or already integral */
50 if (hx > 0x7f800000) /* x is nan */
51 return (x * x); /* + -> * for Cheetah */
52 return (x);
53 } else if (hx == 0) /* x is zero */
54 return (x);
55
56 /* get the rounding mode */
57 rm = fegetround();
58
59 /* flip the sense of directed roundings if x is negative */
60 if (sx && (rm == FE_UPWARD || rm == FE_DOWNWARD))
61 rm = (FE_UPWARD + FE_DOWNWARD) - rm;
62
63 /* handle |x| < 1 */
64 if (hx < 0x3f800000) {
65 if (rm == FE_UPWARD || (rm == FE_TONEAREST && hx > 0x3f000000))
66 xx.i = sx | 0x3f800000;
67 else
68 xx.i = sx;
69 return (xx.f);
70 }
71
72 /* round x at the integer bit */
73 i = 1 << (0x96 - (hx >> 23));
74 frac = hx & (i - 1);
75 if (!frac)
76 return (x);
77
78 hx &= ~(i - 1);
79 if (rm == FE_UPWARD || (rm == FE_TONEAREST && (frac > (i >> 1) ||
80 ((frac == (i >> 1)) && (hx & i)))))
81 xx.i = sx | (hx + i);
82 else
83 xx.i = sx | hx;
84 return (xx.f);
85 }
86
87 #if 0
88
89 /*
90 * Alternate implementations for SPARC, x86, using fp ops. These may
91 * be faster depending on how expensive saving and restoring the fp
92 * modes and status flags is.
93 */
94
95 #include "libm.h"
96 #include "fma.h"
97
98 #if defined(__sparc)
99
100 float
101 __nearbyintf(float x) {
102 union {
103 unsigned i;
104 float f;
105 } xx, yy;
106 float z;
107 unsigned hx, sx, fsr, oldfsr;
108 int rm;
109
110 xx.f = x;
111 sx = xx.i & 0x80000000;
112 hx = xx.i & ~0x80000000;
113
114 /* handle trivial cases */
115 if (hx >= 0x4b000000) /* x is nan, inf, or already integral */
116 return (x + 0.0f);
117 else if (hx == 0) /* x is zero */
118 return (x);
119
120 /* save the fsr */
121 __fenv_getfsr(&oldfsr);
122
123 /* handle |x| < 1 */
124 if (hx < 0x3f800000) {
125 /* flip the sense of directed roundings if x is negative */
126 rm = oldfsr >> 30;
127 if (sx)
128 rm ^= rm >> 1;
129 if (rm == FSR_RP || (rm == FSR_RN && hx > 0x3f000000))
130 xx.i = sx | 0x3f800000;
131 else
132 xx.i = sx;
133 return (xx.f);
134 }
135
136 /* clear the inexact trap */
137 fsr = oldfsr & ~FSR_NXM;
138 __fenv_setfsr(&fsr);
139
140 /* round x at the integer bit */
141 yy.i = sx | 0x4b000000;
142 z = (x + yy.f) - yy.f;
143
144 /* restore the old fsr */
145 __fenv_setfsr(&oldfsr);
146
147 return (z);
148 }
149
150 #elif defined(__x86)
151
152 /* inline template */
153 extern long double frndint(long double);
154
155 float
156 __nearbyintf(float x) {
157 long double z;
158 unsigned oldcwsw, cwsw;
159
160 /* save the control and status words, mask the inexact exception */
161 __fenv_getcwsw(&oldcwsw);
162 cwsw = oldcwsw | 0x00200000;
163 __fenv_setcwsw(&cwsw);
164
165 z = frndint((long double) x);
166
167 /*
168 * restore the control and status words, preserving all but the
169 * inexact flag
170 */
171 __fenv_getcwsw(&cwsw);
172 oldcwsw |= (cwsw & 0x1f);
173 __fenv_setcwsw(&oldcwsw);
174
175 /* note: the value of z is representable in single precision */
176 return (z);
177 }
178
179 #else
180 #error Unknown architecture
181 #endif
182
183 #endif
|
5 * Common Development and Distribution License (the "License").
6 * You may not use this file except in compliance with the License.
7 *
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.
12 *
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]
18 *
19 * CDDL HEADER END
20 */
21
22 /*
23 * Copyright 2011 Nexenta Systems, Inc. All rights reserved.
24 */
25
26 /*
27 * Copyright 2006 Sun Microsystems, Inc. All rights reserved.
28 * Use is subject to license terms.
29 */
30
31 #pragma weak nearbyintf = __nearbyintf
32
33 #include "libm.h"
34 #include <fenv.h>
35
36 float
37 __nearbyintf(float x)
38 {
39 union {
40 unsigned i;
41 float f;
42 } xx;
43
44 unsigned hx, sx, i, frac;
45 int rm;
46
47 xx.f = x;
48 sx = xx.i & 0x80000000;
49 hx = xx.i & ~0x80000000;
50
51 /* handle trivial cases */
52 if (hx >= 0x4b000000) { /* x is nan, inf, or already integral */
53 if (hx > 0x7f800000) /* x is nan */
54 return (x * x); /* + -> * for Cheetah */
55
56 return (x);
57 } else if (hx == 0) { /* x is zero */
58 return (x);
59 }
60
61 /* get the rounding mode */
62 rm = fegetround();
63
64 /* flip the sense of directed roundings if x is negative */
65 if (sx && (rm == FE_UPWARD || rm == FE_DOWNWARD))
66 rm = (FE_UPWARD + FE_DOWNWARD) - rm;
67
68 /* handle |x| < 1 */
69 if (hx < 0x3f800000) {
70 if (rm == FE_UPWARD || (rm == FE_TONEAREST && hx > 0x3f000000))
71 xx.i = sx | 0x3f800000;
72 else
73 xx.i = sx;
74
75 return (xx.f);
76 }
77
78 /* round x at the integer bit */
79 i = 1 << (0x96 - (hx >> 23));
80 frac = hx & (i - 1);
81
82 if (!frac)
83 return (x);
84
85 hx &= ~(i - 1);
86
87 if (rm == FE_UPWARD || (rm == FE_TONEAREST && (frac > (i >> 1) ||
88 ((frac == (i >> 1)) && (hx & i)))))
89 xx.i = sx | (hx + i);
90 else
91 xx.i = sx | hx;
92
93 return (xx.f);
94 }
95
96 #if 0
97 /*
98 * Alternate implementations for SPARC, x86, using fp ops. These may
99 * be faster depending on how expensive saving and restoring the fp
100 * modes and status flags is.
101 */
102
103 #include "libm.h"
104 #include "fma.h"
105
106 #if defined(__sparc)
107 float
108 __nearbyintf(float x)
109 {
110 union {
111 unsigned i;
112 float f;
113 } xx, yy;
114
115 float z;
116 unsigned hx, sx, fsr, oldfsr;
117 int rm;
118
119 xx.f = x;
120 sx = xx.i & 0x80000000;
121 hx = xx.i & ~0x80000000;
122
123 /* handle trivial cases */
124 if (hx >= 0x4b000000) /* x is nan, inf, or already integral */
125 return (x + 0.0f);
126 else if (hx == 0) /* x is zero */
127 return (x);
128
129 /* save the fsr */
130 __fenv_getfsr(&oldfsr);
131
132 /* handle |x| < 1 */
133 if (hx < 0x3f800000) {
134 /* flip the sense of directed roundings if x is negative */
135 rm = oldfsr >> 30;
136
137 if (sx)
138 rm ^= rm >> 1;
139
140 if (rm == FSR_RP || (rm == FSR_RN && hx > 0x3f000000))
141 xx.i = sx | 0x3f800000;
142 else
143 xx.i = sx;
144
145 return (xx.f);
146 }
147
148 /* clear the inexact trap */
149 fsr = oldfsr & ~FSR_NXM;
150 __fenv_setfsr(&fsr);
151
152 /* round x at the integer bit */
153 yy.i = sx | 0x4b000000;
154 z = (x + yy.f) - yy.f;
155
156 /* restore the old fsr */
157 __fenv_setfsr(&oldfsr);
158
159 return (z);
160 }
161 #elif defined(__x86)
162 /* inline template */
163 extern long double frndint(long double);
164
165 float
166 __nearbyintf(float x)
167 {
168 long double z;
169 unsigned oldcwsw, cwsw;
170
171 /* save the control and status words, mask the inexact exception */
172 __fenv_getcwsw(&oldcwsw);
173 cwsw = oldcwsw | 0x00200000;
174 __fenv_setcwsw(&cwsw);
175
176 z = frndint((long double)x);
177
178 /*
179 * restore the control and status words, preserving all but the
180 * inexact flag
181 */
182 __fenv_getcwsw(&cwsw);
183 oldcwsw |= (cwsw & 0x1f);
184 __fenv_setcwsw(&oldcwsw);
185
186 /* note: the value of z is representable in single precision */
187 return (z);
188 }
189 #else
190 #error Unknown architecture
191 #endif
192 #endif
|