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