Embedded Artistry libc
C Standard Library Support for Bare-metal Systems
strtodg.c
Go to the documentation of this file.
1 /****************************************************************
2 
3 The author of this software is David M. Gay.
4 
5 Copyright (C) 1998-2001 by Lucent Technologies
6 All Rights Reserved
7 
8 Permission to use, copy, modify, and distribute this software and
9 its documentation for any purpose and without fee is hereby
10 granted, provided that the above copyright notice appear in all
11 copies and that both that the copyright notice and this
12 permission notice and warranty disclaimer appear in supporting
13 documentation, and that the name of Lucent or any of its entities
14 not be used in advertising or publicity pertaining to
15 distribution of the software without specific, written prior
16 permission.
17 
18 LUCENT DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE,
19 INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS.
20 IN NO EVENT SHALL LUCENT OR ANY OF ITS ENTITIES BE LIABLE FOR ANY
21 SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
22 WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER
23 IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION,
24 ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF
25 THIS SOFTWARE.
26 
27 ****************************************************************/
28 
29 /* Please send bug reports to David M. Gay (dmg at acm dot org,
30  * with " at " changed at "@" and " dot " changed to "."). */
31 
32 #include "gdtoaimp.h"
33 #include <stdint.h>
34 
35 #ifdef USE_LOCALE
36 #include "locale.h"
37 #endif
38 
39 static CONST int fivesbits[] = {0,
40  3,
41  5,
42  7,
43  10,
44  12,
45  14,
46  17,
47  19,
48  21,
49  24,
50  26,
51  28,
52  31,
53  33,
54  35,
55  38,
56  40,
57  42,
58  45,
59  47,
60  49,
61  52
62 #ifdef VAX
63  ,
64  54,
65  56
66 #endif
67 };
68 
69 Bigint*
70 #ifdef KR_headers
71  increment(b) Bigint* b;
72 #else
74 #endif
75 {
76  ULong *x, *xe;
77  Bigint* b1;
78 #ifdef Pack_16
79  ULong carry = 1, y;
80 #endif
81 
82  x = b->x;
83  xe = x + b->wds;
84 #ifdef Pack_32
85  do
86  {
87  if(*x < (ULong)0xffffffffL)
88  {
89  ++*x;
90  return b;
91  }
92  *x++ = 0;
93  } while(x < xe);
94 #else
95  do
96  {
97  y = *x + carry;
98  carry = y >> 16;
99  *x++ = y & 0xffff;
100  if(!carry)
101  return b;
102  } while(x < xe);
103  if(carry)
104 #endif
105  {
106  if(b->wds >= b->maxwds)
107  {
108  b1 = Balloc(b->k + 1);
109  Bcopy(b1, b);
110  Bfree(b);
111  b = b1;
112  }
113  b->x[b->wds++] = 1;
114  }
115  return b;
116 }
117 
118 int
119 #ifdef KR_headers
120  decrement(b) Bigint* b;
121 #else
123 #endif
124 {
125  ULong *x, *xe;
126 #ifdef Pack_16
127  ULong borrow = 1, y;
128 #endif
129 
130  x = b->x;
131  xe = x + b->wds;
132 #ifdef Pack_32
133  do
134  {
135  if(*x)
136  {
137  --*x;
138  break;
139  }
140  *x++ = 0xffffffffL;
141  } while(x < xe);
142 #else
143  do
144  {
145  y = *x - borrow;
146  borrow = (y & 0x10000) >> 16;
147  *x++ = y & 0xffff;
148  } while(borrow && x < xe);
149 #endif
150  return STRTOG_Inexlo;
151 }
152 
153 static int
154 #ifdef KR_headers
155  all_on(b, n) Bigint* b;
156 int n;
157 #else
158  all_on(Bigint* b, int n)
159 #endif
160 {
161  ULong *x, *xe;
162 
163  x = b->x;
164  xe = x + (n >> kshift);
165  while(x < xe)
166  {
167  {
168  if((*x++ & ALL_ON) != ALL_ON)
169  {
170  {
171  return 0;
172  }
173  }
174  }
175  }
176  if(n &= kmask)
177  {
178  {
179  return ((*x | (ALL_ON << n)) & ALL_ON) == ALL_ON;
180  }
181  }
182  return 1;
183 }
184 
185 Bigint*
186 #ifdef KR_headers
187  set_ones(b, n) Bigint* b;
188 int n;
189 #else
190  set_ones(Bigint* b, int n)
191 #endif
192 {
193  int k;
194  ULong *x, *xe;
195 
196  k = (n + ((1 << kshift) - 1)) >> kshift;
197  if(b->k < k)
198  {
199  Bfree(b);
200  b = Balloc(k);
201  }
202  k = n >> kshift;
203  if(n &= kmask)
204  {
205  {
206  k++;
207  }
208  }
209  b->wds = k;
210  x = b->x;
211  xe = x + k;
212  while(x < xe)
213  {
214  {
215  *x++ = ALL_ON;
216  }
217  }
218  if(n)
219  {
220  {
221  x[-1] >>= ULbits - n;
222  }
223  }
224  return b;
225 }
226 
227 static int rvOK
228 #ifdef KR_headers
229  (d, fpi, exp, bits, exact, rd, irv) double d;
230 FPI* fpi;
231 Long* exp;
232 ULong* bits;
233 int exact, rd, *irv;
234 #else
235  (double d, FPI* fpi, Long* exp, ULong* bits, int exact, int rd, int* irv)
236 #endif
237 {
238  Bigint* b;
239  ULong carry, inex, lostbits;
240  int bdif, e, j, k, k1, nb, rv;
241 
242  carry = rv = 0;
243  b = d2b(d, &e, &bdif);
244  bdif -= nb = fpi->nbits;
245  e += bdif;
246  if(bdif <= 0)
247  {
248  if(exact)
249  {
250  {
251  goto trunc;
252  }
253  }
254  goto ret;
255  }
256  if(P == nb)
257  {
258  if(
259 #ifndef IMPRECISE_INEXACT
260  exact &&
261 #endif
262  fpi->rounding ==
263 #ifdef RND_PRODQUOT
265 #else
266  Flt_Rounds
267 #endif
268  )
269  goto trunc;
270  goto ret;
271  }
272  switch(rd)
273  {
274  case 1:
275  goto trunc;
276  case 2:
277  break;
278  default: /* round near */
279  k = bdif - 1;
280  if(k < 0)
281  {
282  {
283  goto trunc;
284  }
285  }
286  if(!k)
287  {
288  if(!exact)
289  {
290  {
291  goto ret;
292  }
293  }
294  if(b->x[0] & 2)
295  {
296  {
297  break;
298  }
299  }
300  goto trunc;
301  }
302  if(b->x[k >> kshift] & ((ULong)1 << (k & kmask)))
303  {
304  {
305  break;
306  }
307  }
308  goto trunc;
309  }
310  /* "break" cases: round up 1 bit, then truncate; bdif > 0 */
311  carry = 1;
312 trunc:
313  inex = lostbits = 0;
314  if(bdif > 0)
315  {
316  if((lostbits = any_on(b, bdif)) != 0)
317  {
318  {
319  inex = STRTOG_Inexlo;
320  }
321  }
322  rshift(b, bdif);
323  if(carry)
324  {
325  inex = STRTOG_Inexhi;
326  b = increment(b);
327  if((j = nb & kmask) != 0)
328  {
329  {
330  j = ULbits - j;
331  }
332  }
333  if(hi0bits(b->x[b->wds - 1]) != j)
334  {
335  if(!lostbits)
336  {
337  {
338  lostbits = b->x[0] & 1;
339  }
340  }
341  rshift(b, 1);
342  e++;
343  }
344  }
345  }
346  else if(bdif < 0)
347  {
348  {
349  b = lshift(b, -bdif);
350  }
351  }
352  if(e < fpi->emin)
353  {
354  k = fpi->emin - e;
355  e = fpi->emin;
356  if(k > nb || fpi->sudden_underflow)
357  {
358  b->wds = inex = 0;
360  }
361  else
362  {
363  k1 = k - 1;
364  if(k1 > 0 && !lostbits)
365  {
366  {
367  lostbits = any_on(b, k1);
368  }
369  }
370  if(!lostbits && !exact)
371  {
372  {
373  goto ret;
374  }
375  }
376  lostbits |= carry = b->x[k1 >> kshift] & (1 << (k1 & kmask));
377  rshift(b, k);
378  *irv = STRTOG_Denormal;
379  if(carry)
380  {
381  b = increment(b);
383  }
384  else if(lostbits)
385  {
386  {
388  }
389  }
390  }
391  }
392  else if(e > fpi->emax)
393  {
394  e = fpi->emax + 1;
396 #ifndef NO_ERRNO
397  errno = ERANGE;
398 #endif
399  b->wds = inex = 0;
400  }
401  *exp = e;
402  copybits(bits, nb, b);
403  *irv |= inex;
404  rv = 1;
405 ret:
406  Bfree(b);
407  return rv;
408 }
409 
410 static int
411 #ifdef KR_headers
412  mantbits(d) double d;
413 #else
414  mantbits(double d)
415 #endif
416 {
417  ULong L;
418 #ifdef VAX
419  L = word1(d) << 16 | word1(d) >> 16;
420  if(L)
421 #else
422  if((L = word1(d)) != 0)
423 #endif
424  {
425  return P - lo0bits(&L);
426  }
427 #ifdef VAX
428  L = word0(d) << 16 | word0(d) >> 16 | Exp_msk11;
429 #else
430  L = word0(d) | Exp_msk1;
431 #endif
432  return P - 32 - lo0bits(&L);
433 }
434 
435 int strtodg
436 #ifdef KR_headers
437  (s00, se, fpi, exp, bits) CONST char* s00;
438 char** se;
439 FPI* fpi;
440 Long* exp;
441 ULong* bits;
442 #else
443  (CONST char* s00, char** se, FPI* fpi, Long* exp, ULong* bits)
444 #endif
445 {
446  int abe, abits, asub;
447  int bb0, bb2, bb5, bbe, bd2, bd5, bbbits, bs2, c, decpt, denorm;
448  int dsign, e, e1, e2, emin, esign, finished, i, inex, irv;
449  int j, k, nbits, nd, nd0, nf, nz, nz0, rd, rvbits, rve, rve1, sign;
450  int sudden_underflow = 0;
451  CONST char *s, *s0, *s1;
452  double adj, adj0, rv, tol;
453  Long L;
454  ULong y, z;
455  Bigint *ab, *bb, *bb1, *bd, *bd0, *bs, *delta, *rvb, *rvb0;
456 
457  irv = STRTOG_Zero;
458  denorm = sign = nz0 = nz = 0;
459  dval(rv) = 0.;
460  rvb = 0;
461  nbits = fpi->nbits;
462  for(s = s00;; s++)
463  {
464  {
465  switch(*s)
466  {
467  case '-':
468  sign = 1;
469  /* no break */
470  case '+':
471  if(*++s)
472  {
473  {
474  goto break2;
475  }
476  }
477  /* no break */
478  case 0:
479  sign = 0;
480  irv = STRTOG_NoNumber;
481  s = s00;
482  goto ret;
483  case '\t':
484  case '\n':
485  case '\v':
486  case '\f':
487  case '\r':
488  case ' ':
489  continue;
490  default:
491  goto break2;
492  }
493  }
494  }
495 break2:
496  if(*s == '0')
497  {
498 #ifndef NO_HEX_FP
499  switch(s[1])
500  {
501  case 'x':
502  case 'X':
503  irv = gethex(&s, fpi, exp, &rvb, sign);
504  if(irv == STRTOG_NoNumber)
505  {
506  s = s00;
507  sign = 0;
508  }
509  goto ret;
510  }
511 #endif
512  nz0 = 1;
513  while(*++s == '0')
514  {
515  {
516  ;
517  }
518  }
519  if(!*s)
520  {
521  {
522  goto ret;
523  }
524  }
525  }
526  sudden_underflow = fpi->sudden_underflow;
527  s0 = s;
528  y = z = 0;
529  for(decpt = nd = nf = 0; (c = *s) >= '0' && c <= '9'; nd++, s++)
530  {
531  {
532  if(nd < 9)
533  {
534  y = 10 * y + (unsigned)c - '0';
535  }
536  else if(nd < 16)
537  {
538  z = 10 * z + (unsigned)c - '0';
539  }
540  }
541  }
542  nd0 = nd;
543 #ifdef USE_LOCALE
544  if(c == *localeconv()->decimal_point)
545 #else
546  if(c == '.')
547 #endif
548  {
549  decpt = 1;
550  c = *++s;
551  if(!nd)
552  {
553  for(; c == '0'; c = *++s)
554  {
555  {
556  nz++;
557  }
558  }
559  if(c > '0' && c <= '9')
560  {
561  s0 = s;
562  nf += nz;
563  nz = 0;
564  goto have_dig;
565  }
566  goto dig_done;
567  }
568  for(; c >= '0' && c <= '9'; c = *++s)
569  {
570  have_dig:
571  nz++;
572  if(c -= '0')
573  {
574  nf += nz;
575  for(i = 1; i < nz; i++)
576  {
577  {
578  if(nd++ < 9)
579  {
580  y *= 10;
581  }
582  else if(nd <= DBL_DIG + 1)
583  {
584  z *= 10;
585  }
586  }
587  }
588  if(nd++ < 9)
589  {
590  y = 10 * y + (unsigned)c;
591  }
592  else if(nd <= DBL_DIG + 1)
593  {
594  z = 10 * z + (unsigned)c;
595  }
596  nz = 0;
597  }
598  }
599  }
600 dig_done:
601  e = 0;
602  if(c == 'e' || c == 'E')
603  {
604  if(!nd && !nz && !nz0)
605  {
606  irv = STRTOG_NoNumber;
607  s = s00;
608  goto ret;
609  }
610  s00 = s;
611  esign = 0;
612  switch(c = *++s)
613  {
614  case '-':
615  esign = 1;
616  case '+':
617  c = *++s;
618  }
619  if(c >= '0' && c <= '9')
620  {
621  while(c == '0')
622  {
623  {
624  c = *++s;
625  }
626  }
627  if(c > '0' && c <= '9')
628  {
629  L = c - '0';
630  s1 = s;
631  while((c = *++s) >= '0' && c <= '9')
632  {
633  {
634  L = 10 * L + c - '0';
635  }
636  }
637  if(s - s1 > 8 || L > 19999)
638  {
639  {
640  /* Avoid confusion from exponents
641  * so large that e might overflow.
642  */
643  e = 19999; /* safe for 16 bit ints */
644  }
645  }
646  else
647  {
648  {
649  e = (int)L;
650  }
651  }
652  if(esign)
653  {
654  {
655  e = -e;
656  }
657  }
658  }
659  else
660  {
661  {
662  e = 0;
663  }
664  }
665  }
666  else
667  {
668  {
669  s = s00;
670  }
671  }
672  }
673  if(!nd)
674  {
675  if(!nz && !nz0)
676  {
677 #ifdef INFNAN_CHECK
678  /* Check for Nan and Infinity */
679  if(!decpt)
680  switch(c)
681  {
682  case 'i':
683  case 'I':
684  if(match(&s, "nf"))
685  {
686  --s;
687  if(!match(&s, "inity"))
688  ++s;
689  irv = STRTOG_Infinite;
690  goto infnanexp;
691  }
692  break;
693  case 'n':
694  case 'N':
695  if(match(&s, "an"))
696  {
697  irv = STRTOG_NaN;
698  *exp = fpi->emax + 1;
699 #ifndef No_Hex_NaN
700  if(*s == '(') /*)*/
701  irv = hexnan(&s, fpi, bits);
702 #endif
703  goto infnanexp;
704  }
705  }
706 #endif /* INFNAN_CHECK */
707  irv = STRTOG_NoNumber;
708  s = s00;
709  }
710  goto ret;
711  }
712 
713  irv = STRTOG_Normal;
714  e1 = e -= nf;
715  rd = 0;
716  switch(fpi->rounding & 3)
717  {
718  case FPI_Round_up:
719  rd = 2 - sign;
720  break;
721  case FPI_Round_zero:
722  rd = 1;
723  break;
724  case FPI_Round_down:
725  rd = 1 + sign;
726  }
727 
728  /* Now we have nd0 digits, starting at s0, followed by a
729  * decimal point, followed by nd-nd0 digits. The number we're
730  * after is the integer represented by those digits times
731  * 10**e */
732 
733  if(!nd0)
734  {
735  {
736  nd0 = nd;
737  }
738  }
739  k = nd < DBL_DIG + 1 ? nd : DBL_DIG + 1;
740  dval(rv) = y;
741  if(k > 9)
742  {
743  {
744  dval(rv) = tens[k - 9] * dval(rv) + z;
745  }
746  }
747  bd0 = 0;
748  if(nbits <= P && nd <= DBL_DIG)
749  {
750  if(!e)
751  {
752  if(rvOK(dval(rv), fpi, exp, bits, 1, rd, &irv))
753  {
754  {
755  goto ret;
756  }
757  }
758  }
759  else if(e > 0)
760  {
761  if(e <= Ten_pmax)
762  {
763 #ifdef VAX
764  goto vax_ovfl_check;
765 #else
766  i = fivesbits[e] + mantbits(dval(rv)) <= P;
767  /* rv = */ rounded_product(dval(rv), tens[e]);
768  if(rvOK(dval(rv), fpi, exp, bits, i, rd, &irv))
769  {
770  {
771  goto ret;
772  }
773  }
774  e1 -= e;
775  goto rv_notOK;
776 #endif
777  }
778  i = DBL_DIG - nd;
779  if(e <= Ten_pmax + i)
780  {
781  /* A fancier test would sometimes let us do
782  * this for larger i values.
783  */
784  e2 = e - i;
785  e1 -= i;
786  dval(rv) *= tens[i];
787 #ifdef VAX
788  /* VAX exponent range is so narrow we must
789  * worry about overflow here...
790  */
791  vax_ovfl_check:
792  dval(adj) = dval(rv);
793  word0(adj) -= P * Exp_msk1;
794  /* adj = */ rounded_product(dval(adj), tens[e2]);
795  if((word0(adj) & Exp_mask) > Exp_msk1 * (DBL_MAX_EXP + Bias - 1 - P))
796  goto rv_notOK;
797  word0(adj) += P * Exp_msk1;
798  dval(rv) = dval(adj);
799 #else
800  /* rv = */ rounded_product(dval(rv), tens[e2]);
801 #endif
802  if(rvOK(dval(rv), fpi, exp, bits, 0, rd, &irv))
803  {
804  {
805  goto ret;
806  }
807  }
808  e1 -= e2;
809  }
810  }
811 #ifndef Inaccurate_Divide
812  else if(e >= -Ten_pmax)
813  {
814  /* rv = */ rounded_quotient(dval(rv), tens[-e]);
815  if(rvOK(dval(rv), fpi, exp, bits, 0, rd, &irv))
816  {
817  {
818  goto ret;
819  }
820  }
821  e1 -= e;
822  }
823 #endif
824  }
825 rv_notOK:
826  e1 += nd - k;
827 
828  /* Get starting approximation = rv * 10**e1 */
829 
830  e2 = 0;
831  if(e1 > 0)
832  {
833  if((i = e1 & 15) != 0)
834  {
835  {
836  dval(rv) *= tens[i];
837  }
838  }
839  if(e1 &= ~15)
840  {
841  e1 >>= 4;
842  while(e1 >= (1 << (n_bigtens - 1)))
843  {
844  e2 += ((word0(rv) & Exp_mask) >> Exp_shift1) - Bias;
845  word0(rv) &= (unsigned)~Exp_mask;
846  word0(rv) |= Bias << Exp_shift1;
847  dval(rv) *= bigtens[n_bigtens - 1];
848  e1 -= 1 << (n_bigtens - 1);
849  }
850  e2 += ((word0(rv) & Exp_mask) >> Exp_shift1) - Bias;
851  word0(rv) &= (unsigned)~Exp_mask;
852  word0(rv) |= Bias << Exp_shift1;
853  for(j = 0; e1 > 0; j++, e1 >>= 1)
854  {
855  {
856  if(e1 & 1)
857  {
858  {
859  dval(rv) *= bigtens[j];
860  }
861  }
862  }
863  }
864  }
865  }
866  else if(e1 < 0)
867  {
868  e1 = -e1;
869  if((i = e1 & 15) != 0)
870  {
871  {
872  dval(rv) /= tens[i];
873  }
874  }
875  if(e1 &= ~15)
876  {
877  e1 >>= 4;
878  while(e1 >= (1 << (n_bigtens - 1)))
879  {
880  e2 += ((word0(rv) & Exp_mask) >> Exp_shift1) - Bias;
881  word0(rv) &= (unsigned)~Exp_mask;
882  word0(rv) |= Bias << Exp_shift1;
883  dval(rv) *= tinytens[n_bigtens - 1];
884  e1 -= 1 << (n_bigtens - 1);
885  }
886  e2 += ((word0(rv) & Exp_mask) >> Exp_shift1) - Bias;
887  word0(rv) &= (unsigned)~Exp_mask;
888  word0(rv) |= Bias << Exp_shift1;
889  for(j = 0; e1 > 0; j++, e1 >>= 1)
890  {
891  {
892  if(e1 & 1)
893  {
894  {
895  dval(rv) *= tinytens[j];
896  }
897  }
898  }
899  }
900  }
901  }
902 #ifdef IBM
903  /* e2 is a correction to the (base 2) exponent of the return
904  * value, reflecting adjustments above to avoid overflow in the
905  * native arithmetic. For native IBM (base 16) arithmetic, we
906  * must multiply e2 by 4 to change from base 16 to 2.
907  */
908  e2 <<= 2;
909 #endif
910  rvb = d2b(dval(rv), &rve, &rvbits); /* rv = rvb * 2^rve */
911  rve += e2;
912  if((j = rvbits - nbits) > 0)
913  {
914  rshift(rvb, j);
915  rvbits = nbits;
916  rve += j;
917  }
918  bb0 = 0; /* trailing zero bits in rvb */
919  e2 = rve + rvbits - nbits;
920  if(e2 > fpi->emax + 1)
921  {
922  {
923  goto huge;
924  }
925  }
926  rve1 = rve + rvbits - nbits;
927  if(e2 < (emin = fpi->emin))
928  {
929  denorm = 1;
930  j = rve - emin;
931  if(j > 0)
932  {
933  rvb = lshift(rvb, j);
934  rvbits += j;
935  }
936  else if(j < 0)
937  {
938  rvbits += j;
939  if(rvbits <= 0)
940  {
941  if(rvbits < -1)
942  {
943  ufl:
944  rvb->wds = 0;
945  rvb->x[0] = 0;
946  *exp = emin;
948  goto ret;
949  }
950  rvb->x[0] = rvb->wds = rvbits = 1;
951  }
952  else
953  {
954  {
955  rshift(rvb, -j);
956  }
957  }
958  }
959  rve = rve1 = emin;
960  if(sudden_underflow && e2 + 1 < emin)
961  {
962  {
963  goto ufl;
964  }
965  }
966  }
967 
968  /* Now the hard part -- adjusting rv to the correct value.*/
969 
970  /* Put digits into bd: true value = bd * 10^e */
971 
972  bd0 = s2b(s0, nd0, nd, y);
973 
974  for(;;)
975  {
976  bd = Balloc(bd0->k);
977  Bcopy(bd, bd0);
978  bb = Balloc(rvb->k);
979  Bcopy(bb, rvb);
980  bbbits = rvbits - bb0;
981  bbe = rve + bb0;
982  bs = i2b(1);
983 
984  if(e >= 0)
985  {
986  bb2 = bb5 = 0;
987  bd2 = bd5 = e;
988  }
989  else
990  {
991  bb2 = bb5 = -e;
992  bd2 = bd5 = 0;
993  }
994  if(bbe >= 0)
995  {
996  {
997  bb2 += bbe;
998  }
999  }
1000  else
1001  {
1002  {
1003  bd2 -= bbe;
1004  }
1005  }
1006  bs2 = bb2;
1007  j = nbits + 1 - bbbits;
1008  i = bbe + bbbits - nbits;
1009  if(i < emin)
1010  {
1011  { /* denormal */
1012  j += i - emin;
1013  }
1014  }
1015  bb2 += j;
1016  bd2 += j;
1017  i = bb2 < bd2 ? bb2 : bd2;
1018  if(i > bs2)
1019  {
1020  {
1021  i = bs2;
1022  }
1023  }
1024  if(i > 0)
1025  {
1026  bb2 -= i;
1027  bd2 -= i;
1028  bs2 -= i;
1029  }
1030  if(bb5 > 0)
1031  {
1032  bs = pow5mult(bs, bb5);
1033  bb1 = mult(bs, bb);
1034  Bfree(bb);
1035  bb = bb1;
1036  }
1037  bb2 -= bb0;
1038  if(bb2 > 0)
1039  {
1040  {
1041  bb = lshift(bb, bb2);
1042  }
1043  }
1044  else if(bb2 < 0)
1045  {
1046  {
1047  rshift(bb, -bb2);
1048  }
1049  }
1050  if(bd5 > 0)
1051  {
1052  {
1053  bd = pow5mult(bd, bd5);
1054  }
1055  }
1056  if(bd2 > 0)
1057  {
1058  {
1059  bd = lshift(bd, bd2);
1060  }
1061  }
1062  if(bs2 > 0)
1063  {
1064  {
1065  bs = lshift(bs, bs2);
1066  }
1067  }
1068  asub = 1;
1069  inex = STRTOG_Inexhi;
1070  delta = diff(bb, bd);
1071  if(delta->wds <= 1 && !delta->x[0])
1072  {
1073  {
1074  break;
1075  }
1076  }
1077  dsign = delta->sign;
1078  delta->sign = finished = 0;
1079  L = 0;
1080  i = cmp(delta, bs);
1081  if(rd && i <= 0)
1082  {
1083  irv = STRTOG_Normal;
1084  if((finished = dsign ^ (rd & 1)) != 0)
1085  {
1086  if(dsign != 0)
1087  {
1088  irv |= STRTOG_Inexhi;
1089  goto adj1;
1090  }
1091  irv |= STRTOG_Inexlo;
1092  if(rve1 == emin)
1093  {
1094  {
1095  goto adj1;
1096  }
1097  }
1098  for(i = 0, j = nbits; j >= ULbits; i++, j -= ULbits)
1099  {
1100  if(rvb->x[i] & ALL_ON)
1101  {
1102  {
1103  goto adj1;
1104  }
1105  }
1106  }
1107  if(j > 1 && lo0bits(rvb->x + i) < j - 1)
1108  {
1109  {
1110  goto adj1;
1111  }
1112  }
1113  rve = rve1 - 1;
1114  rvb = set_ones(rvb, rvbits = nbits);
1115  break;
1116  }
1117  irv |= dsign ? STRTOG_Inexlo : STRTOG_Inexhi;
1118  break;
1119  }
1120  if(i < 0)
1121  {
1122  /* Error is less than half an ulp -- check for
1123  * special case of mantissa a power of two.
1124  */
1126  if(dsign || bbbits > 1 || denorm || rve1 == emin)
1127  {
1128  {
1129  break;
1130  }
1131  }
1132  delta = lshift(delta, 1);
1133  if(cmp(delta, bs) > 0)
1134  {
1135  irv = STRTOG_Normal | STRTOG_Inexlo;
1136  goto drop_down;
1137  }
1138  break;
1139  }
1140  if(i == 0)
1141  {
1142  /* exactly half-way between */
1143  if(dsign)
1144  {
1145  if(denorm && all_on(rvb, rvbits))
1146  {
1147  /*boundary case -- increment exponent*/
1148  rvb->wds = 1;
1149  rvb->x[0] = 1;
1150  rve = emin + nbits - (rvbits = 1);
1151  irv = STRTOG_Normal | STRTOG_Inexhi;
1152  denorm = 0;
1153  break;
1154  }
1155  irv = STRTOG_Normal | STRTOG_Inexlo;
1156  }
1157  else if(bbbits == 1)
1158  {
1159  irv = STRTOG_Normal;
1160  drop_down:
1161  /* boundary case -- decrement exponent */
1162  if(rve1 == emin)
1163  {
1164  irv = STRTOG_Normal | STRTOG_Inexhi;
1165  if(rvb->wds == 1 && rvb->x[0] == 1)
1166  {
1167  {
1168  sudden_underflow = 1;
1169  }
1170  }
1171  break;
1172  }
1173  rve -= nbits;
1174  rvb = set_ones(rvb, rvbits = nbits);
1175  break;
1176  }
1177  else
1178  {
1179  {
1180  irv = STRTOG_Normal | STRTOG_Inexhi;
1181  }
1182  }
1183  if(((bbbits < nbits) && !denorm) || !(rvb->x[0] & 1))
1184  {
1185  {
1186  break;
1187  }
1188  }
1189  if(dsign)
1190  {
1191  rvb = increment(rvb);
1192  j = kmask & (ULbits - (rvbits & kmask));
1193  if(hi0bits(rvb->x[rvb->wds - 1]) != j)
1194  {
1195  {
1196  rvbits++;
1197  }
1198  }
1199  irv = STRTOG_Normal | STRTOG_Inexhi;
1200  }
1201  else
1202  {
1203  if(bbbits == 1)
1204  {
1205  {
1206  goto undfl;
1207  }
1208  }
1209  decrement(rvb);
1210  irv = STRTOG_Normal | STRTOG_Inexlo;
1211  }
1212  break;
1213  }
1214  if((dval(adj) = ratio(delta, bs)) <= 2.)
1215  {
1216  adj1:
1217  inex = STRTOG_Inexlo;
1218  if(dsign)
1219  {
1220  asub = 0;
1221  inex = STRTOG_Inexhi;
1222  }
1223  else if(denorm && bbbits <= 1)
1224  {
1225  undfl:
1226  rvb->wds = 0;
1227  rve = emin;
1229  break;
1230  }
1231  adj0 = dval(adj) = 1.;
1232  }
1233  else
1234  {
1235  adj0 = dval(adj) *= 0.5;
1236  if(dsign)
1237  {
1238  asub = 0;
1239  inex = STRTOG_Inexlo;
1240  }
1241  if(dval(adj) < 2147483647.)
1242  {
1243  L = (int)adj0;
1244  adj0 -= L;
1245  switch(rd)
1246  {
1247  case 0:
1248  if(adj0 >= .5)
1249  {
1250  {
1251  goto inc_L;
1252  }
1253  }
1254  break;
1255  case 1:
1256  if(asub && adj0 > 0.)
1257  {
1258  {
1259  goto inc_L;
1260  }
1261  }
1262  break;
1263  case 2:
1264  if(!asub && adj0 > 0.)
1265  {
1266  inc_L:
1267  L++;
1268  inex = STRTOG_Inexact - inex;
1269  }
1270  }
1271  dval(adj) = L;
1272  }
1273  }
1274  y = (ULong)(rve + rvbits);
1275 
1276  /* adj *= ulp(dval(rv)); */
1277  /* if (asub) rv -= adj; else rv += adj; */
1278 
1279  if(!denorm && rvbits < nbits)
1280  {
1281  rvb = lshift(rvb, j = nbits - rvbits);
1282  rve -= j;
1283  rvbits = nbits;
1284  }
1285  ab = d2b(dval(adj), &abe, &abits);
1286  if(abe < 0)
1287  {
1288  {
1289  rshift(ab, -abe);
1290  }
1291  }
1292  else if(abe > 0)
1293  {
1294  {
1295  ab = lshift(ab, abe);
1296  }
1297  }
1298  rvb0 = rvb;
1299  if(asub)
1300  {
1301  /* rv -= adj; */
1302  j = hi0bits(rvb->x[rvb->wds - 1]);
1303  rvb = diff(rvb, ab);
1304  k = rvb0->wds - 1;
1305  if(denorm)
1306  {
1307  {
1308  /* do nothing */;
1309  }
1310  }
1311  else if(rvb->wds <= k || hi0bits(rvb->x[k]) > hi0bits(rvb0->x[k]))
1312  {
1313  /* unlikely; can only have lost 1 high bit */
1314  if(rve1 == emin)
1315  {
1316  --rvbits;
1317  denorm = 1;
1318  }
1319  else
1320  {
1321  rvb = lshift(rvb, 1);
1322  --rve;
1323  --rve1;
1324  L = finished = 0;
1325  }
1326  }
1327  }
1328  else
1329  {
1330  rvb = sum(rvb, ab);
1331  k = rvb->wds - 1;
1332  if(k >= rvb0->wds || hi0bits(rvb->x[k]) < hi0bits(rvb0->x[k]))
1333  {
1334  if(denorm)
1335  {
1336  if(++rvbits == nbits)
1337  {
1338  {
1339  denorm = 0;
1340  }
1341  }
1342  }
1343  else
1344  {
1345  rshift(rvb, 1);
1346  rve++;
1347  rve1++;
1348  L = 0;
1349  }
1350  }
1351  }
1352  Bfree(ab);
1353  Bfree(rvb0);
1354  if(finished)
1355  {
1356  {
1357  break;
1358  }
1359  }
1360 
1361  z = (ULong)(rve + rvbits);
1362  if(y == z && L)
1363  {
1364  /* Can we stop now? */
1365  tol = dval(adj) * 5e-16; /* > max rel error */
1366  dval(adj) = adj0 - .5;
1367  if(dval(adj) < -tol)
1368  {
1369  if(adj0 > tol)
1370  {
1371  irv |= inex;
1372  break;
1373  }
1374  }
1375  else if(dval(adj) > tol && adj0 < 1. - tol)
1376  {
1377  irv |= inex;
1378  break;
1379  }
1380  }
1381  bb0 = denorm ? 0 : trailz(rvb);
1382  Bfree(bb);
1383  Bfree(bd);
1384  Bfree(bs);
1385  Bfree(delta);
1386  }
1387  if(!denorm && (j = nbits - rvbits))
1388  {
1389  if(j > 0)
1390  {
1391  {
1392  rvb = lshift(rvb, j);
1393  }
1394  }
1395  else
1396  {
1397  {
1398  rshift(rvb, -j);
1399  }
1400  }
1401  rve -= j;
1402  }
1403  *exp = rve;
1404  Bfree(bb);
1405  Bfree(bd);
1406  Bfree(bs);
1407  Bfree(bd0);
1408  Bfree(delta);
1409  if(rve > fpi->emax)
1410  {
1411  huge:
1412  rvb->wds = 0;
1414 #ifndef NO_ERRNO
1415  errno = ERANGE;
1416 #endif
1417 #ifdef INFNAN_CHECK
1418  infnanexp:
1419  *exp = fpi->emax + 1;
1420 #endif
1421  }
1422 ret:
1423  if(denorm)
1424  {
1425  if(sudden_underflow)
1426  {
1427  rvb->wds = 0;
1429  }
1430  else
1431  {
1432  irv = (irv & ~STRTOG_Retmask) | (rvb->wds > 0 ? STRTOG_Denormal : STRTOG_Zero);
1433  if(irv & STRTOG_Inexact)
1434  {
1435  irv |= STRTOG_Underflow;
1436  }
1437  }
1438  }
1439  if(se)
1440  {
1441  *se = (char*)(uintptr_t)s;
1442  }
1443  if(sign)
1444  {
1445  irv |= STRTOG_Neg;
1446  }
1447  if(rvb)
1448  {
1449  copybits(bits, nbits, rvb);
1450  Bfree(rvb);
1451  }
1452  return irv;
1453 }
#define ALL_ON
Definition: gdtoaimp.h:469
int errno
#define Ten_pmax
Definition: gdtoaimp.h:406
int hexnan(CONST char **sp, FPI *fpi, ULong *x0)
Definition: hexnan.c:61
#define bigtens
Definition: gdtoaimp.h:514
#define DBL_MAX_EXP
Definition: float.h:42
#define P
Definition: gdtoaimp.h:399
ULong x[1]
Definition: gdtoaimp.h:488
Bigint * i2b(int i)
Definition: misc.c:271
#define Flt_Rounds
Definition: gdtoaimp.h:393
ULong any_on(Bigint *b, int k)
Definition: smisc.c:210
static int all_on(Bigint *b, int n)
Definition: strtodg.c:158
int wds
Definition: gdtoaimp.h:487
#define Bias
Definition: gdtoaimp.h:400
#define Exp_shift1
Definition: gdtoaimp.h:395
#define CONST
Definition: gdtoa.h:61
int trailz(Bigint *b)
Definition: gmisc.c:92
Bigint * lshift(Bigint *b, int k)
Definition: misc.c:495
#define Exp_msk1
Definition: gdtoaimp.h:396
#define word1(x)
Definition: gdtoaimp.h:305
int sign
Definition: gdtoaimp.h:487
Bigint * diff(Bigint *a, Bigint *b)
Definition: misc.c:617
#define Bcopy(x, y)
Definition: gdtoaimp.h:501
#define tens
Definition: gdtoaimp.h:544
#define match
Definition: gdtoaimp.h:530
#define dval(x)
Definition: gdtoaimp.h:307
#define ERANGE
Definition: errno.h:50
Bigint * increment(Bigint *b)
Definition: strtodg.c:73
int lo0bits(ULong *y)
Definition: misc.c:108
#define word0(x)
Definition: gdtoaimp.h:304
struct lconv * localeconv(void)
#define hi0bits(x)
Definition: gdtoaimp.h:525
int gethex(CONST char **sp, FPI *fpi, Long *exp, Bigint **bp, int sign)
Definition: gethex.c:46
int maxwds
Definition: gdtoaimp.h:487
Bigint * d2b(double d, int *e, int *bits)
Definition: misc.c:798
#define rounded_product(a, b)
Definition: gdtoaimp.h:433
int decrement(Bigint *b)
Definition: strtodg.c:122
#define kshift
Definition: gdtoaimp.h:467
#define rounded_quotient(a, b)
Definition: gdtoaimp.h:434
double ratio(Bigint *a, Bigint *b)
Definition: smisc.c:96
unsigned Long ULong
Definition: gdtoa.h:41
void copybits(ULong *c, int n, Bigint *b)
Definition: smisc.c:171
void Bfree(Bigint *v)
Definition: misc.c:92
#define tinytens
Definition: gdtoaimp.h:546
Bigint * s2b(CONST char *s, int nd0, int nd, ULong y9)
Definition: smisc.c:40
static CONST int fivesbits[]
Definition: strtodg.c:39
void rshift(Bigint *b, int k)
Definition: gmisc.c:39
#define kmask
Definition: gdtoaimp.h:468
#define Exp_mask
Definition: gdtoaimp.h:398
Bigint * pow5mult(Bigint *b, int k)
Definition: misc.c:420
Bigint * set_ones(Bigint *b, int n)
Definition: strtodg.c:190
#define Exp_msk11
Definition: gdtoaimp.h:397
static int rvOK(double d, FPI *fpi, Long *exp, ULong *bits, int exact, int rd, int *irv)
Definition: strtodg.c:235
static int mantbits(double d)
Definition: strtodg.c:414
int cmp(Bigint *a, Bigint *b)
Definition: misc.c:570
Definition: gdtoa.h:86
Bigint * sum(Bigint *a, Bigint *b)
Definition: sum.c:39
Bigint * Balloc(int k)
Definition: misc.c:47
#define DBL_DIG
Definition: float.h:45
Bigint * mult(Bigint *a, Bigint *b)
Definition: misc.c:287
#define Long
Definition: gdtoa.h:38
int k
Definition: gdtoaimp.h:487
#define ULbits
Definition: gdtoaimp.h:466
char * decimal_point
Definition: locale.h:20
int strtodg(CONST char *s00, char **se, FPI *fpi, Long *exp, ULong *bits)
Definition: strtodg.c:443