xref: /illumos-gate/usr/src/lib/libm/common/R/rintf.c (revision b515258426fed6c7311fd3f1dea697cfbd4085c6)
1 /*
2  * CDDL HEADER START
3  *
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.
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 #if defined(ELFOBJ)
31 #pragma weak aintf = __aintf
32 #pragma weak anintf = __anintf
33 #pragma weak irintf = __irintf
34 #pragma weak nintf = __nintf
35 #pragma weak rintf = __rintf
36 #endif
37 
38 /* INDENT OFF */
39 /*
40  * aintf(x)	return x chopped to integral value
41  * anintf(x)	return sign(x)*(|x|+0.5) chopped to integral value
42  * irintf(x)	return rint(x) in integer format
43  * nintf(x)	return anint(x) in integer format
44  * rintf(x)	return x rounded to integral according to the rounding direction
45  *
46  * NOTE: rintf(x), aintf(x) and anintf(x) return results with the same sign as
47  * x's,  including 0.0.
48  */
49 
50 #include "libm.h"
51 
52 static const float xf[] = {
53 /* ZEROF */	0.0f,
54 /* TWO_23F */	8.3886080000e6f,
55 /* MTWO_23F */	-8.3886080000e6f,
56 /* ONEF */	1.0f,
57 /* MONEF */	-1.0f,
58 /* HALFF */	0.5f,
59 /* MHALFF */	-0.5f,
60 /* HUGEF */	1.0e30f,
61 };
62 
63 #define	ZEROF		xf[0]
64 #define	TWO_23F		xf[1]
65 #define	MTWO_23F	xf[2]
66 #define	ONEF		xf[3]
67 #define	MONEF		xf[4]
68 #define	HALFF		xf[5]
69 #define	MHALFF		xf[6]
70 #define	HUGEF		xf[7]
71 /* INDENT ON */
72 
73 float
74 aintf(float x) {
75 	int hx, k;
76 	float y;
77 
78 	hx = *(int *) &x;
79 	k = (hx & ~0x80000000) >> 23;
80 	if (k < 150) {
81 		y = (float) ((int) x);
82 		/*
83 		 * make sure y has the same sign of x when |x|<0.5
84 		 * (i.e., y=0.0)
85 		 */
86 		return (((k - 127) & hx) < 0 ? -y : y);
87 	} else
88 		/* signal invalid if x is a SNaN */
89 		return (x * ONEF);		/* +0 -> *1 for Cheetah */
90 }
91 
92 float
93 anintf(float x) {
94 	volatile float dummy;
95 	int hx, k, j, ix;
96 
97 	hx = *(int *) &x;
98 	ix = hx & ~0x80000000;
99 	k = ix >> 23;
100 	if (((k - 127) ^ (k - 150)) < 0) {
101 		j = 1 << (149 - k);
102 		k = j + j - 1;
103 		if ((k & hx) != 0)
104 			dummy = HUGEF + x;	/* raise inexact */
105 		*(int *) &x = (hx + j) & ~k;
106 		return (x);
107 	} else if (k <= 126) {
108 		dummy = HUGEF + x;
109 		*(int *) &x = (0x3f800000 & ((125 - k) >> 31)) |
110 			(0x80000000 & hx);
111 		return (x);
112 	} else
113 		/* signal invalid if x is a SNaN */
114 		return (x * ONEF);		/* +0 -> *1 for Cheetah */
115 }
116 
117 int
118 irintf(float x) {
119 	float v;
120 	int hx, k;
121 
122 	hx = *(int *) &x;
123 	k = (hx & ~0x80000000) >> 23;
124 	v = xf[((k - 150) >> 31) & (1 - (hx >> 31))];
125 	return ((int) ((float) (x + v) - v));
126 }
127 
128 int
129 nintf(float x) {
130 	int hx, ix, k, j, m;
131 	volatile float dummy;
132 
133 	hx = *(int *) &x;
134 	k = (hx & ~0x80000000) >> 23;
135 	if (((k - 126) ^ (k - 150)) < 0) {
136 		ix = (hx & 0x00ffffff) | 0x800000;
137 		m = 149 - k;
138 		j = 1 << m;
139 		if ((ix & (j + j - 1)) != 0)
140 			dummy = HUGEF + x;
141 		hx = hx >> 31;
142 		return ((((ix + j) >> (m + 1)) ^ hx) - hx);
143 	} else
144 		return ((int) x);
145 }
146 
147 float
148 rintf(float x) {
149 	float w, v;
150 	int hx, k;
151 
152 	hx = *(int *) &x;
153 	k = (hx & ~0x80000000) >> 23;
154 #if defined(FPADD_TRAPS_INCOMPLETE_ON_NAN)
155 	if (k >= 150)
156 		return (x * ONEF);
157 	v = xf[1 - (hx >> 31)];
158 #else
159 	v = xf[((k - 150) >> 31) & (1 - (hx >> 31))];
160 #endif
161 	w = (float) (x + v);
162 	if (k < 127 && w == v)
163 		return (ZEROF * x);
164 	else
165 		return (w - v);
166 }
167