]>
Commit | Line | Data |
---|---|---|
2aa62f2b | 1 | /* $NetBSD: gdtoa.c,v 1.1.1.1.4.1.4.1 2008/04/08 21:10:55 jdc Exp $ */\r |
2 | \r | |
3 | /****************************************************************\r | |
4 | \r | |
5 | The author of this software is David M. Gay.\r | |
6 | \r | |
7 | Copyright (C) 1998, 1999 by Lucent Technologies\r | |
8 | All Rights Reserved\r | |
9 | \r | |
10 | Permission to use, copy, modify, and distribute this software and\r | |
11 | its documentation for any purpose and without fee is hereby\r | |
12 | granted, provided that the above copyright notice appear in all\r | |
13 | copies and that both that the copyright notice and this\r | |
14 | permission notice and warranty disclaimer appear in supporting\r | |
15 | documentation, and that the name of Lucent or any of its entities\r | |
16 | not be used in advertising or publicity pertaining to\r | |
17 | distribution of the software without specific, written prior\r | |
18 | permission.\r | |
19 | \r | |
20 | LUCENT DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE,\r | |
21 | INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS.\r | |
22 | IN NO EVENT SHALL LUCENT OR ANY OF ITS ENTITIES BE LIABLE FOR ANY\r | |
23 | SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES\r | |
24 | WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER\r | |
25 | IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION,\r | |
26 | ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF\r | |
27 | THIS SOFTWARE.\r | |
28 | \r | |
29 | ****************************************************************/\r | |
30 | \r | |
31 | /* Please send bug reports to David M. Gay (dmg at acm dot org,\r | |
32 | * with " at " changed at "@" and " dot " changed to "."). */\r | |
33 | #include <LibConfig.h>\r | |
34 | \r | |
35 | #include "gdtoaimp.h"\r | |
36 | \r | |
37 | #if defined(_MSC_VER)\r | |
38 | /* Disable warnings about conversions to narrower data types. */\r | |
39 | #pragma warning ( disable : 4244 )\r | |
40 | // Squelch bogus warnings about uninitialized variable use.\r | |
41 | #pragma warning ( disable : 4701 )\r | |
42 | #endif\r | |
43 | \r | |
44 | static Bigint *\r | |
45 | bitstob(ULong *bits, int nbits, int *bbits)\r | |
46 | {\r | |
47 | int i, k;\r | |
48 | Bigint *b;\r | |
49 | ULong *be, *x, *x0;\r | |
50 | \r | |
51 | i = ULbits;\r | |
52 | k = 0;\r | |
53 | while(i < nbits) {\r | |
54 | i <<= 1;\r | |
55 | k++;\r | |
56 | }\r | |
57 | #ifndef Pack_32\r | |
58 | if (!k)\r | |
59 | k = 1;\r | |
60 | #endif\r | |
61 | b = Balloc(k);\r | |
62 | if (b == NULL)\r | |
63 | return NULL;\r | |
64 | be = bits + (((unsigned int)nbits - 1) >> kshift);\r | |
65 | x = x0 = b->x;\r | |
66 | do {\r | |
67 | *x++ = *bits & ALL_ON;\r | |
68 | #ifdef Pack_16\r | |
69 | *x++ = (*bits >> 16) & ALL_ON;\r | |
70 | #endif\r | |
71 | } while(++bits <= be);\r | |
72 | i = x - x0;\r | |
73 | while(!x0[--i])\r | |
74 | if (!i) {\r | |
75 | b->wds = 0;\r | |
76 | *bbits = 0;\r | |
77 | goto ret;\r | |
78 | }\r | |
79 | b->wds = i + 1;\r | |
80 | *bbits = i*ULbits + 32 - hi0bits(b->x[i]);\r | |
81 | ret:\r | |
82 | return b;\r | |
83 | }\r | |
84 | \r | |
85 | /* dtoa for IEEE arithmetic (dmg): convert double to ASCII string.\r | |
86 | *\r | |
87 | * Inspired by "How to Print Floating-Point Numbers Accurately" by\r | |
88 | * Guy L. Steele, Jr. and Jon L. White [Proc. ACM SIGPLAN '90, pp. 112-126].\r | |
89 | *\r | |
90 | * Modifications:\r | |
91 | * 1. Rather than iterating, we use a simple numeric overestimate\r | |
92 | * to determine k = floor(log10(d)). We scale relevant\r | |
93 | * quantities using O(log2(k)) rather than O(k) multiplications.\r | |
94 | * 2. For some modes > 2 (corresponding to ecvt and fcvt), we don't\r | |
95 | * try to generate digits strictly left to right. Instead, we\r | |
96 | * compute with fewer bits and propagate the carry if necessary\r | |
97 | * when rounding the final digit up. This is often faster.\r | |
98 | * 3. Under the assumption that input will be rounded nearest,\r | |
99 | * mode 0 renders 1e23 as 1e23 rather than 9.999999999999999e22.\r | |
100 | * That is, we allow equality in stopping tests when the\r | |
101 | * round-nearest rule will give the same floating-point value\r | |
102 | * as would satisfaction of the stopping test with strict\r | |
103 | * inequality.\r | |
104 | * 4. We remove common factors of powers of 2 from relevant\r | |
105 | * quantities.\r | |
106 | * 5. When converting floating-point integers less than 1e16,\r | |
107 | * we use floating-point arithmetic rather than resorting\r | |
108 | * to multiple-precision integers.\r | |
109 | * 6. When asked to produce fewer than 15 digits, we first try\r | |
110 | * to get by with floating-point arithmetic; we resort to\r | |
111 | * multiple-precision integer arithmetic only if we cannot\r | |
112 | * guarantee that the floating-point calculation has given\r | |
113 | * the correctly rounded result. For k requested digits and\r | |
114 | * "uniformly" distributed input, the probability is\r | |
115 | * something like 10^(k-15) that we must resort to the Long\r | |
116 | * calculation.\r | |
117 | */\r | |
118 | \r | |
119 | char *\r | |
120 | gdtoa\r | |
121 | (FPI *fpi, int be, ULong *bits, int *kindp, int mode, int ndigits, int *decpt, char **rve)\r | |
122 | {\r | |
123 | /* Arguments ndigits and decpt are similar to the second and third\r | |
124 | arguments of ecvt and fcvt; trailing zeros are suppressed from\r | |
125 | the returned string. If not null, *rve is set to point\r | |
126 | to the end of the return value. If d is +-Infinity or NaN,\r | |
127 | then *decpt is set to 9999.\r | |
128 | \r | |
129 | mode:\r | |
130 | 0 ==> shortest string that yields d when read in\r | |
131 | and rounded to nearest.\r | |
132 | 1 ==> like 0, but with Steele & White stopping rule;\r | |
133 | e.g. with IEEE P754 arithmetic , mode 0 gives\r | |
134 | 1e23 whereas mode 1 gives 9.999999999999999e22.\r | |
135 | 2 ==> max(1,ndigits) significant digits. This gives a\r | |
136 | return value similar to that of ecvt, except\r | |
137 | that trailing zeros are suppressed.\r | |
138 | 3 ==> through ndigits past the decimal point. This\r | |
139 | gives a return value similar to that from fcvt,\r | |
140 | except that trailing zeros are suppressed, and\r | |
141 | ndigits can be negative.\r | |
142 | 4-9 should give the same return values as 2-3, i.e.,\r | |
143 | 4 <= mode <= 9 ==> same return as mode\r | |
144 | 2 + (mode & 1). These modes are mainly for\r | |
145 | debugging; often they run slower but sometimes\r | |
146 | faster than modes 2-3.\r | |
147 | 4,5,8,9 ==> left-to-right digit generation.\r | |
148 | 6-9 ==> don't try fast floating-point estimate\r | |
149 | (if applicable).\r | |
150 | \r | |
151 | Values of mode other than 0-9 are treated as mode 0.\r | |
152 | \r | |
153 | Sufficient space is allocated to the return value\r | |
154 | to hold the suppressed trailing zeros.\r | |
155 | */\r | |
156 | \r | |
157 | int bbits, b2, b5, be0, dig, i, ieps, ilim = 0, ilim0, ilim1 = 0, inex;\r | |
158 | int j, jj1, k, k0, k_check, kind, leftright, m2, m5, nbits;\r | |
159 | int rdir, s2, s5, spec_case, try_quick;\r | |
160 | Long L;\r | |
161 | Bigint *b, *b1, *delta, *mlo, *mhi, *mhi1, *S;\r | |
162 | double d, d2, ds, eps;\r | |
163 | char *s, *s0;\r | |
164 | \r | |
165 | #ifndef MULTIPLE_THREADS\r | |
166 | if (dtoa_result) {\r | |
167 | freedtoa(dtoa_result);\r | |
168 | dtoa_result = 0;\r | |
169 | }\r | |
170 | #endif\r | |
171 | inex = 0;\r | |
172 | if (*kindp & STRTOG_NoMemory)\r | |
173 | return NULL;\r | |
174 | kind = *kindp &= ~STRTOG_Inexact;\r | |
175 | switch(kind & STRTOG_Retmask) {\r | |
176 | case STRTOG_Zero:\r | |
177 | goto ret_zero;\r | |
178 | case STRTOG_Normal:\r | |
179 | case STRTOG_Denormal:\r | |
180 | break;\r | |
181 | case STRTOG_Infinite:\r | |
182 | *decpt = -32768;\r | |
183 | return nrv_alloc("Infinity", rve, 8);\r | |
184 | case STRTOG_NaN:\r | |
185 | *decpt = -32768;\r | |
186 | return nrv_alloc("NaN", rve, 3);\r | |
187 | default:\r | |
188 | return 0;\r | |
189 | }\r | |
190 | b = bitstob(bits, nbits = fpi->nbits, &bbits);\r | |
191 | if (b == NULL)\r | |
192 | return NULL;\r | |
193 | be0 = be;\r | |
194 | if ( (i = trailz(b)) !=0) {\r | |
195 | rshift(b, i);\r | |
196 | be += i;\r | |
197 | bbits -= i;\r | |
198 | }\r | |
199 | if (!b->wds) {\r | |
200 | Bfree(b);\r | |
201 | ret_zero:\r | |
202 | *decpt = 1;\r | |
203 | return nrv_alloc("0", rve, 1);\r | |
204 | }\r | |
205 | \r | |
206 | dval(d) = b2d(b, &i);\r | |
207 | i = be + bbits - 1;\r | |
208 | word0(d) &= Frac_mask1;\r | |
209 | word0(d) |= Exp_11;\r | |
210 | #ifdef IBM\r | |
211 | if ( (j = 11 - hi0bits(word0(d) & Frac_mask)) !=0)\r | |
212 | dval(d) /= 1 << j;\r | |
213 | #endif\r | |
214 | \r | |
215 | /* log(x) ~=~ log(1.5) + (x-1.5)/1.5\r | |
216 | * log10(x) = log(x) / log(10)\r | |
217 | * ~=~ log(1.5)/log(10) + (x-1.5)/(1.5*log(10))\r | |
218 | * log10(d) = (i-Bias)*log(2)/log(10) + log10(d2)\r | |
219 | *\r | |
220 | * This suggests computing an approximation k to log10(d) by\r | |
221 | *\r | |
222 | * k = (i - Bias)*0.301029995663981\r | |
223 | * + ( (d2-1.5)*0.289529654602168 + 0.176091259055681 );\r | |
224 | *\r | |
225 | * We want k to be too large rather than too small.\r | |
226 | * The error in the first-order Taylor series approximation\r | |
227 | * is in our favor, so we just round up the constant enough\r | |
228 | * to compensate for any error in the multiplication of\r | |
229 | * (i - Bias) by 0.301029995663981; since |i - Bias| <= 1077,\r | |
230 | * and 1077 * 0.30103 * 2^-52 ~=~ 7.2e-14,\r | |
231 | * adding 1e-13 to the constant term more than suffices.\r | |
232 | * Hence we adjust the constant term to 0.1760912590558.\r | |
233 | * (We could get a more accurate k by invoking log10,\r | |
234 | * but this is probably not worthwhile.)\r | |
235 | */\r | |
236 | #ifdef IBM\r | |
237 | i <<= 2;\r | |
238 | i += j;\r | |
239 | #endif\r | |
240 | ds = (dval(d)-1.5)*0.289529654602168 + 0.1760912590558 + i*0.301029995663981;\r | |
241 | \r | |
242 | /* correct assumption about exponent range */\r | |
243 | if ((j = i) < 0)\r | |
244 | j = -j;\r | |
245 | if ((j -= 1077) > 0)\r | |
246 | ds += j * 7e-17;\r | |
247 | \r | |
248 | k = (int)ds;\r | |
249 | if (ds < 0. && ds != k)\r | |
250 | k--; /* want k = floor(ds) */\r | |
251 | k_check = 1;\r | |
252 | #ifdef IBM\r | |
253 | j = be + bbits - 1;\r | |
254 | if ( (jj1 = j & 3) !=0)\r | |
255 | dval(d) *= 1 << jj1;\r | |
256 | word0(d) += j << Exp_shift - 2 & Exp_mask;\r | |
257 | #else\r | |
258 | word0(d) += (be + bbits - 1) << Exp_shift;\r | |
259 | #endif\r | |
260 | if (k >= 0 && k <= Ten_pmax) {\r | |
261 | if (dval(d) < tens[k])\r | |
262 | k--;\r | |
263 | k_check = 0;\r | |
264 | }\r | |
265 | j = bbits - i - 1;\r | |
266 | if (j >= 0) {\r | |
267 | b2 = 0;\r | |
268 | s2 = j;\r | |
269 | }\r | |
270 | else {\r | |
271 | b2 = -j;\r | |
272 | s2 = 0;\r | |
273 | }\r | |
274 | if (k >= 0) {\r | |
275 | b5 = 0;\r | |
276 | s5 = k;\r | |
277 | s2 += k;\r | |
278 | }\r | |
279 | else {\r | |
280 | b2 -= k;\r | |
281 | b5 = -k;\r | |
282 | s5 = 0;\r | |
283 | }\r | |
284 | if (mode < 0 || mode > 9)\r | |
285 | mode = 0;\r | |
286 | try_quick = 1;\r | |
287 | if (mode > 5) {\r | |
288 | mode -= 4;\r | |
289 | try_quick = 0;\r | |
290 | }\r | |
291 | leftright = 1;\r | |
292 | switch(mode) {\r | |
293 | case 0:\r | |
294 | case 1:\r | |
295 | ilim = ilim1 = -1;\r | |
296 | i = (int)(nbits * .30103) + 3;\r | |
297 | ndigits = 0;\r | |
298 | break;\r | |
299 | case 2:\r | |
300 | leftright = 0;\r | |
301 | /*FALLTHROUGH*/\r | |
302 | case 4:\r | |
303 | if (ndigits <= 0)\r | |
304 | ndigits = 1;\r | |
305 | ilim = ilim1 = i = ndigits;\r | |
306 | break;\r | |
307 | case 3:\r | |
308 | leftright = 0;\r | |
309 | /*FALLTHROUGH*/\r | |
310 | case 5:\r | |
311 | i = ndigits + k + 1;\r | |
312 | ilim = i;\r | |
313 | ilim1 = i - 1;\r | |
314 | if (i <= 0)\r | |
315 | i = 1;\r | |
316 | }\r | |
317 | s = s0 = rv_alloc((size_t)i);\r | |
318 | if (s == NULL)\r | |
319 | return NULL;\r | |
320 | \r | |
321 | if ( (rdir = fpi->rounding - 1) !=0) {\r | |
322 | if (rdir < 0)\r | |
323 | rdir = 2;\r | |
324 | if (kind & STRTOG_Neg)\r | |
325 | rdir = 3 - rdir;\r | |
326 | }\r | |
327 | \r | |
328 | /* Now rdir = 0 ==> round near, 1 ==> round up, 2 ==> round down. */\r | |
329 | \r | |
330 | if (ilim >= 0 && ilim <= Quick_max && try_quick && !rdir\r | |
331 | #ifndef IMPRECISE_INEXACT\r | |
332 | && k == 0\r | |
333 | #endif\r | |
334 | ) {\r | |
335 | \r | |
336 | /* Try to get by with floating-point arithmetic. */\r | |
337 | \r | |
338 | i = 0;\r | |
339 | d2 = dval(d);\r | |
340 | #ifdef IBM\r | |
341 | if ( (j = 11 - hi0bits(word0(d) & Frac_mask)) !=0)\r | |
342 | dval(d) /= 1 << j;\r | |
343 | #endif\r | |
344 | k0 = k;\r | |
345 | ilim0 = ilim;\r | |
346 | ieps = 2; /* conservative */\r | |
347 | if (k > 0) {\r | |
348 | ds = tens[k&0xf];\r | |
349 | j = (unsigned int)k >> 4;\r | |
350 | if (j & Bletch) {\r | |
351 | /* prevent overflows */\r | |
352 | j &= Bletch - 1;\r | |
353 | dval(d) /= bigtens[n_bigtens-1];\r | |
354 | ieps++;\r | |
355 | }\r | |
356 | for(; j; j /= 2, i++)\r | |
357 | if (j & 1) {\r | |
358 | ieps++;\r | |
359 | ds *= bigtens[i];\r | |
360 | }\r | |
361 | }\r | |
362 | else {\r | |
363 | ds = 1.;\r | |
364 | if ( (jj1 = -k) !=0) {\r | |
365 | dval(d) *= tens[jj1 & 0xf];\r | |
366 | for(j = jj1 >> 4; j; j >>= 1, i++)\r | |
367 | if (j & 1) {\r | |
368 | ieps++;\r | |
369 | dval(d) *= bigtens[i];\r | |
370 | }\r | |
371 | }\r | |
372 | }\r | |
373 | if (k_check && dval(d) < 1. && ilim > 0) {\r | |
374 | if (ilim1 <= 0)\r | |
375 | goto fast_failed;\r | |
376 | ilim = ilim1;\r | |
377 | k--;\r | |
378 | dval(d) *= 10.;\r | |
379 | ieps++;\r | |
380 | }\r | |
381 | dval(eps) = ieps*dval(d) + 7.;\r | |
382 | word0(eps) -= (P-1)*Exp_msk1;\r | |
383 | if (ilim == 0) {\r | |
384 | S = mhi = 0;\r | |
385 | dval(d) -= 5.;\r | |
386 | if (dval(d) > dval(eps))\r | |
387 | goto one_digit;\r | |
388 | if (dval(d) < -dval(eps))\r | |
389 | goto no_digits;\r | |
390 | goto fast_failed;\r | |
391 | }\r | |
392 | #ifndef No_leftright\r | |
393 | if (leftright) {\r | |
394 | /* Use Steele & White method of only\r | |
395 | * generating digits needed.\r | |
396 | */\r | |
397 | dval(eps) = ds*0.5/tens[ilim-1] - dval(eps);\r | |
398 | for(i = 0;;) {\r | |
399 | L = (Long)(dval(d)/ds);\r | |
400 | dval(d) -= L*ds;\r | |
401 | *s++ = '0' + (int)L;\r | |
402 | if (dval(d) < dval(eps)) {\r | |
403 | if (dval(d))\r | |
404 | inex = STRTOG_Inexlo;\r | |
405 | goto ret1;\r | |
406 | }\r | |
407 | if (ds - dval(d) < dval(eps))\r | |
408 | goto bump_up;\r | |
409 | if (++i >= ilim)\r | |
410 | break;\r | |
411 | dval(eps) *= 10.;\r | |
412 | dval(d) *= 10.;\r | |
413 | }\r | |
414 | }\r | |
415 | else {\r | |
416 | #endif\r | |
417 | /* Generate ilim digits, then fix them up. */\r | |
418 | dval(eps) *= tens[ilim-1];\r | |
419 | for(i = 1;; i++, dval(d) *= 10.) {\r | |
420 | if ( (L = (Long)(dval(d)/ds)) !=0)\r | |
421 | dval(d) -= L*ds;\r | |
422 | *s++ = '0' + (int)L;\r | |
423 | if (i == ilim) {\r | |
424 | ds *= 0.5;\r | |
425 | if (dval(d) > ds + dval(eps))\r | |
426 | goto bump_up;\r | |
427 | else if (dval(d) < ds - dval(eps)) {\r | |
428 | while(*--s == '0'){}\r | |
429 | s++;\r | |
430 | if (dval(d))\r | |
431 | inex = STRTOG_Inexlo;\r | |
432 | goto ret1;\r | |
433 | }\r | |
434 | break;\r | |
435 | }\r | |
436 | }\r | |
437 | #ifndef No_leftright\r | |
438 | }\r | |
439 | #endif\r | |
440 | fast_failed:\r | |
441 | s = s0;\r | |
442 | dval(d) = d2;\r | |
443 | k = k0;\r | |
444 | ilim = ilim0;\r | |
445 | }\r | |
446 | \r | |
447 | /* Do we have a "small" integer? */\r | |
448 | \r | |
449 | if (be >= 0 && k <= Int_max) {\r | |
450 | /* Yes. */\r | |
451 | ds = tens[k];\r | |
452 | if (ndigits < 0 && ilim <= 0) {\r | |
453 | S = mhi = 0;\r | |
454 | if (ilim < 0 || dval(d) <= 5*ds)\r | |
455 | goto no_digits;\r | |
456 | goto one_digit;\r | |
457 | }\r | |
458 | for(i = 1;; i++, dval(d) *= 10.) {\r | |
459 | L = dval(d) / ds;\r | |
460 | dval(d) -= L*ds;\r | |
461 | #ifdef Check_FLT_ROUNDS\r | |
462 | /* If FLT_ROUNDS == 2, L will usually be high by 1 */\r | |
463 | if (dval(d) < 0) {\r | |
464 | L--;\r | |
465 | dval(d) += ds;\r | |
466 | }\r | |
467 | #endif\r | |
468 | *s++ = '0' + (int)L;\r | |
469 | if (dval(d) == 0.)\r | |
470 | break;\r | |
471 | if (i == ilim) {\r | |
472 | if (rdir) {\r | |
473 | if (rdir == 1)\r | |
474 | goto bump_up;\r | |
475 | inex = STRTOG_Inexlo;\r | |
476 | goto ret1;\r | |
477 | }\r | |
478 | dval(d) += dval(d);\r | |
479 | if (dval(d) > ds || (dval(d) == ds && L & 1)) {\r | |
480 | bump_up:\r | |
481 | inex = STRTOG_Inexhi;\r | |
482 | while(*--s == '9')\r | |
483 | if (s == s0) {\r | |
484 | k++;\r | |
485 | *s = '0';\r | |
486 | break;\r | |
487 | }\r | |
488 | ++*s++;\r | |
489 | }\r | |
490 | else\r | |
491 | inex = STRTOG_Inexlo;\r | |
492 | break;\r | |
493 | }\r | |
494 | }\r | |
495 | goto ret1;\r | |
496 | }\r | |
497 | \r | |
498 | m2 = b2;\r | |
499 | m5 = b5;\r | |
500 | mhi = mlo = 0;\r | |
501 | if (leftright) {\r | |
502 | if (mode < 2) {\r | |
503 | i = nbits - bbits;\r | |
504 | if (be - i++ < fpi->emin)\r | |
505 | /* denormal */\r | |
506 | i = be - fpi->emin + 1;\r | |
507 | }\r | |
508 | else {\r | |
509 | j = ilim - 1;\r | |
510 | if (m5 >= j)\r | |
511 | m5 -= j;\r | |
512 | else {\r | |
513 | s5 += j -= m5;\r | |
514 | b5 += j;\r | |
515 | m5 = 0;\r | |
516 | }\r | |
517 | if ((i = ilim) < 0) {\r | |
518 | m2 -= i;\r | |
519 | i = 0;\r | |
520 | }\r | |
521 | }\r | |
522 | b2 += i;\r | |
523 | s2 += i;\r | |
524 | mhi = i2b(1);\r | |
525 | }\r | |
526 | if (m2 > 0 && s2 > 0) {\r | |
527 | i = m2 < s2 ? m2 : s2;\r | |
528 | b2 -= i;\r | |
529 | m2 -= i;\r | |
530 | s2 -= i;\r | |
531 | }\r | |
532 | if (b5 > 0) {\r | |
533 | if (leftright) {\r | |
534 | if (m5 > 0) {\r | |
535 | mhi = pow5mult(mhi, m5);\r | |
536 | if (mhi == NULL)\r | |
537 | return NULL;\r | |
538 | b1 = mult(mhi, b);\r | |
539 | if (b1 == NULL)\r | |
540 | return NULL;\r | |
541 | Bfree(b);\r | |
542 | b = b1;\r | |
543 | }\r | |
544 | if ( (j = b5 - m5) !=0) {\r | |
545 | b = pow5mult(b, j);\r | |
546 | if (b == NULL)\r | |
547 | return NULL;\r | |
548 | }\r | |
549 | }\r | |
550 | else {\r | |
551 | b = pow5mult(b, b5);\r | |
552 | if (b == NULL)\r | |
553 | return NULL;\r | |
554 | }\r | |
555 | }\r | |
556 | S = i2b(1);\r | |
557 | if (S == NULL)\r | |
558 | return NULL;\r | |
559 | if (s5 > 0) {\r | |
560 | S = pow5mult(S, s5);\r | |
561 | if (S == NULL)\r | |
562 | return NULL;\r | |
563 | }\r | |
564 | \r | |
565 | /* Check for special case that d is a normalized power of 2. */\r | |
566 | \r | |
567 | spec_case = 0;\r | |
568 | if (mode < 2) {\r | |
569 | if (bbits == 1 && be0 > fpi->emin + 1) {\r | |
570 | /* The special case */\r | |
571 | b2++;\r | |
572 | s2++;\r | |
573 | spec_case = 1;\r | |
574 | }\r | |
575 | }\r | |
576 | \r | |
577 | /* Arrange for convenient computation of quotients:\r | |
578 | * shift left if necessary so divisor has 4 leading 0 bits.\r | |
579 | *\r | |
580 | * Perhaps we should just compute leading 28 bits of S once\r | |
581 | * and for all and pass them and a shift to quorem, so it\r | |
582 | * can do shifts and ors to compute the numerator for q.\r | |
583 | */\r | |
584 | #ifdef Pack_32\r | |
585 | if ( (i = ((s5 ? 32 - hi0bits(S->x[S->wds-1]) : 1) + s2) & 0x1f) !=0)\r | |
586 | i = 32 - i;\r | |
587 | #else\r | |
588 | if ( (i = ((s5 ? 32 - hi0bits(S->x[S->wds-1]) : 1) + s2) & 0xf) !=0)\r | |
589 | i = 16 - i;\r | |
590 | #endif\r | |
591 | if (i > 4) {\r | |
592 | i -= 4;\r | |
593 | b2 += i;\r | |
594 | m2 += i;\r | |
595 | s2 += i;\r | |
596 | }\r | |
597 | else if (i < 4) {\r | |
598 | i += 28;\r | |
599 | b2 += i;\r | |
600 | m2 += i;\r | |
601 | s2 += i;\r | |
602 | }\r | |
603 | if (b2 > 0)\r | |
604 | b = lshift(b, b2);\r | |
605 | if (s2 > 0)\r | |
606 | S = lshift(S, s2);\r | |
607 | if (k_check) {\r | |
608 | if (cmp(b,S) < 0) {\r | |
609 | k--;\r | |
610 | b = multadd(b, 10, 0); /* we botched the k estimate */\r | |
611 | if (b == NULL)\r | |
612 | return NULL;\r | |
613 | if (leftright) {\r | |
614 | mhi = multadd(mhi, 10, 0);\r | |
615 | if (mhi == NULL)\r | |
616 | return NULL;\r | |
617 | }\r | |
618 | ilim = ilim1;\r | |
619 | }\r | |
620 | }\r | |
621 | if (ilim <= 0 && mode > 2) {\r | |
622 | if (ilim < 0 || cmp(b,S = multadd(S,5,0)) <= 0) {\r | |
623 | /* no digits, fcvt style */\r | |
624 | no_digits:\r | |
625 | k = -1 - ndigits;\r | |
626 | inex = STRTOG_Inexlo;\r | |
627 | goto ret;\r | |
628 | }\r | |
629 | one_digit:\r | |
630 | inex = STRTOG_Inexhi;\r | |
631 | *s++ = '1';\r | |
632 | k++;\r | |
633 | goto ret;\r | |
634 | }\r | |
635 | if (leftright) {\r | |
636 | if (m2 > 0) {\r | |
637 | mhi = lshift(mhi, m2);\r | |
638 | if (mhi == NULL)\r | |
639 | return NULL;\r | |
640 | }\r | |
641 | \r | |
642 | /* Compute mlo -- check for special case\r | |
643 | * that d is a normalized power of 2.\r | |
644 | */\r | |
645 | \r | |
646 | mlo = mhi;\r | |
647 | if (spec_case) {\r | |
648 | mhi = Balloc(mhi->k);\r | |
649 | if (mhi == NULL)\r | |
650 | return NULL;\r | |
651 | Bcopy(mhi, mlo);\r | |
652 | mhi = lshift(mhi, 1);\r | |
653 | if (mhi == NULL)\r | |
654 | return NULL;\r | |
655 | }\r | |
656 | \r | |
657 | for(i = 1;;i++) {\r | |
658 | dig = quorem(b,S) + '0';\r | |
659 | /* Do we yet have the shortest decimal string\r | |
660 | * that will round to d?\r | |
661 | */\r | |
662 | j = cmp(b, mlo);\r | |
663 | delta = diff(S, mhi);\r | |
664 | if (delta == NULL)\r | |
665 | return NULL;\r | |
666 | jj1 = delta->sign ? 1 : cmp(b, delta);\r | |
667 | Bfree(delta);\r | |
668 | #ifndef ROUND_BIASED\r | |
669 | if (jj1 == 0 && !mode && !(bits[0] & 1) && !rdir) {\r | |
670 | if (dig == '9')\r | |
671 | goto round_9_up;\r | |
672 | if (j <= 0) {\r | |
673 | if (b->wds > 1 || b->x[0])\r | |
674 | inex = STRTOG_Inexlo;\r | |
675 | }\r | |
676 | else {\r | |
677 | dig++;\r | |
678 | inex = STRTOG_Inexhi;\r | |
679 | }\r | |
680 | *s++ = dig;\r | |
681 | goto ret;\r | |
682 | }\r | |
683 | #endif\r | |
684 | if (j < 0 || (j == 0 && !mode\r | |
685 | #ifndef ROUND_BIASED\r | |
686 | && !(bits[0] & 1)\r | |
687 | #endif\r | |
688 | )) {\r | |
689 | if (rdir && (b->wds > 1 || b->x[0])) {\r | |
690 | if (rdir == 2) {\r | |
691 | inex = STRTOG_Inexlo;\r | |
692 | goto accept;\r | |
693 | }\r | |
694 | while (cmp(S,mhi) > 0) {\r | |
695 | *s++ = dig;\r | |
696 | mhi1 = multadd(mhi, 10, 0);\r | |
697 | if (mhi1 == NULL)\r | |
698 | return NULL;\r | |
699 | if (mlo == mhi)\r | |
700 | mlo = mhi1;\r | |
701 | mhi = mhi1;\r | |
702 | b = multadd(b, 10, 0);\r | |
703 | if (b == NULL)\r | |
704 | return NULL;\r | |
705 | dig = quorem(b,S) + '0';\r | |
706 | }\r | |
707 | if (dig++ == '9')\r | |
708 | goto round_9_up;\r | |
709 | inex = STRTOG_Inexhi;\r | |
710 | goto accept;\r | |
711 | }\r | |
712 | if (jj1 > 0) {\r | |
713 | b = lshift(b, 1);\r | |
714 | if (b == NULL)\r | |
715 | return NULL;\r | |
716 | jj1 = cmp(b, S);\r | |
717 | if ((jj1 > 0 || (jj1 == 0 && dig & 1))\r | |
718 | && dig++ == '9')\r | |
719 | goto round_9_up;\r | |
720 | inex = STRTOG_Inexhi;\r | |
721 | }\r | |
722 | if (b->wds > 1 || b->x[0])\r | |
723 | inex = STRTOG_Inexlo;\r | |
724 | accept:\r | |
725 | *s++ = dig;\r | |
726 | goto ret;\r | |
727 | }\r | |
728 | if (jj1 > 0 && rdir != 2) {\r | |
729 | if (dig == '9') { /* possible if i == 1 */\r | |
730 | round_9_up:\r | |
731 | *s++ = '9';\r | |
732 | inex = STRTOG_Inexhi;\r | |
733 | goto roundoff;\r | |
734 | }\r | |
735 | inex = STRTOG_Inexhi;\r | |
736 | *s++ = dig + 1;\r | |
737 | goto ret;\r | |
738 | }\r | |
739 | *s++ = dig;\r | |
740 | if (i == ilim)\r | |
741 | break;\r | |
742 | b = multadd(b, 10, 0);\r | |
743 | if (b == NULL)\r | |
744 | return NULL;\r | |
745 | if (mlo == mhi) {\r | |
746 | mlo = mhi = multadd(mhi, 10, 0);\r | |
747 | if (mlo == NULL)\r | |
748 | return NULL;\r | |
749 | }\r | |
750 | else {\r | |
751 | mlo = multadd(mlo, 10, 0);\r | |
752 | if (mlo == NULL)\r | |
753 | return NULL;\r | |
754 | mhi = multadd(mhi, 10, 0);\r | |
755 | if (mhi == NULL)\r | |
756 | return NULL;\r | |
757 | }\r | |
758 | }\r | |
759 | }\r | |
760 | else\r | |
761 | for(i = 1;; i++) {\r | |
762 | *s++ = dig = quorem(b,S) + '0';\r | |
763 | if (i >= ilim)\r | |
764 | break;\r | |
765 | b = multadd(b, 10, 0);\r | |
766 | if (b == NULL)\r | |
767 | return NULL;\r | |
768 | }\r | |
769 | \r | |
770 | /* Round off last digit */\r | |
771 | \r | |
772 | if (rdir) {\r | |
773 | if (rdir == 2 || (b->wds <= 1 && !b->x[0]))\r | |
774 | goto chopzeros;\r | |
775 | goto roundoff;\r | |
776 | }\r | |
777 | b = lshift(b, 1);\r | |
778 | if (b == NULL)\r | |
779 | return NULL;\r | |
780 | j = cmp(b, S);\r | |
781 | if (j > 0 || (j == 0 && dig & 1)) {\r | |
782 | roundoff:\r | |
783 | inex = STRTOG_Inexhi;\r | |
784 | while(*--s == '9')\r | |
785 | if (s == s0) {\r | |
786 | k++;\r | |
787 | *s++ = '1';\r | |
788 | goto ret;\r | |
789 | }\r | |
790 | ++*s++;\r | |
791 | }\r | |
792 | else {\r | |
793 | chopzeros:\r | |
794 | if (b->wds > 1 || b->x[0])\r | |
795 | inex = STRTOG_Inexlo;\r | |
796 | while(*--s == '0'){}\r | |
797 | s++;\r | |
798 | }\r | |
799 | ret:\r | |
800 | Bfree(S);\r | |
801 | if (mhi) {\r | |
802 | if (mlo && mlo != mhi)\r | |
803 | Bfree(mlo);\r | |
804 | Bfree(mhi);\r | |
805 | }\r | |
806 | ret1:\r | |
807 | Bfree(b);\r | |
808 | *s = 0;\r | |
809 | *decpt = k + 1;\r | |
810 | if (rve)\r | |
811 | *rve = s;\r | |
812 | *kindp |= inex;\r | |
813 | return s0;\r | |
814 | }\r |