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 2008 Sun Microsystems, Inc. All rights reserved.
24 * Use is subject to license terms.
25 */
26
27 /*
28 * Short cut for conversion from double precision to decimal
29 * floating point
30 */
31
32 #include "lint.h"
33 #include <sys/types.h>
34 #include <sys/isa_defs.h>
35 #include "base_conversion.h"
36
37 /*
38 * Powers of ten rounded up. If i is the largest index such that
39 * tbl_decade[i] <= x, then:
40 *
41 * if i == 0 then x < 10^-49
42 * else if i == TBL_DECADE_MAX then x >= 10^67
43 * else 10^(i-TBL_DECADE_OFFSET) <= x < 10^(i-TBL_DECADE_OFFSET+1)
44 */
45
46 #define TBL_DECADE_OFFSET 50
47 #define TBL_DECADE_MAX 117
48
49 static const double tbl_decade[TBL_DECADE_MAX + 1] = {
50 0.0,
51 1.00000000000000012631e-49, 1.00000000000000012631e-48,
52 1.00000000000000009593e-47, 1.00000000000000002300e-46,
53 1.00000000000000013968e-45, 1.00000000000000007745e-44,
54 1.00000000000000007745e-43, 1.00000000000000003762e-42,
55 1.00000000000000000576e-41, 1.00000000000000013321e-40,
56 1.00000000000000009243e-39, 1.00000000000000009243e-38,
57 1.00000000000000006632e-37, 1.00000000000000010809e-36,
58 1.00000000000000000786e-35, 1.00000000000000014150e-34,
59 1.00000000000000005597e-33, 1.00000000000000005597e-32,
60 1.00000000000000008334e-31, 1.00000000000000008334e-30,
61 1.00000000000000008334e-29, 1.00000000000000008334e-28,
62 1.00000000000000003849e-27, 1.00000000000000003849e-26,
63 1.00000000000000003849e-25, 1.00000000000000010737e-24,
64 1.00000000000000010737e-23, 1.00000000000000004860e-22,
65 1.00000000000000009562e-21, 1.00000000000000009562e-20,
66 1.00000000000000009562e-19, 1.00000000000000007154e-18,
67 1.00000000000000007154e-17, 1.00000000000000010236e-16,
68 1.00000000000000007771e-15, 1.00000000000000015659e-14,
69 1.00000000000000003037e-13, 1.00000000000000018184e-12,
70 1.00000000000000010106e-11, 1.00000000000000003643e-10,
71 1.00000000000000006228e-09, 1.00000000000000002092e-08,
72 1.00000000000000008710e-07, 1.00000000000000016651e-06,
73 1.00000000000000008180e-05, 1.00000000000000004792e-04,
74 1.00000000000000002082e-03, 1.00000000000000002082e-02,
75 1.00000000000000005551e-01, 1.00000000000000000000e+00,
76 1.00000000000000000000e+01, 1.00000000000000000000e+02,
77 1.00000000000000000000e+03, 1.00000000000000000000e+04,
78 1.00000000000000000000e+05, 1.00000000000000000000e+06,
79 1.00000000000000000000e+07, 1.00000000000000000000e+08,
80 1.00000000000000000000e+09, 1.00000000000000000000e+10,
81 1.00000000000000000000e+11, 1.00000000000000000000e+12,
82 1.00000000000000000000e+13, 1.00000000000000000000e+14,
83 1.00000000000000000000e+15, 1.00000000000000000000e+16,
84 1.00000000000000000000e+17, 1.00000000000000000000e+18,
85 1.00000000000000000000e+19, 1.00000000000000000000e+20,
86 1.00000000000000000000e+21, 1.00000000000000000000e+22,
87 1.00000000000000008389e+23, 1.00000000000000011744e+24,
88 1.00000000000000009060e+25, 1.00000000000000004765e+26,
89 1.00000000000000001329e+27, 1.00000000000000017821e+28,
90 1.00000000000000009025e+29, 1.00000000000000001988e+30,
91 1.00000000000000007618e+31, 1.00000000000000005366e+32,
92 1.00000000000000008969e+33, 1.00000000000000006087e+34,
93 1.00000000000000015310e+35, 1.00000000000000004242e+36,
94 1.00000000000000007194e+37, 1.00000000000000016638e+38,
95 1.00000000000000009082e+39, 1.00000000000000003038e+40,
96 1.00000000000000000620e+41, 1.00000000000000004489e+42,
97 1.00000000000000001394e+43, 1.00000000000000008821e+44,
98 1.00000000000000008821e+45, 1.00000000000000011990e+46,
99 1.00000000000000004385e+47, 1.00000000000000004385e+48,
100 1.00000000000000007630e+49, 1.00000000000000007630e+50,
101 1.00000000000000015937e+51, 1.00000000000000012614e+52,
102 1.00000000000000020590e+53, 1.00000000000000007829e+54,
103 1.00000000000000001024e+55, 1.00000000000000009190e+56,
104 1.00000000000000004835e+57, 1.00000000000000008319e+58,
105 1.00000000000000008319e+59, 1.00000000000000012779e+60,
106 1.00000000000000009211e+61, 1.00000000000000003502e+62,
107 1.00000000000000005786e+63, 1.00000000000000002132e+64,
108 1.00000000000000010901e+65, 1.00000000000000013239e+66,
109 1.00000000000000013239e+67
110 };
111
112 /*
113 * Convert a positive double precision integer x <= 2147483647999999744
114 * (the largest double less than 2^31 * 10^9; this implementation works
115 * up to the largest double less than 2^25 * 10^12) to a string of ASCII
116 * decimal digits, adding leading zeroes so that the result has at least
117 * n digits. The string is terminated by a null byte, and its length
118 * is returned.
119 *
120 * This routine assumes round-to-nearest mode is in effect and any
121 * exceptions raised will be ignored.
122 */
123
124 #define tenm4 tbl_decade[TBL_DECADE_OFFSET - 4]
125 #define ten4 tbl_decade[TBL_DECADE_OFFSET + 4]
126 #define tenm12 tbl_decade[TBL_DECADE_OFFSET - 12]
127 #define ten12 tbl_decade[TBL_DECADE_OFFSET + 12]
128 #define one tbl_decade[TBL_DECADE_OFFSET]
129
130 static int
__double_to_digits(double x,char * s,int n)131 __double_to_digits(double x, char *s, int n)
132 {
133 double y;
134 int d[5], i, j;
135 char *ss, tmp[4];
136
137 /* decompose x into four-digit chunks */
138 y = (int)(x * tenm12);
139 x -= y * ten12;
140 if (x < 0.0) {
141 y -= one;
142 x += ten12;
143 }
144 d[0] = (int)(y * tenm4);
145 d[1] = (int)(y - d[0] * ten4);
146 y = (int)(x * tenm4);
147 d[4] = (int)(x - y * ten4);
148 d[2] = (int)(y * tenm4);
149 d[3] = (int)(y - d[2] * ten4);
150
151 /*
152 * Find the first nonzero chunk or the point at which to start
153 * converting so we have n digits, whichever comes first.
154 */
155 ss = s;
156 if (n > 20) {
157 for (j = 0; j < n - 20; j++)
158 *ss++ = '0';
159 i = 0;
160 } else {
161 for (i = 0; d[i] == 0 && n <= ((4 - i) << 2); i++)
162 ;
163 __four_digits_quick(d[i], tmp);
164 for (j = 0; j < 4 && tmp[j] == '0' &&
165 n <= ((4 - i) << 2) + 3 - j; j++)
166 ;
167 while (j < 4)
168 *ss++ = tmp[j++];
169 i++;
170 }
171
172 /* continue converting four-digit chunks */
173 while (i < 5) {
174 __four_digits_quick(d[i], ss);
175 ss += 4;
176 i++;
177 }
178
179 *ss = '\0';
180 return (ss - s);
181 }
182
183 /*
184 * Round a positive double precision number *x to the nearest integer,
185 * returning the result and passing back an indication of accuracy in
186 * *pe. On entry, nrx is the number of rounding errors already com-
187 * mitted in forming *x. On exit, *pe is 0 if *x was already integral
188 * and exact, 1 if the result is the correctly rounded integer value
189 * but not exact, and 2 if error in *x precludes determining the cor-
190 * rectly rounded integer value (i.e., the error might be larger than
191 * 1/2 - |*x - rx|, where rx is the nearest integer to *x).
192 */
193
194 static union {
195 unsigned int i[2];
196 double d;
197 } C[] = {
198 #ifdef _LITTLE_ENDIAN
199 { 0x00000000, 0x43300000 },
200 { 0x00000000, 0x3ca00000 },
201 { 0x00000000, 0x3fe00000 },
202 { 0xffffffff, 0x3fdfffff },
203 #else
204 { 0x43300000, 0x00000000 },
205 { 0x3ca00000, 0x00000000 },
206 { 0x3fe00000, 0x00000000 },
207 { 0x3fdfffff, 0xffffffff }, /* nextafter(1/2, 0) */
208 #endif
209 };
210
211 #define two52 C[0].d
212 #define twom53 C[1].d
213 #define half C[2].d
214 #define halfdec C[3].d
215
216 static double
__arint_set_n(double * x,int nrx,int * pe)217 __arint_set_n(double *x, int nrx, int *pe)
218 {
219 int hx;
220 double rx, rmx;
221
222 #ifdef _LITTLE_ENDIAN
223 hx = *(1+(int *)x);
224 #else
225 hx = *(int *)x;
226 #endif
227 if (hx >= 0x43300000) {
228 /* x >= 2^52, so it's already integral */
229 if (nrx == 0)
230 *pe = 0;
231 else if (nrx == 1 && hx < 0x43400000)
232 *pe = 1;
233 else
234 *pe = 2;
235 return (*x);
236 } else if (hx < 0x3fe00000) {
237 /* x < 1/2 */
238 if (nrx > 1 && hx == 0x3fdfffff)
239 *pe = (*x == halfdec)? 2 : 1;
240 else
241 *pe = 1;
242 return (0.0);
243 }
244
245 rx = (*x + two52) - two52;
246 if (nrx == 0) {
247 *pe = (rx == *x)? 0 : 1;
248 } else {
249 rmx = rx - *x;
250 if (rmx < 0.0)
251 rmx = -rmx;
252 *pe = (nrx * twom53 * *x < half - rmx)? 1 : 2;
253 }
254 return (rx);
255 }
256
257 /*
258 * Attempt to convert dd to a decimal record *pd according to the
259 * modes in *pm using double precision floating point. Return zero
260 * and sets *ps to reflect any exceptions incurred if successful.
261 * Return a nonzero value if unsuccessful.
262 */
263 int
__fast_double_to_decimal(double * dd,decimal_mode * pm,decimal_record * pd,fp_exception_field_type * ps)264 __fast_double_to_decimal(double *dd, decimal_mode *pm, decimal_record *pd,
265 fp_exception_field_type *ps)
266 {
267 int i, is, esum, eround, hd;
268 double dds;
269 __ieee_flags_type fb;
270
271 if (pm->rd != fp_nearest)
272 return (1);
273
274 if (pm->df == fixed_form) {
275 /* F format */
276 if (pm->ndigits < 0 || pm->ndigits > __TBL_TENS_MAX)
277 return (1);
278 __get_ieee_flags(&fb);
279 dds = __dabs(dd);
280 esum = 0;
281 if (pm->ndigits) {
282 /* scale by a positive power of ten */
283 if (pm->ndigits > __TBL_TENS_EXACT) {
284 dds *= __tbl_tens[pm->ndigits];
285 esum = 2;
286 } else {
287 dds = __mul_set(dds, __tbl_tens[pm->ndigits],
288 &eround);
289 esum = eround;
290 }
291 }
292 if (dds > 2147483647999999744.0) {
293 __set_ieee_flags(&fb);
294 return (1);
295 }
296 dds = __arint_set_n(&dds, esum, &eround);
297 if (eround == 2) {
298 /* error is too large to round reliably; punt */
299 __set_ieee_flags(&fb);
300 return (1);
301 }
302 if (dds == 0.0) {
303 is = (pm->ndigits > 0)? pm->ndigits : 1;
304 for (i = 0; i < is; i++)
305 pd->ds[i] = '0';
306 pd->ds[is] = '\0';
307 eround++;
308 } else {
309 is = __double_to_digits(dds, pd->ds, pm->ndigits);
310 }
311 pd->ndigits = is;
312 pd->exponent = -pm->ndigits;
313 } else {
314 /* E format */
315 if (pm->ndigits < 1 || pm->ndigits > 18)
316 return (1);
317 __get_ieee_flags(&fb);
318 dds = __dabs(dd);
319 /* find the decade containing dds */
320 #ifdef _LITTLE_ENDIAN
321 hd = *(1+(int *)dd);
322 #else
323 hd = *(int *)dd;
324 #endif
325 hd = (hd >> 20) & 0x7ff;
326 if (hd >= 0x400) {
327 if (hd > 0x4e0)
328 i = TBL_DECADE_MAX;
329 else
330 i = TBL_DECADE_MAX - ((0x4e0 - hd) >> 2);
331 } else {
332 if (hd < 0x358)
333 i = 0;
334 else
335 i = TBL_DECADE_OFFSET - ((0x3ff - hd) >> 2);
336 }
337 while (dds < tbl_decade[i])
338 i--;
339 /* determine the power of ten by which to scale */
340 i = pm->ndigits - 1 - (i - TBL_DECADE_OFFSET);
341 esum = 0;
342 if (i > 0) {
343 /* scale by a positive power of ten */
344 if (i > __TBL_TENS_EXACT) {
345 if (i > __TBL_TENS_MAX) {
346 __set_ieee_flags(&fb);
347 return (1);
348 }
349 dds *= __tbl_tens[i];
350 esum = 2;
351 } else {
352 dds = __mul_set(dds, __tbl_tens[i], &eround);
353 esum = eround;
354 }
355 } else if (i < 0) {
356 /* scale by a negative power of ten */
357 if (-i > __TBL_TENS_EXACT) {
358 if (-i > __TBL_TENS_MAX) {
359 __set_ieee_flags(&fb);
360 return (1);
361 }
362 dds /= __tbl_tens[-i];
363 esum = 2;
364 } else {
365 dds = __div_set(dds, __tbl_tens[-i], &eround);
366 esum = eround;
367 }
368 }
369 dds = __arint_set_n(&dds, esum, &eround);
370 if (eround == 2) {
371 /* error is too large to round reliably; punt */
372 __set_ieee_flags(&fb);
373 return (1);
374 }
375 is = __double_to_digits(dds, pd->ds, 1);
376 if (is > pm->ndigits) {
377 /*
378 * The result rounded up to the next larger power
379 * of ten; just discard the last zero and adjust
380 * the exponent.
381 */
382 pd->ds[--is] = '\0';
383 i--;
384 }
385 pd->ndigits = is;
386 pd->exponent = -i;
387 }
388 *ps = (eround == 0)? 0 : (1 << fp_inexact);
389 __set_ieee_flags(&fb);
390 return (0);
391 }
392