← Back to branch summary

~drizzle-trunk/drizzle/development

« back to all changes in this revision

Viewing changes to mystrings/dtoa.c

  • Committer: Mark Atwood
  • Date: 2008-10-16 11:33:16 UTC
  • mto: (520.1.13 drizzle)
  • mto: This revision was merged to the branch mainline in revision 530.
  • Revision ID: mark@fallenpegasus.com-20081016113316-ff6jdt31ck90sjdh
an implemention of the errmsg plugin

Show diffs side-by-side

added added

removed removed

Lines of Context:
mystrings/dtoa.c
11
11
 
12
12
   You should have received a copy of the GNU General Public License
13
13
   along with this program; if not, write to the Free Software
14
 
   Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301 USA */
 
14
   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA */
15
15
 
16
16
/****************************************************************
17
17
 
20
20
 
21
21
  The author of this software is David M. Gay.
22
22
 
23
 
  Copyright (C) 1991, 2000, 2001 by Lucent Technologies.
 
23
  Copyright (c) 1991, 2000, 2001 by Lucent Technologies.
24
24
 
25
25
  Permission to use, copy, modify, and distribute this software for any
26
26
  purpose without fee is hereby granted, provided that this entire notice
35
35
 
36
36
 ***************************************************************/
37
37
 
38
 
#include <config.h>
39
 
 
40
 
#include <drizzled/internal/m_string.h>  /* for memcpy and NOT_FIXED_DEC */
41
 
 
42
 
#include <float.h>
43
 
 
44
 
#include <cstdlib>
45
 
#include <cerrno>
46
 
#include <algorithm>
47
 
 
48
 
using namespace std;
49
 
 
50
 
namespace drizzled
51
 
{
52
 
namespace internal
53
 
{
 
38
#include <m_string.h>  /* for memcpy and NOT_FIXED_DEC */
 
39
 
 
40
/**
 
41
   Appears to suffice to not call malloc() in most cases.
 
42
   @todo
 
43
     see if it is possible to get rid of malloc().
 
44
     this constant is sufficient to avoid malloc() on all inputs I have tried.
 
45
*/
 
46
#define DTOA_BUFF_SIZE (420 * sizeof(void *))
54
47
 
55
48
/* Magic value returned by dtoa() to indicate overflow */
56
49
#define DTOA_OVERFLOW 9999
57
50
 
58
 
static double my_strtod_int(const char *, char **, int *);
59
 
static char *dtoa(double, int, int, int *, int *, char **);
 
51
static double my_strtod_int(const char *, char **, int *, char *, size_t);
 
52
static char *dtoa(double, int, int, int *, int *, char **, char *, size_t);
 
53
static void dtoa_free(char *, char *, size_t);
60
54
 
61
55
/**
62
56
   @brief
93
87
{
94
88
  int decpt, sign, len, i;
95
89
  char *res, *src, *end, *dst= to;
 
90
  char buf[DTOA_BUFF_SIZE];
96
91
  assert(precision >= 0 && precision < NOT_FIXED_DEC && to != NULL);
97
 
 
98
 
  res= dtoa(x, 5, precision, &decpt, &sign, &end);
 
92
  
 
93
  res= dtoa(x, 5, precision, &decpt, &sign, &end, buf, sizeof(buf));
99
94
 
100
95
  if (decpt == DTOA_OVERFLOW)
101
96
  {
102
 
    free(res);
 
97
    dtoa_free(res, buf, sizeof(buf));
103
98
    *to++= '0';
104
99
    *to= '\0';
105
100
    if (error != NULL)
134
129
  {
135
130
    if (len <= decpt)
136
131
      *dst++= '.';
137
 
 
138
 
    for (i= precision - max(0, (len - decpt)); i > 0; i--)
 
132
    
 
133
    for (i= precision - cmax(0, (len - decpt)); i > 0; i--)
139
134
      *dst++= '0';
140
135
  }
141
 
 
 
136
  
142
137
  *dst= '\0';
143
138
  if (error != NULL)
144
139
    *error= false;
145
140
 
146
 
  free(res);
 
141
  dtoa_free(res, buf, sizeof(buf));
147
142
 
148
143
  return dst - to;
149
144
}
201
196
     my_gcvt(55, ..., 1, ...);
202
197
 
203
198
   We do our best to minimize such cases by:
204
 
 
 
199
   
205
200
   - passing to dtoa() the field width as the number of significant digits
206
 
 
 
201
   
207
202
   - removing the sign of the number early (and decreasing the width before
208
203
     passing it to dtoa())
209
 
 
 
204
   
210
205
   - choosing the proper format to preserve the most number of significant
211
206
     digits.
212
207
*/
216
211
{
217
212
  int decpt, sign, len, exp_len;
218
213
  char *res, *src, *end, *dst= to, *dend= dst + width;
 
214
  char buf[DTOA_BUFF_SIZE];
219
215
  bool have_space, force_e_format;
220
216
  assert(width > 0 && to != NULL);
221
 
 
 
217
  
222
218
  /* We want to remove '-' from equations early */
223
219
  if (x < 0.)
224
220
    width--;
225
221
 
226
 
  res= dtoa(x, 4, type == MY_GCVT_ARG_DOUBLE ? min(width, DBL_DIG) : 
227
 
            min(width, FLT_DIG), &decpt, &sign, &end);
228
 
 
 
222
  res= dtoa(x, 4, type == MY_GCVT_ARG_DOUBLE ? width : cmin(width, FLT_DIG),
 
223
            &decpt, &sign, &end, buf, sizeof(buf));
229
224
  if (decpt == DTOA_OVERFLOW)
230
225
  {
231
 
    free(res);
 
226
    dtoa_free(res, buf, sizeof(buf));
232
227
    *to++= '0';
233
228
    *to= '\0';
234
229
    if (error != NULL)
248
243
     to count it here.
249
244
   */
250
245
  exp_len= 1 + (decpt >= 101 || decpt <= -99) + (decpt >= 11 || decpt <= -9);
251
 
 
 
246
  
252
247
  /*
253
248
     Do we have enough space for all digits in the 'f' format?
254
249
     Let 'len' be the number of significant digits returned by dtoa,
301
296
       ((decpt <= width && (decpt >= -1 || (decpt == -2 &&
302
297
                                            (len > 1 || !force_e_format)))) &&
303
298
         !force_e_format)) &&
304
 
 
 
299
      
305
300
       /*
306
301
         Use the 'e' format in some cases even if we have enough space for the
307
302
         'f' one. See comment for MAX_DECPT_FOR_F_FORMAT.
323
318
          *error= true;
324
319
        width= decpt;
325
320
      }
326
 
 
 
321
      
327
322
      /*
328
323
        We want to truncate (len - width) least significant digits after the
329
324
        decimal point. For this we are calling dtoa with mode=5, passing the
330
325
        number of significant digits = (len-decpt) - (len-width) = width-decpt
331
326
      */
332
 
      free(res);
333
 
      res= dtoa(x, 5, width - decpt, &decpt, &sign, &end);
 
327
      dtoa_free(res, buf, sizeof(buf));
 
328
      res= dtoa(x, 5, width - decpt, &decpt, &sign, &end, buf, sizeof(buf));
334
329
      src= res;
335
330
      len= end - res;
336
331
    }
341
336
      *dst++= '0';
342
337
      goto end;
343
338
    }
344
 
 
 
339
    
345
340
    /*
346
341
      At this point we are sure we have enough space to put all digits
347
342
      returned by dtoa
390
385
        *error= true;
391
386
      width= 0;
392
387
    }
393
 
 
 
388
      
394
389
    /* Do we have to truncate any digits? */
395
390
    if (width < len)
396
391
    {
397
392
      /* Yes, re-convert with a smaller width */
398
 
      free(res);
399
 
      res= dtoa(x, 4, width, &decpt, &sign, &end);
 
393
      dtoa_free(res, buf, sizeof(buf));
 
394
      res= dtoa(x, 4, width, &decpt, &sign, &end, buf, sizeof(buf));
400
395
      src= res;
401
396
      len= end - res;
402
397
      if (--decpt < 0)
436
431
  }
437
432
 
438
433
end:
439
 
  free(res);
 
434
  dtoa_free(res, buf, sizeof(buf));
440
435
  *dst= '\0';
441
436
 
442
437
  return dst - to;
455
450
                  rejected character.
456
451
   @param error   Upon return is set to EOVERFLOW in case of underflow or
457
452
                  overflow.
458
 
 
 
453
   
459
454
   @return        The resulting double value. In case of underflow, 0.0 is
460
455
                  returned. In case overflow, signed DBL_MAX is returned.
461
456
*/
462
457
 
463
458
double my_strtod(const char *str, char **end, int *error)
464
459
{
 
460
  char buf[DTOA_BUFF_SIZE];
465
461
  double res;
466
462
  assert(str != NULL && end != NULL && *end != NULL && error != NULL);
467
463
 
468
 
  res= my_strtod_int(str, end, error);
 
464
  res= my_strtod_int(str, end, error, buf, sizeof(buf));
469
465
  return (*error == 0) ? res : (res < 0 ? -DBL_MAX : DBL_MAX);
470
466
}
471
467
 
482
478
 *
483
479
 * The author of this software is David M. Gay.
484
480
 *
485
 
 * Copyright (C) 1991, 2000, 2001 by Lucent Technologies.
 
481
 * Copyright (c) 1991, 2000, 2001 by Lucent Technologies.
486
482
 *
487
483
 * Permission to use, copy, modify, and distribute this software for any
488
484
 * purpose without fee is hereby granted, provided that this entire notice
532
528
  file.
533
529
*/
534
530
 
 
531
/*
 
532
  #define Honor_FLT_ROUNDS if FLT_ROUNDS can assume the values 2 or 3
 
533
       and dtoa should round accordingly.
 
534
  #define Check_FLT_ROUNDS if FLT_ROUNDS can assume the values 2 or 3
 
535
       and Honor_FLT_ROUNDS is not #defined.
 
536
 
 
537
  TODO: check if we can get rid of the above two
 
538
*/
 
539
 
535
540
typedef int32_t Long;
536
541
typedef uint32_t ULong;
537
542
typedef int64_t LLong;
587
592
#endif
588
593
#endif /*Flt_Rounds*/
589
594
 
 
595
#ifdef Honor_FLT_ROUNDS
 
596
#define Rounding rounding
 
597
#undef Check_FLT_ROUNDS
 
598
#define Check_FLT_ROUNDS
 
599
#else
 
600
#define Rounding Flt_Rounds
 
601
#endif
 
602
 
590
603
#define rounded_product(a,b) a*= b
591
604
#define rounded_quotient(a,b) a/= b
592
605
 
594
607
#define Big1 0xffffffff
595
608
#define FFFFFFFF 0xffffffffUL
596
609
 
 
610
/* This is tested to be enough for dtoa */
 
611
 
 
612
#define Kmax 15
 
613
 
 
614
#define Bcopy(x,y) memcpy(&x->sign, &y->sign, \
 
615
                          2*sizeof(int) + y->wds*sizeof(ULong))
 
616
 
597
617
/* Arbitrary-length integer */
598
618
 
599
619
typedef struct Bigint
608
628
  int wds;                 /* current length in 32-bit words */
609
629
} Bigint;
610
630
 
611
 
static Bigint *Bcopy(Bigint* dst, Bigint* src)
 
631
 
 
632
/* A simple stack-memory based allocator for Bigints */
 
633
 
 
634
typedef struct Stack_alloc
612
635
{
613
 
  dst->sign= src->sign;
614
 
  dst->wds= src->wds;
615
 
 
616
 
  assert(dst->maxwds >= src->wds);
617
 
 
618
 
  memcpy(dst->p.x, src->p.x, src->wds*sizeof(ULong));
619
 
 
620
 
  return dst;
621
 
}
622
 
 
623
 
static Bigint *Balloc(int k)
 
636
  char *begin;
 
637
  char *free;
 
638
  char *end;
 
639
  /*
 
640
    Having list of free blocks lets us reduce maximum required amount
 
641
    of memory from ~4000 bytes to < 1680 (tested on x86).
 
642
  */
 
643
  Bigint *freelist[Kmax+1];
 
644
} Stack_alloc;
 
645
 
 
646
 
 
647
/*
 
648
  Try to allocate object on stack, and resort to malloc if all
 
649
  stack memory is used.
 
650
*/
 
651
 
 
652
static Bigint *Balloc(int k, Stack_alloc *alloc)
624
653
{
625
654
  Bigint *rv;
626
 
 
627
 
  /* TODO: some malloc failure checking */
628
 
 
629
 
  int x= 1 << k;
630
 
  rv= (Bigint*) malloc(sizeof(Bigint));
631
 
 
632
 
  rv->p.x= (ULong*)malloc(x * sizeof(ULong));
633
 
 
634
 
  rv->k= k;
635
 
  rv->maxwds= x;
636
 
 
 
655
  if (k <= Kmax &&  alloc->freelist[k])
 
656
  {
 
657
    rv= alloc->freelist[k];
 
658
    alloc->freelist[k]= rv->p.next;
 
659
  }
 
660
  else
 
661
  {
 
662
    int x, len;
 
663
 
 
664
    x= 1 << k;
 
665
    len= sizeof(Bigint) + x * sizeof(ULong);
 
666
 
 
667
    if (alloc->free + len <= alloc->end)
 
668
    {
 
669
      rv= (Bigint*) alloc->free;
 
670
      alloc->free+= len;
 
671
    }
 
672
    else
 
673
      rv= (Bigint*) malloc(len);
 
674
 
 
675
    rv->k= k;
 
676
    rv->maxwds= x;
 
677
  }
637
678
  rv->sign= rv->wds= 0;
638
 
 
 
679
  rv->p.x= (ULong*) (rv + 1);
639
680
  return rv;
640
681
}
641
682
 
645
686
  list. Otherwise call free().
646
687
*/
647
688
 
648
 
static void Bfree(Bigint *v)
649
 
{
650
 
  if(!v)
651
 
    return;
652
 
  free(v->p.x);
653
 
  free(v);
654
 
}
 
689
static void Bfree(Bigint *v, Stack_alloc *alloc)
 
690
{
 
691
  char *gptr= (char*) v;                       /* generic pointer */
 
692
  if (gptr < alloc->begin || gptr >= alloc->end)
 
693
    free(gptr);
 
694
  else if (v->k <= Kmax)
 
695
  {
 
696
    /*
 
697
      Maintain free lists only for stack objects: this way we don't
 
698
      have to bother with freeing lists in the end of dtoa;
 
699
      heap should not be used normally anyway.
 
700
    */
 
701
    v->p.next= alloc->freelist[v->k];
 
702
    alloc->freelist[v->k]= v;
 
703
  }
 
704
}
 
705
 
 
706
 
 
707
/*
 
708
  This is to place return value of dtoa in: tries to use stack
 
709
  as well, but passes by free lists management and just aligns len by
 
710
  sizeof(ULong).
 
711
*/
 
712
 
 
713
static char *dtoa_alloc(int i, Stack_alloc *alloc)
 
714
{
 
715
  char *rv;
 
716
  int aligned_size= (i + sizeof(ULong) - 1) / sizeof(ULong) * sizeof(ULong);
 
717
  if (alloc->free + aligned_size <= alloc->end)
 
718
  {
 
719
    rv= alloc->free;
 
720
    alloc->free+= aligned_size;
 
721
  }
 
722
  else
 
723
    rv= malloc(i);
 
724
  return rv;
 
725
}
 
726
 
 
727
 
 
728
/*
 
729
  dtoa_free() must be used to free values s returned by dtoa()
 
730
  This is the counterpart of dtoa_alloc()
 
731
*/
 
732
 
 
733
static void dtoa_free(char *gptr, char *buf, size_t buf_size)
 
734
{
 
735
  if (gptr < buf || gptr >= buf + buf_size)
 
736
    free(gptr);
 
737
}
 
738
 
655
739
 
656
740
/* Bigint arithmetic functions */
657
741
 
658
742
/* Multiply by m and add a */
659
743
 
660
 
static Bigint *multadd(Bigint *b, int m, int a)
 
744
static Bigint *multadd(Bigint *b, int m, int a, Stack_alloc *alloc)
661
745
{
662
746
  int i, wds;
663
747
  ULong *x;
679
763
  {
680
764
    if (wds >= b->maxwds)
681
765
    {
682
 
      b1= Balloc(b->k+1);
 
766
      b1= Balloc(b->k+1, alloc);
683
767
      Bcopy(b1, b);
684
 
      Bfree(b);
 
768
      Bfree(b, alloc);
685
769
      b= b1;
686
770
    }
687
771
    b->p.x[wds++]= (ULong) carry;
691
775
}
692
776
 
693
777
 
694
 
static Bigint *s2b(const char *s, int nd0, int nd, ULong y9)
 
778
static Bigint *s2b(const char *s, int nd0, int nd, ULong y9, Stack_alloc *alloc)
695
779
{
696
780
  Bigint *b;
697
781
  int i, k;
699
783
 
700
784
  x= (nd + 8) / 9;
701
785
  for (k= 0, y= 1; x > y; y <<= 1, k++) ;
702
 
  b= Balloc(k);
 
786
  b= Balloc(k, alloc);
703
787
  b->p.x[0]= y9;
704
788
  b->wds= 1;
705
 
 
 
789
  
706
790
  i= 9;
707
791
  if (9 < nd0)
708
792
  {
709
793
    s+= 9;
710
794
    do
711
 
      b= multadd(b, 10, *s++ - '0');
 
795
      b= multadd(b, 10, *s++ - '0', alloc);
712
796
    while (++i < nd0);
713
797
    s++;
714
798
  }
715
799
  else
716
800
    s+= 10;
717
801
  for(; i < nd; i++)
718
 
    b= multadd(b, 10, *s++ - '0');
 
802
    b= multadd(b, 10, *s++ - '0', alloc);
719
803
  return b;
720
804
}
721
805
 
722
806
 
723
 
static int hi0bits(ULong x)
 
807
static int hi0bits(register ULong x)
724
808
{
725
 
  int k= 0;
 
809
  register int k= 0;
726
810
 
727
811
  if (!(x & 0xffff0000))
728
812
  {
756
840
 
757
841
static int lo0bits(ULong *y)
758
842
{
759
 
  int k;
760
 
  ULong x= *y;
 
843
  register int k;
 
844
  register ULong x= *y;
761
845
 
762
846
  if (x & 7)
763
847
  {
806
890
 
807
891
/* Convert integer to Bigint number */
808
892
 
809
 
static Bigint *i2b(int i)
 
893
static Bigint *i2b(int i, Stack_alloc *alloc)
810
894
{
811
895
  Bigint *b;
812
896
 
813
 
  b= Balloc(1);
 
897
  b= Balloc(1, alloc);
814
898
  b->p.x[0]= i;
815
899
  b->wds= 1;
816
900
  return b;
819
903
 
820
904
/* Multiply two Bigint numbers */
821
905
 
822
 
static Bigint *mult(Bigint *a, Bigint *b)
 
906
static Bigint *mult(Bigint *a, Bigint *b, Stack_alloc *alloc)
823
907
{
824
908
  Bigint *c;
825
909
  int k, wa, wb, wc;
839
923
  wc= wa + wb;
840
924
  if (wc > a->maxwds)
841
925
    k++;
842
 
  c= Balloc(k);
 
926
  c= Balloc(k, alloc);
843
927
  for (x= c->p.x, xa= x + wc; x < xa; x++)
844
928
    *x= 0;
845
929
  xa= a->p.x;
911
995
 
912
996
#define P5A_MAX (sizeof(p5_a)/sizeof(*p5_a) - 1)
913
997
 
914
 
static Bigint *pow5mult(Bigint *b, int k)
 
998
static Bigint *pow5mult(Bigint *b, int k, Stack_alloc *alloc)
915
999
{
916
1000
  Bigint *b1, *p5, *p51;
917
1001
  int i;
918
1002
  static int p05[3]= { 5, 25, 125 };
919
1003
 
920
1004
  if ((i= k & 3))
921
 
    b= multadd(b, p05[i-1], 0);
 
1005
    b= multadd(b, p05[i-1], 0, alloc);
922
1006
 
923
1007
  if (!(k>>= 2))
924
1008
    return b;
927
1011
  {
928
1012
    if (k & 1)
929
1013
    {
930
 
      b1= mult(b, p5);
931
 
      Bfree(b);
 
1014
      b1= mult(b, p5, alloc);
 
1015
      Bfree(b, alloc);
932
1016
      b= b1;
933
1017
    }
934
1018
    if (!(k>>= 1))
937
1021
    if (p5 < p5_a + P5A_MAX)
938
1022
      ++p5;
939
1023
    else if (p5 == p5_a + P5A_MAX)
940
 
      p5= mult(p5, p5);
 
1024
      p5= mult(p5, p5, alloc);
941
1025
    else
942
1026
    {
943
 
      p51= mult(p5, p5);
944
 
      Bfree(p5);
 
1027
      p51= mult(p5, p5, alloc);
 
1028
      Bfree(p5, alloc);
945
1029
      p5= p51;
946
1030
    }
947
1031
  }
949
1033
}
950
1034
 
951
1035
 
952
 
static Bigint *lshift(Bigint *b, int k)
 
1036
static Bigint *lshift(Bigint *b, int k, Stack_alloc *alloc)
953
1037
{
954
1038
  int i, k1, n, n1;
955
1039
  Bigint *b1;
960
1044
  n1= n + b->wds + 1;
961
1045
  for (i= b->maxwds; n1 > i; i<<= 1)
962
1046
    k1++;
963
 
  b1= Balloc(k1);
 
1047
  b1= Balloc(k1, alloc);
964
1048
  x1= b1->p.x;
965
1049
  for (i= 0; i < n; i++)
966
1050
    *x1++= 0;
984
1068
      *x1++= *x++;
985
1069
    while (x < xe);
986
1070
  b1->wds= n1 - 1;
987
 
  Bfree(b);
 
1071
  Bfree(b, alloc);
988
1072
  return b1;
989
1073
}
990
1074
 
1013
1097
}
1014
1098
 
1015
1099
 
1016
 
static Bigint *diff(Bigint *a, Bigint *b)
 
1100
static Bigint *diff(Bigint *a, Bigint *b, Stack_alloc *alloc)
1017
1101
{
1018
1102
  Bigint *c;
1019
1103
  int i, wa, wb;
1023
1107
  i= cmp(a,b);
1024
1108
  if (!i)
1025
1109
  {
1026
 
    c= Balloc(0);
 
1110
    c= Balloc(0, alloc);
1027
1111
    c->wds= 1;
1028
1112
    c->p.x[0]= 0;
1029
1113
    return c;
1037
1121
  }
1038
1122
  else
1039
1123
    i= 0;
1040
 
  c= Balloc(a->k);
 
1124
  c= Balloc(a->k, alloc);
1041
1125
  c->sign= i;
1042
1126
  wa= a->wds;
1043
1127
  xa= a->p.x;
1069
1153
 
1070
1154
static double ulp(double x)
1071
1155
{
1072
 
  Long L;
 
1156
  register Long L;
1073
1157
  double a;
1074
1158
 
1075
1159
  L= (word0(x) & Exp_mask) - (P - 1)*Exp_msk1;
1118
1202
}
1119
1203
 
1120
1204
 
1121
 
static Bigint *d2b(double d, int *e, int *bits)
 
1205
static Bigint *d2b(double d, int *e, int *bits, Stack_alloc *alloc)
1122
1206
{
1123
1207
  Bigint *b;
1124
1208
  int de, k;
1127
1211
#define d0 word0(d)
1128
1212
#define d1 word1(d)
1129
1213
 
1130
 
  b= Balloc(1);
 
1214
  b= Balloc(1, alloc);
1131
1215
  x= b->p.x;
1132
1216
 
1133
1217
  z= d0 & Frac_mask;
1199
1283
  9007199254740992.*9007199254740992.e-256 /* = 2^106 * 1e-53 */
1200
1284
};
1201
1285
/*
1202
 
  The factor of 2^53 in tinytens[4] helps us avoid setting the underflow
 
1286
  The factor of 2^53 in tinytens[4] helps us avoid setting the underflow 
1203
1287
  flag unnecessarily.  It leads to a song and dance at the end of strtod.
1204
1288
*/
1205
1289
#define Scale_Bit 0x10
1207
1291
 
1208
1292
/*
1209
1293
  strtod for IEEE--arithmetic machines.
1210
 
 
 
1294
 
1211
1295
  This strtod returns a nearest machine number to the input decimal
1212
1296
  string (or sets errno to EOVERFLOW). Ties are broken by the IEEE round-even
1213
1297
  rule.
1214
 
 
 
1298
 
1215
1299
  Inspired loosely by William D. Clinger's paper "How to Read Floating
1216
1300
  Point Numbers Accurately" [Proc. ACM SIGPLAN '90, pp. 92-101].
1217
 
 
 
1301
 
1218
1302
  Modifications:
1219
 
 
 
1303
 
1220
1304
   1. We only require IEEE (not IEEE double-extended).
1221
1305
   2. We get by with floating-point arithmetic in a case that
1222
1306
     Clinger missed -- when we're computing d * 10^n
1234
1318
     for 0 <= k <= 22).
1235
1319
*/
1236
1320
 
1237
 
static double my_strtod_int(const char *s00, char **se, int *error)
 
1321
static double my_strtod_int(const char *s00, char **se, int *error, char *buf, size_t buf_size)
1238
1322
{
1239
1323
  int scale;
1240
1324
  int bb2, bb5, bbe, bd2, bd5, bbbits, bs2, c, dsign,
1247
1331
#ifdef SET_INEXACT
1248
1332
  int inexact, oldinexact;
1249
1333
#endif
 
1334
#ifdef Honor_FLT_ROUNDS
 
1335
  int rounding;
 
1336
#endif
 
1337
  Stack_alloc alloc;
1250
1338
 
1251
1339
  c= 0;
1252
1340
  *error= 0;
1253
1341
 
 
1342
  alloc.begin= alloc.free= buf;
 
1343
  alloc.end= buf + buf_size;
 
1344
  memset(alloc.freelist, 0, sizeof(alloc.freelist));
 
1345
 
1254
1346
  sign= nz0= nz= 0;
1255
1347
  dval(rv)= 0.;
1256
1348
  for (s= s00; s < end; s++)
1274
1366
 break2:
1275
1367
  if (s >= end)
1276
1368
    goto ret0;
1277
 
 
 
1369
  
1278
1370
  if (*s == '0')
1279
1371
  {
1280
1372
    nz0= 1;
1399
1491
  }
1400
1492
  bd0= 0;
1401
1493
  if (nd <= DBL_DIG
 
1494
#ifndef Honor_FLT_ROUNDS
 
1495
    && Flt_Rounds == 1
 
1496
#endif
1402
1497
      )
1403
1498
  {
1404
1499
    if (!e)
1407
1502
    {
1408
1503
      if (e <= Ten_pmax)
1409
1504
      {
 
1505
#ifdef Honor_FLT_ROUNDS
 
1506
        /* round correctly FLT_ROUNDS = 2 or 3 */
 
1507
        if (sign)
 
1508
        {
 
1509
          rv= -rv;
 
1510
          sign= 0;
 
1511
        }
 
1512
#endif
1410
1513
        /* rv = */ rounded_product(dval(rv), tens[e]);
1411
1514
        goto ret;
1412
1515
      }
1417
1520
          A fancier test would sometimes let us do
1418
1521
          this for larger i values.
1419
1522
        */
 
1523
#ifdef Honor_FLT_ROUNDS
 
1524
        /* round correctly FLT_ROUNDS = 2 or 3 */
 
1525
        if (sign)
 
1526
        {
 
1527
          rv= -rv;
 
1528
          sign= 0;
 
1529
        }
 
1530
#endif
1420
1531
        e-= i;
1421
1532
        dval(rv)*= tens[i];
1422
1533
        /* rv = */ rounded_product(dval(rv), tens[e]);
1426
1537
#ifndef Inaccurate_Divide
1427
1538
    else if (e >= -Ten_pmax)
1428
1539
    {
 
1540
#ifdef Honor_FLT_ROUNDS
 
1541
      /* round correctly FLT_ROUNDS = 2 or 3 */
 
1542
      if (sign)
 
1543
      {
 
1544
        rv= -rv;
 
1545
        sign= 0;
 
1546
      }
 
1547
#endif
1429
1548
      /* rv = */ rounded_quotient(dval(rv), tens[-e]);
1430
1549
      goto ret;
1431
1550
    }
1439
1558
    oldinexact= get_inexact();
1440
1559
#endif
1441
1560
  scale= 0;
 
1561
#ifdef Honor_FLT_ROUNDS
 
1562
  if ((rounding= Flt_Rounds) >= 2)
 
1563
  {
 
1564
    if (sign)
 
1565
      rounding= rounding == 2 ? 0 : 2;
 
1566
    else
 
1567
      if (rounding != 2)
 
1568
        rounding= 0;
 
1569
  }
 
1570
#endif
1442
1571
 
1443
1572
  /* Get starting approximation = rv * 10**e1 */
1444
1573
 
1453
1582
 ovfl:
1454
1583
        *error= EOVERFLOW;
1455
1584
        /* Can't trust HUGE_VAL */
 
1585
#ifdef Honor_FLT_ROUNDS
 
1586
        switch (rounding)
 
1587
        {
 
1588
        case 0: /* toward 0 */
 
1589
        case 3: /* toward -infinity */
 
1590
          word0(rv)= Big0;
 
1591
          word1(rv)= Big1;
 
1592
          break;
 
1593
        default:
 
1594
          word0(rv)= Exp_mask;
 
1595
          word1(rv)= 0;
 
1596
        }
 
1597
#else /*Honor_FLT_ROUNDS*/
1456
1598
        word0(rv)= Exp_mask;
1457
1599
        word1(rv)= 0;
 
1600
#endif /*Honor_FLT_ROUNDS*/
1458
1601
#ifdef SET_INEXACT
1459
1602
        /* set overflow bit */
1460
1603
        dval(rv0)= 1e300;
1526
1669
 
1527
1670
  /* Put digits into bd: true value = bd * 10^e */
1528
1671
 
1529
 
  bd0= s2b(s0, nd0, nd, y);
 
1672
  bd0= s2b(s0, nd0, nd, y, &alloc);
1530
1673
 
1531
1674
  for(;;)
1532
1675
  {
1533
 
    bd= Balloc(bd0->k);
 
1676
    bd= Balloc(bd0->k, &alloc);
1534
1677
    Bcopy(bd, bd0);
1535
 
    bb= d2b(dval(rv), &bbe, &bbbits);  /* rv = bb * 2^bbe */
1536
 
    bs= i2b(1);
 
1678
    bb= d2b(dval(rv), &bbe, &bbbits, &alloc);  /* rv = bb * 2^bbe */
 
1679
    bs= i2b(1, &alloc);
1537
1680
 
1538
1681
    if (e >= 0)
1539
1682
    {
1550
1693
    else
1551
1694
      bd2-= bbe;
1552
1695
    bs2= bb2;
 
1696
#ifdef Honor_FLT_ROUNDS
 
1697
    if (rounding != 1)
 
1698
      bs2++;
 
1699
#endif
1553
1700
    j= bbe - scale;
1554
1701
    i= j + bbbits - 1;  /* logb(rv) */
1555
1702
    if (i < Emin)  /* denormal */
1570
1717
    }
1571
1718
    if (bb5 > 0)
1572
1719
    {
1573
 
      bs= pow5mult(bs, bb5);
1574
 
      bb1= mult(bs, bb);
1575
 
      Bfree(bb);
 
1720
      bs= pow5mult(bs, bb5, &alloc);
 
1721
      bb1= mult(bs, bb, &alloc);
 
1722
      Bfree(bb, &alloc);
1576
1723
      bb= bb1;
1577
1724
    }
1578
1725
    if (bb2 > 0)
1579
 
      bb= lshift(bb, bb2);
 
1726
      bb= lshift(bb, bb2, &alloc);
1580
1727
    if (bd5 > 0)
1581
 
      bd= pow5mult(bd, bd5);
 
1728
      bd= pow5mult(bd, bd5, &alloc);
1582
1729
    if (bd2 > 0)
1583
 
      bd= lshift(bd, bd2);
 
1730
      bd= lshift(bd, bd2, &alloc);
1584
1731
    if (bs2 > 0)
1585
 
      bs= lshift(bs, bs2);
1586
 
    delta= diff(bb, bd);
 
1732
      bs= lshift(bs, bs2, &alloc);
 
1733
    delta= diff(bb, bd, &alloc);
1587
1734
    dsign= delta->sign;
1588
1735
    delta->sign= 0;
1589
1736
    i= cmp(delta, bs);
 
1737
#ifdef Honor_FLT_ROUNDS
 
1738
    if (rounding != 1)
 
1739
    {
 
1740
      if (i < 0)
 
1741
      {
 
1742
        /* Error is less than an ulp */
 
1743
        if (!delta->x[0] && delta->wds <= 1)
 
1744
        {
 
1745
          /* exact */
 
1746
#ifdef SET_INEXACT
 
1747
          inexact= 0;
 
1748
#endif
 
1749
          break;
 
1750
        }
 
1751
        if (rounding)
 
1752
        {
 
1753
          if (dsign)
 
1754
          {
 
1755
            adj= 1.;
 
1756
            goto apply_adj;
 
1757
          }
 
1758
        }
 
1759
        else if (!dsign)
 
1760
        {
 
1761
          adj= -1.;
 
1762
          if (!word1(rv) && !(word0(rv) & Frac_mask))
 
1763
          {
 
1764
            y= word0(rv) & Exp_mask;
 
1765
            if (!scale || y > 2*P*Exp_msk1)
 
1766
            {
 
1767
              delta= lshift(delta,Log2P);
 
1768
              if (cmp(delta, bs) <= 0)
 
1769
              adj= -0.5;
 
1770
            }
 
1771
          }
 
1772
 apply_adj:
 
1773
          if (scale && (y= word0(rv) & Exp_mask) <= 2 * P * Exp_msk1)
 
1774
            word0(adj)+= (2 * P + 1) * Exp_msk1 - y;
 
1775
          dval(rv)+= adj * ulp(dval(rv));
 
1776
        }
 
1777
        break;
 
1778
      }
 
1779
      adj= ratio(delta, bs);
 
1780
      if (adj < 1.)
 
1781
        adj= 1.;
 
1782
      if (adj <= 0x7ffffffe)
 
1783
      {
 
1784
        /* adj = rounding ? ceil(adj) : floor(adj); */
 
1785
        y= adj;
 
1786
        if (y != adj)
 
1787
        {
 
1788
          if (!((rounding >> 1) ^ dsign))
 
1789
            y++;
 
1790
          adj= y;
 
1791
        }
 
1792
      }
 
1793
      if (scale && (y= word0(rv) & Exp_mask) <= 2 * P * Exp_msk1)
 
1794
        word0(adj)+= (2 * P + 1) * Exp_msk1 - y;
 
1795
      adj*= ulp(dval(rv));
 
1796
      if (dsign)
 
1797
        dval(rv)+= adj;
 
1798
      else
 
1799
        dval(rv)-= adj;
 
1800
      goto cont;
 
1801
    }
 
1802
#endif /*Honor_FLT_ROUNDS*/
1590
1803
 
1591
1804
    if (i < 0)
1592
1805
    {
1611
1824
#endif
1612
1825
        break;
1613
1826
      }
1614
 
      delta= lshift(delta, Log2P);
 
1827
      delta= lshift(delta, Log2P, &alloc);
1615
1828
      if (cmp(delta, bs) > 0)
1616
1829
        goto drop_down;
1617
1830
      break;
1693
1906
    {
1694
1907
      aadj*= 0.5;
1695
1908
      aadj1= dsign ? aadj : -aadj;
 
1909
#ifdef Check_FLT_ROUNDS
 
1910
      switch (Rounding)
 
1911
      {
 
1912
      case 2: /* towards +infinity */
 
1913
        aadj1-= 0.5;
 
1914
        break;
 
1915
      case 0: /* towards 0 */
 
1916
      case 3: /* towards -infinity */
 
1917
        aadj1+= 0.5;
 
1918
      }
 
1919
#else
1696
1920
      if (Flt_Rounds == 0)
1697
1921
        aadj1+= 0.5;
 
1922
#endif /*Check_FLT_ROUNDS*/
1698
1923
    }
1699
1924
    y= word0(rv) & Exp_mask;
1700
1925
 
1752
1977
      }
1753
1978
#endif
1754
1979
 cont:
1755
 
    Bfree(bb);
1756
 
    Bfree(bd);
1757
 
    Bfree(bs);
1758
 
    Bfree(delta);
 
1980
    Bfree(bb, &alloc);
 
1981
    Bfree(bd, &alloc);
 
1982
    Bfree(bs, &alloc);
 
1983
    Bfree(delta, &alloc);
1759
1984
  }
1760
1985
#ifdef SET_INEXACT
1761
1986
  if (inexact)
1785
2010
  }
1786
2011
#endif
1787
2012
 retfree:
1788
 
  Bfree(bb);
1789
 
  Bfree(bd);
1790
 
  Bfree(bs);
1791
 
  Bfree(bd0);
1792
 
  Bfree(delta);
 
2013
  Bfree(bb, &alloc);
 
2014
  Bfree(bd, &alloc);
 
2015
  Bfree(bs, &alloc);
 
2016
  Bfree(bd0, &alloc);
 
2017
  Bfree(delta, &alloc);
1793
2018
 ret:
1794
2019
  *se= (char *)s;
1795
2020
  return sign ? -dval(rv) : dval(rv);
1896
2121
 */
1897
2122
 
1898
2123
static char *dtoa(double d, int mode, int ndigits, int *decpt, int *sign,
1899
 
                  char **rve)
 
2124
                  char **rve, char *buf, size_t buf_size)
1900
2125
{
1901
2126
  /*
1902
2127
    Arguments ndigits, decpt, sign are similar to those
1908
2133
    mode:
1909
2134
          0 ==> shortest string that yields d when read in
1910
2135
                and rounded to nearest.
1911
 
 
1912
2136
          1 ==> like 0, but with Steele & White stopping rule;
1913
2137
                e.g. with IEEE P754 arithmetic , mode 0 gives
1914
2138
                1e23 whereas mode 1 gives 9.999999999999999e22.
1922
2146
          4,5 ==> similar to 2 and 3, respectively, but (in
1923
2147
                round-nearest mode) with the tests of mode 0 to
1924
2148
                possibly return a shorter string that rounds to d.
 
2149
                With IEEE arithmetic and compilation with
 
2150
                -DHonor_FLT_ROUNDS, modes 4 and 5 behave the same
 
2151
                as modes 2 and 3 when FLT_ROUNDS != 1.
1925
2152
          6-9 ==> Debugging modes similar to mode - 4:  don't try
1926
2153
                fast floating-point estimate (if applicable).
1927
2154
 
1940
2167
  Bigint *b, *b1, *delta, *mlo = NULL, *mhi, *S;
1941
2168
  double d2, ds, eps;
1942
2169
  char *s, *s0;
 
2170
#ifdef Honor_FLT_ROUNDS
 
2171
  int rounding;
 
2172
#endif
 
2173
  Stack_alloc alloc;
 
2174
  
 
2175
  alloc.begin= alloc.free= buf;
 
2176
  alloc.end= buf + buf_size;
 
2177
  memset(alloc.freelist, 0, sizeof(alloc.freelist));
1943
2178
 
1944
2179
  if (word0(d) & Sign_bit)
1945
2180
  {
1951
2186
    *sign= 0;
1952
2187
 
1953
2188
  /* If infinity, set decpt to DTOA_OVERFLOW, if 0 set it to 1 */
1954
 
  if ((((word0(d) & Exp_mask) == Exp_mask) && ((*decpt= DTOA_OVERFLOW) != 0)) ||
1955
 
      (!dval(d) && ((*decpt= 1) != 0)))
 
2189
  if (((word0(d) & Exp_mask) == Exp_mask && (*decpt= DTOA_OVERFLOW)) ||
 
2190
      (!dval(d) && (*decpt= 1)))
1956
2191
  {
1957
2192
    /* Infinity, NaN, 0 */
1958
 
    char *res= (char*) malloc(2);
 
2193
    char *res= (char*) dtoa_alloc(2, &alloc);
1959
2194
    res[0]= '0';
1960
2195
    res[1]= '\0';
1961
2196
    if (rve)
1962
2197
      *rve= res + 1;
1963
2198
    return res;
1964
2199
  }
1965
 
 
1966
 
 
1967
 
  b= d2b(dval(d), &be, &bbits);
 
2200
  
 
2201
#ifdef Honor_FLT_ROUNDS
 
2202
  if ((rounding= Flt_Rounds) >= 2)
 
2203
  {
 
2204
    if (*sign)
 
2205
      rounding= rounding == 2 ? 0 : 2;
 
2206
    else
 
2207
      if (rounding != 2)
 
2208
        rounding= 0;
 
2209
  }
 
2210
#endif
 
2211
 
 
2212
  b= d2b(dval(d), &be, &bbits, &alloc);
1968
2213
  if ((i= (int)(word0(d) >> Exp_shift1 & (Exp_mask>>Exp_shift1))))
1969
2214
  {
1970
2215
    dval(d2)= dval(d);
1976
2221
      log10(x)      =  log(x) / log(10)
1977
2222
                   ~=~ log(1.5)/log(10) + (x-1.5)/(1.5*log(10))
1978
2223
      log10(d)= (i-Bias)*log(2)/log(10) + log10(d2)
1979
 
 
 
2224
     
1980
2225
      This suggests computing an approximation k to log10(d) by
1981
 
 
 
2226
     
1982
2227
      k= (i - Bias)*0.301029995663981
1983
2228
           + ( (d2-1.5)*0.289529654602168 + 0.176091259055681 );
1984
 
 
 
2229
     
1985
2230
      We want k to be too large rather than too small.
1986
2231
      The error in the first-order Taylor series approximation
1987
2232
      is in our favor, so we just round up the constant enough
2046
2291
  if (mode < 0 || mode > 9)
2047
2292
    mode= 0;
2048
2293
 
 
2294
#ifdef Check_FLT_ROUNDS
 
2295
  try_quick= Rounding == 1;
 
2296
#else
2049
2297
  try_quick= 1;
 
2298
#endif
2050
2299
 
2051
2300
  if (mode > 5)
2052
2301
  {
2079
2328
    if (i <= 0)
2080
2329
      i= 1;
2081
2330
  }
2082
 
  s= s0= (char*)malloc(i+1); /* +1 for trailing '\0' appended
2083
 
                                at end of function */
 
2331
  s= s0= dtoa_alloc(i, &alloc);
 
2332
 
 
2333
#ifdef Honor_FLT_ROUNDS
 
2334
  if (mode > 1 && rounding != 1)
 
2335
    leftright= 0;
 
2336
#endif
2084
2337
 
2085
2338
  if (ilim >= 0 && ilim <= Quick_max && try_quick)
2086
2339
  {
2179
2432
            goto bump_up;
2180
2433
          else if (dval(d) < 0.5 - dval(eps))
2181
2434
          {
2182
 
            while (*--s == '0') {}
 
2435
            while (*--s == '0');
2183
2436
            s++;
2184
2437
            goto ret1;
2185
2438
          }
2211
2464
    {
2212
2465
      L= (Long)(dval(d) / ds);
2213
2466
      dval(d)-= L*ds;
 
2467
#ifdef Check_FLT_ROUNDS
 
2468
      /* If FLT_ROUNDS == 2, L will usually be high by 1 */
 
2469
      if (dval(d) < 0)
 
2470
      {
 
2471
        L--;
 
2472
        dval(d)+= ds;
 
2473
      }
 
2474
#endif
2214
2475
      *s++= '0' + (int)L;
2215
2476
      if (!dval(d))
2216
2477
      {
2218
2479
      }
2219
2480
      if (i == ilim)
2220
2481
      {
 
2482
#ifdef Honor_FLT_ROUNDS
 
2483
        if (mode > 1)
 
2484
        {
 
2485
          switch (rounding) {
 
2486
          case 0: goto ret1;
 
2487
          case 2: goto bump_up;
 
2488
          }
 
2489
        }
 
2490
#endif
2221
2491
        dval(d)+= dval(d);
2222
2492
        if (dval(d) > ds || (dval(d) == ds && L & 1))
2223
2493
        {
2245
2515
    i = denorm ? be + (Bias + (P-1) - 1 + 1) : 1 + P - bbits;
2246
2516
    b2+= i;
2247
2517
    s2+= i;
2248
 
    mhi= i2b(1);
 
2518
    mhi= i2b(1, &alloc);
2249
2519
  }
2250
2520
  if (m2 > 0 && s2 > 0)
2251
2521
  {
2260
2530
    {
2261
2531
      if (m5 > 0)
2262
2532
      {
2263
 
        mhi= pow5mult(mhi, m5);
2264
 
        b1= mult(mhi, b);
2265
 
        Bfree(b);
 
2533
        mhi= pow5mult(mhi, m5, &alloc);
 
2534
        b1= mult(mhi, b, &alloc);
 
2535
        Bfree(b, &alloc);
2266
2536
        b= b1;
2267
2537
      }
2268
2538
      if ((j= b5 - m5))
2269
 
        b= pow5mult(b, j);
 
2539
        b= pow5mult(b, j, &alloc);
2270
2540
    }
2271
2541
    else
2272
 
      b= pow5mult(b, b5);
 
2542
      b= pow5mult(b, b5, &alloc);
2273
2543
  }
2274
 
  S= i2b(1);
 
2544
  S= i2b(1, &alloc);
2275
2545
  if (s5 > 0)
2276
 
    S= pow5mult(S, s5);
 
2546
    S= pow5mult(S, s5, &alloc);
2277
2547
 
2278
2548
  /* Check for special case that d is a normalized power of 2. */
2279
2549
 
2280
2550
  spec_case= 0;
2281
2551
  if ((mode < 2 || leftright)
 
2552
#ifdef Honor_FLT_ROUNDS
 
2553
      && rounding == 1
 
2554
#endif
2282
2555
     )
2283
2556
  {
2284
2557
    if (!word1(d) && !(word0(d) & Bndry_mask) &&
2295
2568
  /*
2296
2569
    Arrange for convenient computation of quotients:
2297
2570
    shift left if necessary so divisor has 4 leading 0 bits.
2298
 
 
 
2571
    
2299
2572
    Perhaps we should just compute leading 28 bits of S once
2300
2573
    a nd for all and pass them and a shift to quorem, so it
2301
2574
    can do shifts and ors to compute the numerator for q.
2317
2590
    s2+= i;
2318
2591
  }
2319
2592
  if (b2 > 0)
2320
 
    b= lshift(b, b2);
 
2593
    b= lshift(b, b2, &alloc);
2321
2594
  if (s2 > 0)
2322
 
    S= lshift(S, s2);
 
2595
    S= lshift(S, s2, &alloc);
2323
2596
  if (k_check)
2324
2597
  {
2325
2598
    if (cmp(b,S) < 0)
2326
2599
    {
2327
2600
      k--;
2328
2601
      /* we botched the k estimate */
2329
 
      b= multadd(b, 10, 0);
 
2602
      b= multadd(b, 10, 0, &alloc);
2330
2603
      if (leftright)
2331
 
        mhi= multadd(mhi, 10, 0);
 
2604
        mhi= multadd(mhi, 10, 0, &alloc);
2332
2605
      ilim= ilim1;
2333
2606
    }
2334
2607
  }
2335
2608
  if (ilim <= 0 && (mode == 3 || mode == 5))
2336
2609
  {
2337
 
    if (ilim < 0 || cmp(b,S= multadd(S,5,0)) <= 0)
 
2610
    if (ilim < 0 || cmp(b,S= multadd(S,5,0, &alloc)) <= 0)
2338
2611
    {
2339
2612
      /* no digits, fcvt style */
2340
2613
no_digits:
2349
2622
  if (leftright)
2350
2623
  {
2351
2624
    if (m2 > 0)
2352
 
      mhi= lshift(mhi, m2);
 
2625
      mhi= lshift(mhi, m2, &alloc);
2353
2626
 
2354
2627
    /*
2355
2628
      Compute mlo -- check for special case that d is a normalized power of 2.
2358
2631
    mlo= mhi;
2359
2632
    if (spec_case)
2360
2633
    {
2361
 
      mhi= Balloc(mhi->k);
 
2634
      mhi= Balloc(mhi->k, &alloc);
2362
2635
      Bcopy(mhi, mlo);
2363
 
      mhi= lshift(mhi, Log2P);
 
2636
      mhi= lshift(mhi, Log2P, &alloc);
2364
2637
    }
2365
2638
 
2366
2639
    for (i= 1;;i++)
2368
2641
      dig= quorem(b,S) + '0';
2369
2642
      /* Do we yet have the shortest decimal string that will round to d? */
2370
2643
      j= cmp(b, mlo);
2371
 
      delta= diff(S, mhi);
 
2644
      delta= diff(S, mhi, &alloc);
2372
2645
      j1= delta->sign ? 1 : cmp(b, delta);
2373
 
      Bfree(delta);
 
2646
      Bfree(delta, &alloc);
2374
2647
      if (j1 == 0 && mode != 1 && !(word1(d) & 1)
 
2648
#ifdef Honor_FLT_ROUNDS
 
2649
          && rounding >= 1
 
2650
#endif
2375
2651
         )
2376
2652
      {
2377
2653
        if (dig == '9')
2387
2663
        {
2388
2664
          goto accept_dig;
2389
2665
        }
 
2666
#ifdef Honor_FLT_ROUNDS
 
2667
        if (mode > 1)
 
2668
          switch (rounding) {
 
2669
          case 0: goto accept_dig;
 
2670
          case 2: goto keep_dig;
 
2671
          }
 
2672
#endif /*Honor_FLT_ROUNDS*/
2390
2673
        if (j1 > 0)
2391
2674
        {
2392
 
          b= lshift(b, 1);
 
2675
          b= lshift(b, 1, &alloc);
2393
2676
          j1= cmp(b, S);
2394
2677
          if ((j1 > 0 || (j1 == 0 && dig & 1))
2395
2678
              && dig++ == '9')
2401
2684
      }
2402
2685
      if (j1 > 0)
2403
2686
      {
 
2687
#ifdef Honor_FLT_ROUNDS
 
2688
        if (!rounding)
 
2689
          goto accept_dig;
 
2690
#endif
2404
2691
        if (dig == '9')
2405
2692
        { /* possible if i == 1 */
2406
2693
round_9_up:
2410
2697
        *s++= dig + 1;
2411
2698
        goto ret;
2412
2699
      }
 
2700
#ifdef Honor_FLT_ROUNDS
 
2701
keep_dig:
 
2702
#endif
2413
2703
      *s++= dig;
2414
2704
      if (i == ilim)
2415
2705
        break;
2416
 
      b= multadd(b, 10, 0);
 
2706
      b= multadd(b, 10, 0, &alloc);
2417
2707
      if (mlo == mhi)
2418
 
        mlo= mhi= multadd(mhi, 10, 0);
 
2708
        mlo= mhi= multadd(mhi, 10, 0, &alloc);
2419
2709
      else
2420
2710
      {
2421
 
        mlo= multadd(mlo, 10, 0);
2422
 
        mhi= multadd(mhi, 10, 0);
 
2711
        mlo= multadd(mlo, 10, 0, &alloc);
 
2712
        mhi= multadd(mhi, 10, 0, &alloc);
2423
2713
      }
2424
2714
    }
2425
2715
  }
2433
2723
      }
2434
2724
      if (i >= ilim)
2435
2725
        break;
2436
 
      b= multadd(b, 10, 0);
 
2726
      b= multadd(b, 10, 0, &alloc);
2437
2727
    }
2438
2728
 
2439
2729
  /* Round off last digit */
2440
2730
 
2441
 
  b= lshift(b, 1);
 
2731
#ifdef Honor_FLT_ROUNDS
 
2732
  switch (rounding) {
 
2733
  case 0: goto trimzeros;
 
2734
  case 2: goto roundoff;
 
2735
  }
 
2736
#endif
 
2737
  b= lshift(b, 1, &alloc);
2442
2738
  j= cmp(b, S);
2443
2739
  if (j > 0 || (j == 0 && dig & 1))
2444
2740
  {
2454
2750
  }
2455
2751
  else
2456
2752
  {
2457
 
    while (*--s == '0') {}
 
2753
#ifdef Honor_FLT_ROUNDS
 
2754
trimzeros:
 
2755
#endif
 
2756
    while (*--s == '0');
2458
2757
    s++;
2459
2758
  }
2460
2759
ret:
2461
 
  Bfree(S);
 
2760
  Bfree(S, &alloc);
2462
2761
  if (mhi)
2463
2762
  {
2464
2763
    if (mlo && mlo != mhi)
2465
 
      Bfree(mlo);
2466
 
    Bfree(mhi);
 
2764
      Bfree(mlo, &alloc);
 
2765
    Bfree(mhi, &alloc);
2467
2766
  }
2468
2767
ret1:
2469
 
  Bfree(b);
 
2768
  Bfree(b, &alloc);
2470
2769
  *s= 0;
2471
2770
  *decpt= k + 1;
2472
2771
  if (rve)
2473
2772
    *rve= s;
2474
2773
  return s0;
2475
2774
}
2476
 
 
2477
 
} /* namespace internal */
2478
 
} /* namespace drizzled */