xref: /illumos-gate/usr/src/lib/libm/common/m9x/nexttoward.c (revision 1ec68d336ba97cd53f46053ac10401d16014d075)
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 #pragma weak nexttoward = __nexttoward
31 
32 /*
33  * nexttoward(x, y) delivers the next representable number after x
34  * in the direction of y.  If x and y are both zero, the result is
35  * zero with the same sign as y.  If either x or y is NaN, the result
36  * is NaN.
37  *
38  * If x != y and the result is infinite, overflow is raised; if
39  * x != y and the result is subnormal or zero, underflow is raised.
40  * (This is wrong, but it's what C99 apparently wants.)
41  */
42 
43 #include "libm.h"
44 
45 #if defined(__sparc)
46 
47 static union {
48 	unsigned i[2];
49 	double d;
50 } C[] = {
51 	0x00100000, 0,
52 	0x7fe00000, 0,
53 	0x7fffffff, 0xffffffff
54 };
55 
56 #define	tiny	C[0].d
57 #define	huge	C[1].d
58 #define	qnan	C[2].d
59 
60 enum fcc_type {
61 	fcc_equal = 0,
62 	fcc_less = 1,
63 	fcc_greater = 2,
64 	fcc_unordered = 3
65 };
66 
67 #ifdef __sparcv9
68 #define	_Q_cmp	_Qp_cmp
69 #endif
70 
71 extern enum fcc_type _Q_cmp(const long double *, const long double *);
72 
73 double
__nexttoward(double x,long double y)74 __nexttoward(double x, long double y) {
75 	union {
76 		unsigned i[2];
77 		double d;
78 	} xx;
79 	union {
80 		unsigned i[4];
81 		long double q;
82 	} yy;
83 	long double lx;
84 	unsigned hx;
85 	volatile double	dummy;
86 	enum fcc_type rel;
87 
88 	/*
89 	 * It would be somewhat more efficient to check for NaN and
90 	 * zero operands before converting x to long double and then
91 	 * to code the comparison in line rather than calling _Q_cmp.
92 	 * However, since this code probably won't get used much,
93 	 * I'm opting in favor of simplicity instead.
94 	 */
95 	lx = xx.d = x;
96 	hx = (xx.i[0] & ~0x80000000) | xx.i[1];
97 
98 	/* check for each of four possible orderings */
99 	rel = _Q_cmp(&lx, &y);
100 	if (rel == fcc_unordered)
101 		return (qnan);
102 
103 	if (rel == fcc_equal) {
104 		if (hx == 0) {	/* x is zero; return zero with y's sign */
105 			yy.q = y;
106 			xx.i[0] = yy.i[0];
107 			return (xx.d);
108 		}
109 		return (x);
110 	}
111 
112 	if (rel == fcc_less) {
113 		if (hx == 0) {	/* x is zero */
114 			xx.i[0] = 0;
115 			xx.i[1] = 0x00000001;
116 		} else if ((int)xx.i[0] >= 0) {	/* x is positive */
117 			if (++xx.i[1] == 0)
118 				xx.i[0]++;
119 		} else {
120 			if (xx.i[1]-- == 0)
121 				xx.i[0]--;
122 		}
123 	} else {
124 		if (hx == 0) {	/* x is zero */
125 			xx.i[0] = 0x80000000;
126 			xx.i[1] = 0x00000001;
127 		} else if ((int)xx.i[0] >= 0) {	/* x is positive */
128 			if (xx.i[1]-- == 0)
129 				xx.i[0]--;
130 		} else {
131 			if (++xx.i[1] == 0)
132 				xx.i[0]++;
133 		}
134 	}
135 
136 	/* raise exceptions as needed */
137 	hx = xx.i[0] & ~0x80000000;
138 	if (hx == 0x7ff00000) {
139 		dummy = huge;
140 		dummy *= huge;
141 	} else if (hx < 0x00100000) {
142 		dummy = tiny;
143 		dummy *= tiny;
144 	}
145 
146 	return (xx.d);
147 }
148 
149 #elif defined(__x86)
150 
151 static union {
152 	unsigned i[2];
153 	double d;
154 } C[] = {
155 	0, 0x00100000,
156 	0, 0x7fe00000,
157 };
158 
159 #define	tiny	C[0].d
160 #define	huge	C[1].d
161 
162 double
__nexttoward(double x,long double y)163 __nexttoward(double x, long double y) {
164 	union {
165 		unsigned i[2];
166 		double d;
167 	} xx;
168 	unsigned hx;
169 	long double lx;
170 	volatile double	dummy;
171 
172 	lx = xx.d = x;
173 	hx = (xx.i[1] & ~0x80000000) | xx.i[0];
174 
175 	/* check for each of four possible orderings */
176 	if (isunordered(lx, y))
177 		return ((double) (lx + y));
178 
179 	if (lx == y)
180 		return ((double) y);
181 
182 	if (lx < y) {
183 		if (hx == 0) {	/* x is zero */
184 			xx.i[0] = 0x00000001;
185 			xx.i[1] = 0;
186 		} else if ((int)xx.i[1] >= 0) {	/* x is positive */
187 			if (++xx.i[0] == 0)
188 				xx.i[1]++;
189 		} else {
190 			if (xx.i[0]-- == 0)
191 				xx.i[1]--;
192 		}
193 	} else {
194 		if (hx == 0) {	/* x is zero */
195 			xx.i[0] = 0x00000001;
196 			xx.i[1] = 0x80000000;
197 		} else if ((int)xx.i[1] >= 0) {	/* x is positive */
198 			if (xx.i[0]-- == 0)
199 				xx.i[1]--;
200 		} else {
201 			if (++xx.i[0] == 0)
202 				xx.i[1]++;
203 		}
204 	}
205 
206 	/* raise exceptions as needed */
207 	hx = xx.i[1] & ~0x80000000;
208 	if (hx == 0x7ff00000) {
209 		dummy = huge;
210 		dummy *= huge;
211 	} else if (hx < 0x00100000) {
212 		dummy = tiny;
213 		dummy *= tiny;
214 	}
215 
216 	return (xx.d);
217 }
218 
219 #else
220 #error Unknown architecture
221 #endif
222