← Back to branch summary

~drizzle-trunk/drizzle/development

« back to all changes in this revision

Viewing changes to strings/dtoa.c

  • Committer: Brian Aker
  • Date: 2008-07-20 05:41:52 UTC
  • Revision ID: brian@tangent.org-20080720054152-5laf6plsb0o7h6ss
Documentation cleanup

Show diffs side-by-side

added added

removed removed

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