OpenCores
URL https://opencores.org/ocsvn/openrisc_me/openrisc_me/trunk

Subversion Repositories openrisc_me

[/] [openrisc/] [trunk/] [gnu-src/] [newlib-1.18.0/] [newlib/] [libc/] [stdlib/] [dtoa.c] - Blame information for rev 308

Go to most recent revision | Details | Compare with Previous | View Log

Line No. Rev Author Line
1 207 jeremybenn
/****************************************************************
2
 *
3
 * The author of this software is David M. Gay.
4
 *
5
 * Copyright (c) 1991 by AT&T.
6
 *
7
 * Permission to use, copy, modify, and distribute this software for any
8
 * purpose without fee is hereby granted, provided that this entire notice
9
 * is included in all copies of any software which is or includes a copy
10
 * or modification of this software and in all copies of the supporting
11
 * documentation for such software.
12
 *
13
 * THIS SOFTWARE IS BEING PROVIDED "AS IS", WITHOUT ANY EXPRESS OR IMPLIED
14
 * WARRANTY.  IN PARTICULAR, NEITHER THE AUTHOR NOR AT&T MAKES ANY
15
 * REPRESENTATION OR WARRANTY OF ANY KIND CONCERNING THE MERCHANTABILITY
16
 * OF THIS SOFTWARE OR ITS FITNESS FOR ANY PARTICULAR PURPOSE.
17
 *
18
 ***************************************************************/
19
 
20
/* Please send bug reports to
21
        David M. Gay
22
        AT&T Bell Laboratories, Room 2C-463
23
        600 Mountain Avenue
24
        Murray Hill, NJ 07974-2070
25
        U.S.A.
26
        dmg@research.att.com or research!dmg
27
 */
28
 
29
#include <_ansi.h>
30
#include <stdlib.h>
31
#include <reent.h>
32
#include <string.h>
33
#include "mprec.h"
34
 
35
static int
36
_DEFUN (quorem,
37
        (b, S),
38
        _Bigint * b _AND _Bigint * S)
39
{
40
  int n;
41
  __Long borrow, y;
42
  __ULong carry, q, ys;
43
  __ULong *bx, *bxe, *sx, *sxe;
44
#ifdef Pack_32
45
  __Long z;
46
  __ULong si, zs;
47
#endif
48
 
49
  n = S->_wds;
50
#ifdef DEBUG
51
  /*debug*/ if (b->_wds > n)
52
    /*debug*/ Bug ("oversize b in quorem");
53
#endif
54
  if (b->_wds < n)
55
    return 0;
56
  sx = S->_x;
57
  sxe = sx + --n;
58
  bx = b->_x;
59
  bxe = bx + n;
60
  q = *bxe / (*sxe + 1);        /* ensure q <= true quotient */
61
#ifdef DEBUG
62
  /*debug*/ if (q > 9)
63
    /*debug*/ Bug ("oversized quotient in quorem");
64
#endif
65
  if (q)
66
    {
67
      borrow = 0;
68
      carry = 0;
69
      do
70
        {
71
#ifdef Pack_32
72
          si = *sx++;
73
          ys = (si & 0xffff) * q + carry;
74
          zs = (si >> 16) * q + (ys >> 16);
75
          carry = zs >> 16;
76
          y = (*bx & 0xffff) - (ys & 0xffff) + borrow;
77
          borrow = y >> 16;
78
          Sign_Extend (borrow, y);
79
          z = (*bx >> 16) - (zs & 0xffff) + borrow;
80
          borrow = z >> 16;
81
          Sign_Extend (borrow, z);
82
          Storeinc (bx, z, y);
83
#else
84
          ys = *sx++ * q + carry;
85
          carry = ys >> 16;
86
          y = *bx - (ys & 0xffff) + borrow;
87
          borrow = y >> 16;
88
          Sign_Extend (borrow, y);
89
          *bx++ = y & 0xffff;
90
#endif
91
        }
92
      while (sx <= sxe);
93
      if (!*bxe)
94
        {
95
          bx = b->_x;
96
          while (--bxe > bx && !*bxe)
97
            --n;
98
          b->_wds = n;
99
        }
100
    }
101
  if (cmp (b, S) >= 0)
102
    {
103
      q++;
104
      borrow = 0;
105
      carry = 0;
106
      bx = b->_x;
107
      sx = S->_x;
108
      do
109
        {
110
#ifdef Pack_32
111
          si = *sx++;
112
          ys = (si & 0xffff) + carry;
113
          zs = (si >> 16) + (ys >> 16);
114
          carry = zs >> 16;
115
          y = (*bx & 0xffff) - (ys & 0xffff) + borrow;
116
          borrow = y >> 16;
117
          Sign_Extend (borrow, y);
118
          z = (*bx >> 16) - (zs & 0xffff) + borrow;
119
          borrow = z >> 16;
120
          Sign_Extend (borrow, z);
121
          Storeinc (bx, z, y);
122
#else
123
          ys = *sx++ + carry;
124
          carry = ys >> 16;
125
          y = *bx - (ys & 0xffff) + borrow;
126
          borrow = y >> 16;
127
          Sign_Extend (borrow, y);
128
          *bx++ = y & 0xffff;
129
#endif
130
        }
131
      while (sx <= sxe);
132
      bx = b->_x;
133
      bxe = bx + n;
134
      if (!*bxe)
135
        {
136
          while (--bxe > bx && !*bxe)
137
            --n;
138
          b->_wds = n;
139
        }
140
    }
141
  return q;
142
}
143
 
144
/* dtoa for IEEE arithmetic (dmg): convert double to ASCII string.
145
 *
146
 * Inspired by "How to Print Floating-Point Numbers Accurately" by
147
 * Guy L. Steele, Jr. and Jon L. White [Proc. ACM SIGPLAN '90, pp. 92-101].
148
 *
149
 * Modifications:
150
 *      1. Rather than iterating, we use a simple numeric overestimate
151
 *         to determine k = floor(log10(d)).  We scale relevant
152
 *         quantities using O(log2(k)) rather than O(k) multiplications.
153
 *      2. For some modes > 2 (corresponding to ecvt and fcvt), we don't
154
 *         try to generate digits strictly left to right.  Instead, we
155
 *         compute with fewer bits and propagate the carry if necessary
156
 *         when rounding the final digit up.  This is often faster.
157
 *      3. Under the assumption that input will be rounded nearest,
158
 *         mode 0 renders 1e23 as 1e23 rather than 9.999999999999999e22.
159
 *         That is, we allow equality in stopping tests when the
160
 *         round-nearest rule will give the same floating-point value
161
 *         as would satisfaction of the stopping test with strict
162
 *         inequality.
163
 *      4. We remove common factors of powers of 2 from relevant
164
 *         quantities.
165
 *      5. When converting floating-point integers less than 1e16,
166
 *         we use floating-point arithmetic rather than resorting
167
 *         to multiple-precision integers.
168
 *      6. When asked to produce fewer than 15 digits, we first try
169
 *         to get by with floating-point arithmetic; we resort to
170
 *         multiple-precision integer arithmetic only if we cannot
171
 *         guarantee that the floating-point calculation has given
172
 *         the correctly rounded result.  For k requested digits and
173
 *         "uniformly" distributed input, the probability is
174
 *         something like 10^(k-15) that we must resort to the long
175
 *         calculation.
176
 */
177
 
178
 
179
char *
180
_DEFUN (_dtoa_r,
181
        (ptr, _d, mode, ndigits, decpt, sign, rve),
182
        struct _reent *ptr _AND
183
        double _d _AND
184
        int mode _AND
185
        int ndigits _AND
186
        int *decpt _AND
187
        int *sign _AND
188
        char **rve)
189
{
190
  /*    Arguments ndigits, decpt, sign are similar to those
191
        of ecvt and fcvt; trailing zeros are suppressed from
192
        the returned string.  If not null, *rve is set to point
193
        to the end of the return value.  If d is +-Infinity or NaN,
194
        then *decpt is set to 9999.
195
 
196
        mode:
197
 
198
                        and rounded to nearest.
199
                1 ==> like 0, but with Steele & White stopping rule;
200
                        e.g. with IEEE P754 arithmetic , mode 0 gives
201
                        1e23 whereas mode 1 gives 9.999999999999999e22.
202
                2 ==> max(1,ndigits) significant digits.  This gives a
203
                        return value similar to that of ecvt, except
204
                        that trailing zeros are suppressed.
205
                3 ==> through ndigits past the decimal point.  This
206
                        gives a return value similar to that from fcvt,
207
                        except that trailing zeros are suppressed, and
208
                        ndigits can be negative.
209
                4-9 should give the same return values as 2-3, i.e.,
210
                        4 <= mode <= 9 ==> same return as mode
211
                        2 + (mode & 1).  These modes are mainly for
212
                        debugging; often they run slower but sometimes
213
                        faster than modes 2-3.
214
                4,5,8,9 ==> left-to-right digit generation.
215
                6-9 ==> don't try fast floating-point estimate
216
                        (if applicable).
217
 
218
                Values of mode other than 0-9 are treated as mode 0.
219
 
220
                Sufficient space is allocated to the return value
221
                to hold the suppressed trailing zeros.
222
        */
223
 
224
  int bbits, b2, b5, be, dig, i, ieps, ilim, ilim0, ilim1, j, j1, k, k0,
225
    k_check, leftright, m2, m5, s2, s5, spec_case, try_quick;
226
  union double_union d, d2, eps;
227
  __Long L;
228
#ifndef Sudden_Underflow
229
  int denorm;
230
  __ULong x;
231
#endif
232
  _Bigint *b, *b1, *delta, *mlo = NULL, *mhi, *S;
233
  double ds;
234
  char *s, *s0;
235
 
236
  d.d = _d;
237
 
238
  _REENT_CHECK_MP(ptr);
239
  if (_REENT_MP_RESULT(ptr))
240
    {
241
      _REENT_MP_RESULT(ptr)->_k = _REENT_MP_RESULT_K(ptr);
242
      _REENT_MP_RESULT(ptr)->_maxwds = 1 << _REENT_MP_RESULT_K(ptr);
243
      Bfree (ptr, _REENT_MP_RESULT(ptr));
244
      _REENT_MP_RESULT(ptr) = 0;
245
    }
246
 
247
  if (word0 (d) & Sign_bit)
248
    {
249
      /* set sign for everything, including 0's and NaNs */
250
      *sign = 1;
251
      word0 (d) &= ~Sign_bit;   /* clear sign bit */
252
    }
253
  else
254
    *sign = 0;
255
 
256
#if defined(IEEE_Arith) + defined(VAX)
257
#ifdef IEEE_Arith
258
  if ((word0 (d) & Exp_mask) == Exp_mask)
259
#else
260
  if (word0 (d) == 0x8000)
261
#endif
262
    {
263
      /* Infinity or NaN */
264
      *decpt = 9999;
265
      s =
266
#ifdef IEEE_Arith
267
        !word1 (d) && !(word0 (d) & 0xfffff) ? "Infinity" :
268
#endif
269
        "NaN";
270
      if (rve)
271
        *rve =
272
#ifdef IEEE_Arith
273
          s[3] ? s + 8 :
274
#endif
275
          s + 3;
276
      return s;
277
    }
278
#endif
279
#ifdef IBM
280
  d.d += 0;                      /* normalize */
281
#endif
282
  if (!d.d)
283
    {
284
      *decpt = 1;
285
      s = "0";
286
      if (rve)
287
        *rve = s + 1;
288
      return s;
289
    }
290
 
291
  b = d2b (ptr, d.d, &be, &bbits);
292
#ifdef Sudden_Underflow
293
  i = (int) (word0 (d) >> Exp_shift1 & (Exp_mask >> Exp_shift1));
294
#else
295
  if ((i = (int) (word0 (d) >> Exp_shift1 & (Exp_mask >> Exp_shift1))) != 0)
296
    {
297
#endif
298
      d2.d = d.d;
299
      word0 (d2) &= Frac_mask1;
300
      word0 (d2) |= Exp_11;
301
#ifdef IBM
302
      if (j = 11 - hi0bits (word0 (d2) & Frac_mask))
303
        d2.d /= 1 << j;
304
#endif
305
 
306
      /* log(x) ~=~ log(1.5) + (x-1.5)/1.5
307
                 * log10(x)      =  log(x) / log(10)
308
                 *              ~=~ log(1.5)/log(10) + (x-1.5)/(1.5*log(10))
309
                 * log10(d) = (i-Bias)*log(2)/log(10) + log10(d2)
310
                 *
311
                 * This suggests computing an approximation k to log10(d) by
312
                 *
313
                 * k = (i - Bias)*0.301029995663981
314
                 *      + ( (d2-1.5)*0.289529654602168 + 0.176091259055681 );
315
                 *
316
                 * We want k to be too large rather than too small.
317
                 * The error in the first-order Taylor series approximation
318
                 * is in our favor, so we just round up the constant enough
319
                 * to compensate for any error in the multiplication of
320
                 * (i - Bias) by 0.301029995663981; since |i - Bias| <= 1077,
321
                 * and 1077 * 0.30103 * 2^-52 ~=~ 7.2e-14,
322
                 * adding 1e-13 to the constant term more than suffices.
323
                 * Hence we adjust the constant term to 0.1760912590558.
324
                 * (We could get a more accurate k by invoking log10,
325
                 *  but this is probably not worthwhile.)
326
                 */
327
 
328
      i -= Bias;
329
#ifdef IBM
330
      i <<= 2;
331
      i += j;
332
#endif
333
#ifndef Sudden_Underflow
334
      denorm = 0;
335
    }
336
  else
337
    {
338
      /* d is denormalized */
339
 
340
      i = bbits + be + (Bias + (P - 1) - 1);
341
#if defined (_DOUBLE_IS_32BITS)
342
      x = word0 (d) << (32 - i);
343
#else
344
      x = (i > 32) ? (word0 (d) << (64 - i)) | (word1 (d) >> (i - 32))
345
       : (word1 (d) << (32 - i));
346
#endif
347
      d2.d = x;
348
      word0 (d2) -= 31 * Exp_msk1;      /* adjust exponent */
349
      i -= (Bias + (P - 1) - 1) + 1;
350
      denorm = 1;
351
    }
352
#endif
353
#if defined (_DOUBLE_IS_32BITS)
354
  ds = (d2.d - 1.5) * 0.289529651 + 0.176091269 + i * 0.30103001;
355
#else
356
  ds = (d2.d - 1.5) * 0.289529654602168 + 0.1760912590558 + i * 0.301029995663981;
357
#endif
358
  k = (int) ds;
359
  if (ds < 0. && ds != k)
360
    k--;                        /* want k = floor(ds) */
361
  k_check = 1;
362
  if (k >= 0 && k <= Ten_pmax)
363
    {
364
      if (d.d < tens[k])
365
        k--;
366
      k_check = 0;
367
    }
368
  j = bbits - i - 1;
369
  if (j >= 0)
370
    {
371
      b2 = 0;
372
      s2 = j;
373
    }
374
  else
375
    {
376
      b2 = -j;
377
      s2 = 0;
378
    }
379
  if (k >= 0)
380
    {
381
      b5 = 0;
382
      s5 = k;
383
      s2 += k;
384
    }
385
  else
386
    {
387
      b2 -= k;
388
      b5 = -k;
389
      s5 = 0;
390
    }
391
  if (mode < 0 || mode > 9)
392
    mode = 0;
393
  try_quick = 1;
394
  if (mode > 5)
395
    {
396
      mode -= 4;
397
      try_quick = 0;
398
    }
399
  leftright = 1;
400
  ilim = ilim1 = -1;
401
  switch (mode)
402
    {
403
    case 0:
404
    case 1:
405
      i = 18;
406
      ndigits = 0;
407
      break;
408
    case 2:
409
      leftright = 0;
410
      /* no break */
411
    case 4:
412
      if (ndigits <= 0)
413
        ndigits = 1;
414
      ilim = ilim1 = i = ndigits;
415
      break;
416
    case 3:
417
      leftright = 0;
418
      /* no break */
419
    case 5:
420
      i = ndigits + k + 1;
421
      ilim = i;
422
      ilim1 = i - 1;
423
      if (i <= 0)
424
        i = 1;
425
    }
426
  j = sizeof (__ULong);
427
  for (_REENT_MP_RESULT_K(ptr) = 0; sizeof (_Bigint) - sizeof (__ULong) + j <= i;
428
       j <<= 1)
429
    _REENT_MP_RESULT_K(ptr)++;
430
  _REENT_MP_RESULT(ptr) = Balloc (ptr, _REENT_MP_RESULT_K(ptr));
431
  s = s0 = (char *) _REENT_MP_RESULT(ptr);
432
 
433
  if (ilim >= 0 && ilim <= Quick_max && try_quick)
434
    {
435
      /* Try to get by with floating-point arithmetic. */
436
 
437
      i = 0;
438
      d2.d = d.d;
439
      k0 = k;
440
      ilim0 = ilim;
441
      ieps = 2;                 /* conservative */
442
      if (k > 0)
443
        {
444
          ds = tens[k & 0xf];
445
          j = k >> 4;
446
          if (j & Bletch)
447
            {
448
              /* prevent overflows */
449
              j &= Bletch - 1;
450
              d.d /= bigtens[n_bigtens - 1];
451
              ieps++;
452
            }
453
          for (; j; j >>= 1, i++)
454
            if (j & 1)
455
              {
456
                ieps++;
457
                ds *= bigtens[i];
458
              }
459
          d.d /= ds;
460
        }
461
      else if ((j1 = -k) != 0)
462
        {
463
          d.d *= tens[j1 & 0xf];
464
          for (j = j1 >> 4; j; j >>= 1, i++)
465
            if (j & 1)
466
              {
467
                ieps++;
468
                d.d *= bigtens[i];
469
              }
470
        }
471
      if (k_check && d.d < 1. && ilim > 0)
472
        {
473
          if (ilim1 <= 0)
474
            goto fast_failed;
475
          ilim = ilim1;
476
          k--;
477
          d.d *= 10.;
478
          ieps++;
479
        }
480
      eps.d = ieps * d.d + 7.;
481
      word0 (eps) -= (P - 1) * Exp_msk1;
482
      if (ilim == 0)
483
        {
484
          S = mhi = 0;
485
          d.d -= 5.;
486
          if (d.d > eps.d)
487
            goto one_digit;
488
          if (d.d < -eps.d)
489
            goto no_digits;
490
          goto fast_failed;
491
        }
492
#ifndef No_leftright
493
      if (leftright)
494
        {
495
          /* Use Steele & White method of only
496
           * generating digits needed.
497
           */
498
          eps.d = 0.5 / tens[ilim - 1] - eps.d;
499
          for (i = 0;;)
500
            {
501
              L = d.d;
502
              d.d -= L;
503
              *s++ = '0' + (int) L;
504
              if (d.d < eps.d)
505
                goto ret1;
506
              if (1. - d.d < eps.d)
507
                goto bump_up;
508
              if (++i >= ilim)
509
                break;
510
              eps.d *= 10.;
511
              d.d *= 10.;
512
            }
513
        }
514
      else
515
        {
516
#endif
517
          /* Generate ilim digits, then fix them up. */
518
          eps.d *= tens[ilim - 1];
519
          for (i = 1;; i++, d.d *= 10.)
520
            {
521
              L = d.d;
522
              d.d -= L;
523
              *s++ = '0' + (int) L;
524
              if (i == ilim)
525
                {
526
                  if (d.d > 0.5 + eps.d)
527
                    goto bump_up;
528
                  else if (d.d < 0.5 - eps.d)
529
                    {
530
                      while (*--s == '0');
531
                      s++;
532
                      goto ret1;
533
                    }
534
                  break;
535
                }
536
            }
537
#ifndef No_leftright
538
        }
539
#endif
540
    fast_failed:
541
      s = s0;
542
      d.d = d2.d;
543
      k = k0;
544
      ilim = ilim0;
545
    }
546
 
547
  /* Do we have a "small" integer? */
548
 
549
  if (be >= 0 && k <= Int_max)
550
    {
551
      /* Yes. */
552
      ds = tens[k];
553
      if (ndigits < 0 && ilim <= 0)
554
        {
555
          S = mhi = 0;
556
          if (ilim < 0 || d.d <= 5 * ds)
557
            goto no_digits;
558
          goto one_digit;
559
        }
560
      for (i = 1;; i++)
561
        {
562
          L = d.d / ds;
563
          d.d -= L * ds;
564
#ifdef Check_FLT_ROUNDS
565
          /* If FLT_ROUNDS == 2, L will usually be high by 1 */
566
          if (d.d < 0)
567
            {
568
              L--;
569
              d.d += ds;
570
            }
571
#endif
572
          *s++ = '0' + (int) L;
573
          if (i == ilim)
574
            {
575
              d.d += d.d;
576
             if ((d.d > ds) || ((d.d == ds) && (L & 1)))
577
                {
578
                bump_up:
579
                  while (*--s == '9')
580
                    if (s == s0)
581
                      {
582
                        k++;
583
                        *s = '0';
584
                        break;
585
                      }
586
                  ++*s++;
587
                }
588
              break;
589
            }
590
          if (!(d.d *= 10.))
591
            break;
592
        }
593
      goto ret1;
594
    }
595
 
596
  m2 = b2;
597
  m5 = b5;
598
  mhi = mlo = 0;
599
  if (leftright)
600
    {
601
      if (mode < 2)
602
        {
603
          i =
604
#ifndef Sudden_Underflow
605
            denorm ? be + (Bias + (P - 1) - 1 + 1) :
606
#endif
607
#ifdef IBM
608
            1 + 4 * P - 3 - bbits + ((bbits + be - 1) & 3);
609
#else
610
            1 + P - bbits;
611
#endif
612
        }
613
      else
614
        {
615
          j = ilim - 1;
616
          if (m5 >= j)
617
            m5 -= j;
618
          else
619
            {
620
              s5 += j -= m5;
621
              b5 += j;
622
              m5 = 0;
623
            }
624
          if ((i = ilim) < 0)
625
            {
626
              m2 -= i;
627
              i = 0;
628
            }
629
        }
630
      b2 += i;
631
      s2 += i;
632
      mhi = i2b (ptr, 1);
633
    }
634
  if (m2 > 0 && s2 > 0)
635
    {
636
      i = m2 < s2 ? m2 : s2;
637
      b2 -= i;
638
      m2 -= i;
639
      s2 -= i;
640
    }
641
  if (b5 > 0)
642
    {
643
      if (leftright)
644
        {
645
          if (m5 > 0)
646
            {
647
              mhi = pow5mult (ptr, mhi, m5);
648
              b1 = mult (ptr, mhi, b);
649
              Bfree (ptr, b);
650
              b = b1;
651
            }
652
         if ((j = b5 - m5) != 0)
653
            b = pow5mult (ptr, b, j);
654
        }
655
      else
656
        b = pow5mult (ptr, b, b5);
657
    }
658
  S = i2b (ptr, 1);
659
  if (s5 > 0)
660
    S = pow5mult (ptr, S, s5);
661
 
662
  /* Check for special case that d is a normalized power of 2. */
663
 
664
  spec_case = 0;
665
  if (mode < 2)
666
    {
667
      if (!word1 (d) && !(word0 (d) & Bndry_mask)
668
#ifndef Sudden_Underflow
669
          && word0 (d) & Exp_mask
670
#endif
671
        )
672
        {
673
          /* The special case */
674
          b2 += Log2P;
675
          s2 += Log2P;
676
          spec_case = 1;
677
        }
678
    }
679
 
680
  /* Arrange for convenient computation of quotients:
681
   * shift left if necessary so divisor has 4 leading 0 bits.
682
   *
683
   * Perhaps we should just compute leading 28 bits of S once
684
   * and for all and pass them and a shift to quorem, so it
685
   * can do shifts and ors to compute the numerator for q.
686
   */
687
 
688
#ifdef Pack_32
689
  if ((i = ((s5 ? 32 - hi0bits (S->_x[S->_wds - 1]) : 1) + s2) & 0x1f) != 0)
690
    i = 32 - i;
691
#else
692
  if ((i = ((s5 ? 32 - hi0bits (S->_x[S->_wds - 1]) : 1) + s2) & 0xf) != 0)
693
    i = 16 - i;
694
#endif
695
  if (i > 4)
696
    {
697
      i -= 4;
698
      b2 += i;
699
      m2 += i;
700
      s2 += i;
701
    }
702
  else if (i < 4)
703
    {
704
      i += 28;
705
      b2 += i;
706
      m2 += i;
707
      s2 += i;
708
    }
709
  if (b2 > 0)
710
    b = lshift (ptr, b, b2);
711
  if (s2 > 0)
712
    S = lshift (ptr, S, s2);
713
  if (k_check)
714
    {
715
      if (cmp (b, S) < 0)
716
        {
717
          k--;
718
          b = multadd (ptr, b, 10, 0);   /* we botched the k estimate */
719
          if (leftright)
720
            mhi = multadd (ptr, mhi, 10, 0);
721
          ilim = ilim1;
722
        }
723
    }
724
  if (ilim <= 0 && mode > 2)
725
    {
726
      if (ilim < 0 || cmp (b, S = multadd (ptr, S, 5, 0)) <= 0)
727
        {
728
          /* no digits, fcvt style */
729
        no_digits:
730
          k = -1 - ndigits;
731
          goto ret;
732
        }
733
    one_digit:
734
      *s++ = '1';
735
      k++;
736
      goto ret;
737
    }
738
  if (leftright)
739
    {
740
      if (m2 > 0)
741
        mhi = lshift (ptr, mhi, m2);
742
 
743
      /* Compute mlo -- check for special case
744
       * that d is a normalized power of 2.
745
       */
746
 
747
      mlo = mhi;
748
      if (spec_case)
749
        {
750
          mhi = Balloc (ptr, mhi->_k);
751
          Bcopy (mhi, mlo);
752
          mhi = lshift (ptr, mhi, Log2P);
753
        }
754
 
755
      for (i = 1;; i++)
756
        {
757
          dig = quorem (b, S) + '0';
758
          /* Do we yet have the shortest decimal string
759
           * that will round to d?
760
           */
761
          j = cmp (b, mlo);
762
          delta = diff (ptr, S, mhi);
763
          j1 = delta->_sign ? 1 : cmp (b, delta);
764
          Bfree (ptr, delta);
765
#ifndef ROUND_BIASED
766
          if (j1 == 0 && !mode && !(word1 (d) & 1))
767
            {
768
              if (dig == '9')
769
                goto round_9_up;
770
              if (j > 0)
771
                dig++;
772
              *s++ = dig;
773
              goto ret;
774
            }
775
#endif
776
         if ((j < 0) || ((j == 0) && !mode
777
#ifndef ROUND_BIASED
778
              && !(word1 (d) & 1)
779
#endif
780
           ))
781
            {
782
              if (j1 > 0)
783
                {
784
                  b = lshift (ptr, b, 1);
785
                  j1 = cmp (b, S);
786
                 if (((j1 > 0) || ((j1 == 0) && (dig & 1)))
787
                      && dig++ == '9')
788
                    goto round_9_up;
789
                }
790
              *s++ = dig;
791
              goto ret;
792
            }
793
          if (j1 > 0)
794
            {
795
              if (dig == '9')
796
                {               /* possible if i == 1 */
797
                round_9_up:
798
                  *s++ = '9';
799
                  goto roundoff;
800
                }
801
              *s++ = dig + 1;
802
              goto ret;
803
            }
804
          *s++ = dig;
805
          if (i == ilim)
806
            break;
807
          b = multadd (ptr, b, 10, 0);
808
          if (mlo == mhi)
809
            mlo = mhi = multadd (ptr, mhi, 10, 0);
810
          else
811
            {
812
              mlo = multadd (ptr, mlo, 10, 0);
813
              mhi = multadd (ptr, mhi, 10, 0);
814
            }
815
        }
816
    }
817
  else
818
    for (i = 1;; i++)
819
      {
820
        *s++ = dig = quorem (b, S) + '0';
821
        if (i >= ilim)
822
          break;
823
        b = multadd (ptr, b, 10, 0);
824
      }
825
 
826
  /* Round off last digit */
827
 
828
  b = lshift (ptr, b, 1);
829
  j = cmp (b, S);
830
  if ((j > 0) || ((j == 0) && (dig & 1)))
831
    {
832
    roundoff:
833
      while (*--s == '9')
834
        if (s == s0)
835
          {
836
            k++;
837
            *s++ = '1';
838
            goto ret;
839
          }
840
      ++*s++;
841
    }
842
  else
843
    {
844
      while (*--s == '0');
845
      s++;
846
    }
847
ret:
848
  Bfree (ptr, S);
849
  if (mhi)
850
    {
851
      if (mlo && mlo != mhi)
852
        Bfree (ptr, mlo);
853
      Bfree (ptr, mhi);
854
    }
855
ret1:
856
  Bfree (ptr, b);
857
  *s = 0;
858
  *decpt = k + 1;
859
  if (rve)
860
    *rve = s;
861
  return s0;
862
}

powered by: WebSVN 2.1.0

© copyright 1999-2024 OpenCores.org, equivalent to Oliscience, all rights reserved. OpenCores®, registered trademark.