1 /* Floating point output for `printf'.
2 Copyright (C) 1995-1999, 2000, 2001 Free Software Foundation, Inc.
3 This file is part of the GNU C Library.
4 Written by Ulrich Drepper <drepper@gnu.ai.mit.edu>, 1995.
5
6 The GNU C Library is free software; you can redistribute it and/or
7 modify it under the terms of the GNU Lesser General Public
8 License as published by the Free Software Foundation; either
9 version 2.1 of the License, or (at your option) any later version.
10
11 The GNU C Library is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 Lesser General Public License for more details.
15
16 You should have received a copy of the GNU Lesser General Public
17 License along with the GNU C Library; if not, write to the Free
18 Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
19 02111-1307 USA. */
20
21 /* The gmp headers need some configuration frobs. */
22 #define HAVE_ALLOCA 1
23
24 #ifdef USE_IN_LIBIO
25 # include <libioP.h>
26 #else
27 # include <stdio.h>
28 #endif
29 #include <alloca.h>
30 #include <ctype.h>
31 #include <float.h>
32 #include <gmp-mparam.h>
33 #include <gmp.h>
34 #include <stdlib/gmp-impl.h>
35 #include <stdlib/longlong.h>
36 #include <stdlib/fpioconst.h>
37 #include <locale/localeinfo.h>
38 #include <limits.h>
39 #include <math.h>
40 #include <printf.h>
41 #include <string.h>
42 #include <unistd.h>
43 #include <stdlib.h>
44 #include <wchar.h>
45
46 #ifndef NDEBUG
47 # define NDEBUG /* Undefine this for debugging assertions. */
48 #endif
49 #include <assert.h>
50
51 /* This defines make it possible to use the same code for GNU C library and
52 the GNU I/O library. */
53 #ifdef USE_IN_LIBIO
54 # define PUT(f, s, n) _IO_sputn (f, s, n)
55 # define PAD(f, c, n) (wide ? _IO_wpadn (f, c, n) : _IO_padn (f, c, n))
56 /* We use this file GNU C library and GNU I/O library. So make
57 names equal. */
58 # undef putc
59 # define putc(c, f) (wide \
60 ? (int)_IO_putwc_unlocked (c, f) : _IO_putc_unlocked (c, f))
61 # define size_t _IO_size_t
62 # define FILE _IO_FILE
63 #else /* ! USE_IN_LIBIO */
64 # define PUT(f, s, n) fwrite (s, 1, n, f)
65 # define PAD(f, c, n) __printf_pad (f, c, n)
66 ssize_t __printf_pad __P ((FILE *, char pad, int n)); /* In vfprintf.c. */
67 #endif /* USE_IN_LIBIO */
68
69 /* Macros for doing the actual output. */
70
71 #define outchar(ch) \
72 do \
73 { \
74 register const int outc = (ch); \
75 if (putc (outc, fp) == EOF) \
76 return -1; \
77 ++done; \
78 } while (0)
79
80 #define PRINT(ptr, wptr, len) \
81 do \
82 { \
83 register size_t outlen = (len); \
84 if (len > 20) \
85 { \
86 if (PUT (fp, wide ? (const char *) wptr : ptr, outlen) != outlen) \
87 return -1; \
88 ptr += outlen; \
89 done += outlen; \
90 } \
91 else \
92 { \
93 if (wide) \
94 while (outlen-- > 0) \
95 outchar (*wptr++); \
96 else \
97 while (outlen-- > 0) \
98 outchar (*ptr++); \
99 } \
100 } while (0)
101
102 #define PADN(ch, len) \
103 do \
104 { \
105 if (PAD (fp, ch, len) != len) \
106 return -1; \
107 done += len; \
108 } \
109 while (0)
110
111 /* We use the GNU MP library to handle large numbers.
112
113 An MP variable occupies a varying number of entries in its array. We keep
114 track of this number for efficiency reasons. Otherwise we would always
115 have to process the whole array. */
116 #define MPN_VAR(name) mp_limb_t *name; mp_size_t name##size
117
118 #define MPN_ASSIGN(dst,src) \
119 memcpy (dst, src, (dst##size = src##size) * sizeof (mp_limb_t))
120 #define MPN_GE(u,v) \
121 (u##size > v##size || (u##size == v##size && __mpn_cmp (u, v, u##size) >= 0))
122
123 extern mp_size_t __mpn_extract_double (mp_ptr res_ptr, mp_size_t size,
124 int *expt, int *is_neg,
125 double value);
126 extern mp_size_t __mpn_extract_long_double (mp_ptr res_ptr, mp_size_t size,
127 int *expt, int *is_neg,
128 long double value);
129 extern unsigned int __guess_grouping (unsigned int intdig_max,
130 const char *grouping);
131
132
133 static wchar_t *group_number (wchar_t *buf, wchar_t *bufend,
134 unsigned int intdig_no, const char *grouping,
135 wchar_t thousands_sep, int ngroups)
136 internal_function;
137
138
139 int
__printf_fp(FILE * fp,const struct printf_info * info,const void * const * args)140 __printf_fp (FILE *fp,
141 const struct printf_info *info,
142 const void *const *args)
143 {
144 /* The floating-point value to output. */
145 union
146 {
147 double dbl;
148 __long_double_t ldbl;
149 }
150 fpnum;
151
152 /* Locale-dependent representation of decimal point. */
153 const char *decimal;
154 wchar_t decimalwc;
155
156 /* Locale-dependent thousands separator and grouping specification. */
157 const char *thousands_sep = NULL;
158 wchar_t thousands_sepwc = 0;
159 const char *grouping;
160
161 /* "NaN" or "Inf" for the special cases. */
162 const char *special = NULL;
163 const wchar_t *wspecial = NULL;
164
165 /* We need just a few limbs for the input before shifting to the right
166 position. */
167 mp_limb_t fp_input[(LDBL_MANT_DIG + BITS_PER_MP_LIMB - 1) / BITS_PER_MP_LIMB];
168 /* We need to shift the contents of fp_input by this amount of bits. */
169 int to_shift = 0;
170
171 /* The fraction of the floting-point value in question */
172 MPN_VAR(frac);
173 /* and the exponent. */
174 int exponent;
175 /* Sign of the exponent. */
176 int expsign = 0;
177 /* Sign of float number. */
178 int is_neg = 0;
179
180 /* Scaling factor. */
181 MPN_VAR(scale);
182
183 /* Temporary bignum value. */
184 MPN_VAR(tmp);
185
186 /* Digit which is result of last hack_digit() call. */
187 wchar_t digit;
188
189 /* The type of output format that will be used: 'e'/'E' or 'f'/'F'. */
190 int type;
191
192 /* Counter for number of written characters. */
193 int done = 0;
194
195 /* General helper (carry limb). */
196 mp_limb_t cy;
197
198 /* Nonzero if this is output on a wide character stream. */
199 int wide = info->wide;
200
201 wchar_t hack_digit_ret;
202 int hack_digit_callee;
203
204 while (0)
205 {
206 mp_limb_t hi;
207
208 hack_digit:
209 if (expsign != 0 && _tolower (type) == 'f' && exponent-- > 0)
210 hi = 0;
211 else if (scalesize == 0)
212 {
213 hi = frac[fracsize - 1];
214 cy = __mpn_mul_1 (frac, frac, fracsize - 1, 10);
215 frac[fracsize - 1] = cy;
216 }
217 else
218 {
219 if (fracsize < scalesize)
220 hi = 0;
221 else
222 {
223 hi = mpn_divmod (tmp, frac, fracsize, scale, scalesize);
224 tmp[fracsize - scalesize] = hi;
225 hi = tmp[0];
226
227 fracsize = scalesize;
228 while (fracsize != 0 && frac[fracsize - 1] == 0)
229 --fracsize;
230 if (fracsize == 0)
231 {
232 /* We're not prepared for an mpn variable with zero
233 limbs. */
234 fracsize = 1;
235 hack_digit_ret = L'0' + hi;
236 goto hack_digit_end;
237 }
238 }
239
240 cy = __mpn_mul_1 (frac, frac, fracsize, 10);
241 if (cy != 0)
242 frac[fracsize++] = cy;
243 }
244
245 hack_digit_ret = L'0' + hi;
246 hack_digit_end:
247 switch (hack_digit_callee)
248 {
249 case 1: goto hack_digit_callee1;
250 case 2: goto hack_digit_callee2;
251 case 3: goto hack_digit_callee3;
252 default: abort();
253 }
254 }
255
256
257 /* Figure out the decimal point character. */
258 if (info->extra == 0)
259 {
260 decimal = _NL_CURRENT (LC_NUMERIC, DECIMAL_POINT);
261 decimalwc = _NL_CURRENT_WORD (LC_NUMERIC, _NL_NUMERIC_DECIMAL_POINT_WC);
262 }
263 else
264 {
265 decimal = _NL_CURRENT (LC_MONETARY, MON_DECIMAL_POINT);
266 if (*decimal == '\0')
267 decimal = _NL_CURRENT (LC_NUMERIC, DECIMAL_POINT);
268 decimalwc = _NL_CURRENT_WORD (LC_MONETARY,
269 _NL_MONETARY_DECIMAL_POINT_WC);
270 if (decimalwc == L'\0')
271 decimalwc = _NL_CURRENT_WORD (LC_NUMERIC,
272 _NL_NUMERIC_DECIMAL_POINT_WC);
273 }
274 /* The decimal point character must not be zero. */
275 assert (*decimal != '\0');
276 assert (decimalwc != L'\0');
277
278 if (info->group)
279 {
280 if (info->extra == 0)
281 grouping = _NL_CURRENT (LC_NUMERIC, GROUPING);
282 else
283 grouping = _NL_CURRENT (LC_MONETARY, MON_GROUPING);
284
285 if (*grouping <= 0 || *grouping == CHAR_MAX)
286 grouping = NULL;
287 else
288 {
289 /* Figure out the thousands separator character. */
290 if (wide)
291 {
292 if (info->extra == 0)
293 thousands_sepwc =
294 _NL_CURRENT_WORD (LC_NUMERIC, _NL_NUMERIC_THOUSANDS_SEP_WC);
295 else
296 thousands_sepwc =
297 _NL_CURRENT_WORD (LC_MONETARY,
298 _NL_MONETARY_THOUSANDS_SEP_WC);
299 }
300 else
301 {
302 if (info->extra == 0)
303 thousands_sep = _NL_CURRENT (LC_NUMERIC, THOUSANDS_SEP);
304 else
305 thousands_sep = _NL_CURRENT (LC_MONETARY, MON_THOUSANDS_SEP);
306 }
307
308 if ((wide && thousands_sepwc == L'\0')
309 || (! wide && *thousands_sep == '\0'))
310 grouping = NULL;
311 else if (thousands_sepwc == L'\0')
312 /* If we are printing multibyte characters and there is a
313 multibyte representation for the thousands separator,
314 we must ensure the wide character thousands separator
315 is available, even if it is fake. */
316 thousands_sepwc = 0xfffffffe;
317 }
318 }
319 else
320 grouping = NULL;
321
322 /* Fetch the argument value. */
323 #ifndef __NO_LONG_DOUBLE_MATH
324 if (info->is_long_double && sizeof (long double) > sizeof (double))
325 {
326 fpnum.ldbl = *(const long double *) args[0];
327
328 /* Check for special values: not a number or infinity. */
329 if (isnan (fpnum.ldbl))
330 {
331 if (isupper (info->spec))
332 {
333 special = "NAN";
334 wspecial = L"NAN";
335 }
336 else
337 {
338 special = "nan";
339 wspecial = L"nan";
340 }
341 is_neg = 0;
342 }
343 else if (isinf (fpnum.ldbl))
344 {
345 if (isupper (info->spec))
346 {
347 special = "INF";
348 wspecial = L"INF";
349 }
350 else
351 {
352 special = "inf";
353 wspecial = L"inf";
354 }
355 is_neg = fpnum.ldbl < 0;
356 }
357 else
358 {
359 fracsize = __mpn_extract_long_double (fp_input,
360 (sizeof (fp_input) /
361 sizeof (fp_input[0])),
362 &exponent, &is_neg,
363 fpnum.ldbl);
364 to_shift = 1 + fracsize * BITS_PER_MP_LIMB - LDBL_MANT_DIG;
365 }
366 }
367 else
368 #endif /* no long double */
369 {
370 fpnum.dbl = *(const double *) args[0];
371
372 /* Check for special values: not a number or infinity. */
373 if (isnan (fpnum.dbl))
374 {
375 if (isupper (info->spec))
376 {
377 special = "NAN";
378 wspecial = L"NAN";
379 }
380 else
381 {
382 special = "nan";
383 wspecial = L"nan";
384 }
385 is_neg = 0;
386 }
387 else if (isinf (fpnum.dbl))
388 {
389 if (isupper (info->spec))
390 {
391 special = "INF";
392 wspecial = L"INF";
393 }
394 else
395 {
396 special = "inf";
397 wspecial = L"inf";
398 }
399 is_neg = fpnum.dbl < 0;
400 }
401 else
402 {
403 fracsize = __mpn_extract_double (fp_input,
404 (sizeof (fp_input)
405 / sizeof (fp_input[0])),
406 &exponent, &is_neg, fpnum.dbl);
407 to_shift = 1 + fracsize * BITS_PER_MP_LIMB - DBL_MANT_DIG;
408 }
409 }
410
411 if (special)
412 {
413 int width = info->width;
414
415 if (is_neg || info->showsign || info->space)
416 --width;
417 width -= 3;
418
419 if (!info->left && width > 0)
420 PADN (' ', width);
421
422 if (is_neg)
423 outchar ('-');
424 else if (info->showsign)
425 outchar ('+');
426 else if (info->space)
427 outchar (' ');
428
429 PRINT (special, wspecial, 3);
430
431 if (info->left && width > 0)
432 PADN (' ', width);
433
434 return done;
435 }
436
437
438 /* We need three multiprecision variables. Now that we have the exponent
439 of the number we can allocate the needed memory. It would be more
440 efficient to use variables of the fixed maximum size but because this
441 would be really big it could lead to memory problems. */
442 {
443 mp_size_t bignum_size = ((ABS (exponent) + BITS_PER_MP_LIMB - 1)
444 / BITS_PER_MP_LIMB + 4) * sizeof (mp_limb_t);
445 frac = (mp_limb_t *) alloca (bignum_size);
446 tmp = (mp_limb_t *) alloca (bignum_size);
447 scale = (mp_limb_t *) alloca (bignum_size);
448 }
449
450 /* We now have to distinguish between numbers with positive and negative
451 exponents because the method used for the one is not applicable/efficient
452 for the other. */
453 scalesize = 0;
454 if (exponent > 2)
455 {
456 /* |FP| >= 8.0. */
457 int scaleexpo = 0;
458 int explog = LDBL_MAX_10_EXP_LOG;
459 int exp10 = 0;
460 const struct mp_power *powers = &_fpioconst_pow10[explog + 1];
461 int cnt_h, cnt_l, i;
462
463 if ((exponent + to_shift) % BITS_PER_MP_LIMB == 0)
464 {
465 MPN_COPY_DECR (frac + (exponent + to_shift) / BITS_PER_MP_LIMB,
466 fp_input, fracsize);
467 fracsize += (exponent + to_shift) / BITS_PER_MP_LIMB;
468 }
469 else
470 {
471 cy = __mpn_lshift (frac + (exponent + to_shift) / BITS_PER_MP_LIMB,
472 fp_input, fracsize,
473 (exponent + to_shift) % BITS_PER_MP_LIMB);
474 fracsize += (exponent + to_shift) / BITS_PER_MP_LIMB;
475 if (cy)
476 frac[fracsize++] = cy;
477 }
478 MPN_ZERO (frac, (exponent + to_shift) / BITS_PER_MP_LIMB);
479
480 assert (powers > &_fpioconst_pow10[0]);
481 do
482 {
483 --powers;
484
485 /* The number of the product of two binary numbers with n and m
486 bits respectively has m+n or m+n-1 bits. */
487 if (exponent >= scaleexpo + powers->p_expo - 1)
488 {
489 if (scalesize == 0)
490 {
491 #ifndef __NO_LONG_DOUBLE_MATH
492 if (LDBL_MANT_DIG > _FPIO_CONST_OFFSET * BITS_PER_MP_LIMB
493 && info->is_long_double)
494 {
495 #define _FPIO_CONST_SHIFT \
496 (((LDBL_MANT_DIG + BITS_PER_MP_LIMB - 1) / BITS_PER_MP_LIMB) \
497 - _FPIO_CONST_OFFSET)
498 /* 64bit const offset is not enough for
499 IEEE quad long double. */
500 tmpsize = powers->arraysize + _FPIO_CONST_SHIFT;
501 memcpy (tmp + _FPIO_CONST_SHIFT,
502 &__tens[powers->arrayoff],
503 tmpsize * sizeof (mp_limb_t));
504 MPN_ZERO (tmp, _FPIO_CONST_SHIFT);
505 }
506 else
507 #endif
508 {
509 tmpsize = powers->arraysize;
510 memcpy (tmp, &__tens[powers->arrayoff],
511 tmpsize * sizeof (mp_limb_t));
512 }
513 }
514 else
515 {
516 cy = __mpn_mul (tmp, scale, scalesize,
517 &__tens[powers->arrayoff
518 + _FPIO_CONST_OFFSET],
519 powers->arraysize - _FPIO_CONST_OFFSET);
520 tmpsize = scalesize + powers->arraysize - _FPIO_CONST_OFFSET;
521 if (cy == 0)
522 --tmpsize;
523 }
524
525 if (MPN_GE (frac, tmp))
526 {
527 int cnt;
528 MPN_ASSIGN (scale, tmp);
529 count_leading_zeros (cnt, scale[scalesize - 1]);
530 scaleexpo = (scalesize - 2) * BITS_PER_MP_LIMB - cnt - 1;
531 exp10 |= 1 << explog;
532 }
533 }
534 --explog;
535 }
536 while (powers > &_fpioconst_pow10[0]);
537 exponent = exp10;
538
539 /* Optimize number representations. We want to represent the numbers
540 with the lowest number of bytes possible without losing any
541 bytes. Also the highest bit in the scaling factor has to be set
542 (this is a requirement of the MPN division routines). */
543 if (scalesize > 0)
544 {
545 /* Determine minimum number of zero bits at the end of
546 both numbers. */
547 for (i = 0; scale[i] == 0 && frac[i] == 0; i++)
548 ;
549
550 /* Determine number of bits the scaling factor is misplaced. */
551 count_leading_zeros (cnt_h, scale[scalesize - 1]);
552
553 if (cnt_h == 0)
554 {
555 /* The highest bit of the scaling factor is already set. So
556 we only have to remove the trailing empty limbs. */
557 if (i > 0)
558 {
559 MPN_COPY_INCR (scale, scale + i, scalesize - i);
560 scalesize -= i;
561 MPN_COPY_INCR (frac, frac + i, fracsize - i);
562 fracsize -= i;
563 }
564 }
565 else
566 {
567 if (scale[i] != 0)
568 {
569 count_trailing_zeros (cnt_l, scale[i]);
570 if (frac[i] != 0)
571 {
572 int cnt_l2;
573 count_trailing_zeros (cnt_l2, frac[i]);
574 if (cnt_l2 < cnt_l)
575 cnt_l = cnt_l2;
576 }
577 }
578 else
579 count_trailing_zeros (cnt_l, frac[i]);
580
581 /* Now shift the numbers to their optimal position. */
582 if (i == 0 && BITS_PER_MP_LIMB - cnt_h > cnt_l)
583 {
584 /* We cannot save any memory. So just roll both numbers
585 so that the scaling factor has its highest bit set. */
586
587 (void) __mpn_lshift (scale, scale, scalesize, cnt_h);
588 cy = __mpn_lshift (frac, frac, fracsize, cnt_h);
589 if (cy != 0)
590 frac[fracsize++] = cy;
591 }
592 else if (BITS_PER_MP_LIMB - cnt_h <= cnt_l)
593 {
594 /* We can save memory by removing the trailing zero limbs
595 and by packing the non-zero limbs which gain another
596 free one. */
597
598 (void) __mpn_rshift (scale, scale + i, scalesize - i,
599 BITS_PER_MP_LIMB - cnt_h);
600 scalesize -= i + 1;
601 (void) __mpn_rshift (frac, frac + i, fracsize - i,
602 BITS_PER_MP_LIMB - cnt_h);
603 fracsize -= frac[fracsize - i - 1] == 0 ? i + 1 : i;
604 }
605 else
606 {
607 /* We can only save the memory of the limbs which are zero.
608 The non-zero parts occupy the same number of limbs. */
609
610 (void) __mpn_rshift (scale, scale + (i - 1),
611 scalesize - (i - 1),
612 BITS_PER_MP_LIMB - cnt_h);
613 scalesize -= i;
614 (void) __mpn_rshift (frac, frac + (i - 1),
615 fracsize - (i - 1),
616 BITS_PER_MP_LIMB - cnt_h);
617 fracsize -= frac[fracsize - (i - 1) - 1] == 0 ? i : i - 1;
618 }
619 }
620 }
621 }
622 else if (exponent < 0)
623 {
624 /* |FP| < 1.0. */
625 int exp10 = 0;
626 int explog = LDBL_MAX_10_EXP_LOG;
627 const struct mp_power *powers = &_fpioconst_pow10[explog + 1];
628 mp_size_t used_limbs = fracsize - 1;
629
630 /* Now shift the input value to its right place. */
631 cy = __mpn_lshift (frac, fp_input, fracsize, to_shift);
632 frac[fracsize++] = cy;
633 assert (cy == 1 || (frac[fracsize - 2] == 0 && frac[0] == 0));
634
635 expsign = 1;
636 exponent = -exponent;
637
638 assert (powers != &_fpioconst_pow10[0]);
639 do
640 {
641 --powers;
642
643 if (exponent >= powers->m_expo)
644 {
645 int i, incr, cnt_h, cnt_l;
646 mp_limb_t topval[2];
647
648 /* The __mpn_mul function expects the first argument to be
649 bigger than the second. */
650 if (fracsize < powers->arraysize - _FPIO_CONST_OFFSET)
651 cy = __mpn_mul (tmp, &__tens[powers->arrayoff
652 + _FPIO_CONST_OFFSET],
653 powers->arraysize - _FPIO_CONST_OFFSET,
654 frac, fracsize);
655 else
656 cy = __mpn_mul (tmp, frac, fracsize,
657 &__tens[powers->arrayoff + _FPIO_CONST_OFFSET],
658 powers->arraysize - _FPIO_CONST_OFFSET);
659 tmpsize = fracsize + powers->arraysize - _FPIO_CONST_OFFSET;
660 if (cy == 0)
661 --tmpsize;
662
663 count_leading_zeros (cnt_h, tmp[tmpsize - 1]);
664 incr = (tmpsize - fracsize) * BITS_PER_MP_LIMB
665 + BITS_PER_MP_LIMB - 1 - cnt_h;
666
667 assert (incr <= powers->p_expo);
668
669 /* If we increased the exponent by exactly 3 we have to test
670 for overflow. This is done by comparing with 10 shifted
671 to the right position. */
672 if (incr == exponent + 3)
673 {
674 if (cnt_h <= BITS_PER_MP_LIMB - 4)
675 {
676 topval[0] = 0;
677 topval[1]
678 = ((mp_limb_t) 10) << (BITS_PER_MP_LIMB - 4 - cnt_h);
679 }
680 else
681 {
682 topval[0] = ((mp_limb_t) 10) << (BITS_PER_MP_LIMB - 4);
683 topval[1] = 0;
684 (void) __mpn_lshift (topval, topval, 2,
685 BITS_PER_MP_LIMB - cnt_h);
686 }
687 }
688
689 /* We have to be careful when multiplying the last factor.
690 If the result is greater than 1.0 be have to test it
691 against 10.0. If it is greater or equal to 10.0 the
692 multiplication was not valid. This is because we cannot
693 determine the number of bits in the result in advance. */
694 if (incr < exponent + 3
695 || (incr == exponent + 3 &&
696 (tmp[tmpsize - 1] < topval[1]
697 || (tmp[tmpsize - 1] == topval[1]
698 && tmp[tmpsize - 2] < topval[0]))))
699 {
700 /* The factor is right. Adapt binary and decimal
701 exponents. */
702 exponent -= incr;
703 exp10 |= 1 << explog;
704
705 /* If this factor yields a number greater or equal to
706 1.0, we must not shift the non-fractional digits down. */
707 if (exponent < 0)
708 cnt_h += -exponent;
709
710 /* Now we optimize the number representation. */
711 for (i = 0; tmp[i] == 0; ++i);
712 if (cnt_h == BITS_PER_MP_LIMB - 1)
713 {
714 MPN_COPY (frac, tmp + i, tmpsize - i);
715 fracsize = tmpsize - i;
716 }
717 else
718 {
719 count_trailing_zeros (cnt_l, tmp[i]);
720
721 /* Now shift the numbers to their optimal position. */
722 if (i == 0 && BITS_PER_MP_LIMB - 1 - cnt_h > cnt_l)
723 {
724 /* We cannot save any memory. Just roll the
725 number so that the leading digit is in a
726 separate limb. */
727
728 cy = __mpn_lshift (frac, tmp, tmpsize, cnt_h + 1);
729 fracsize = tmpsize + 1;
730 frac[fracsize - 1] = cy;
731 }
732 else if (BITS_PER_MP_LIMB - 1 - cnt_h <= cnt_l)
733 {
734 (void) __mpn_rshift (frac, tmp + i, tmpsize - i,
735 BITS_PER_MP_LIMB - 1 - cnt_h);
736 fracsize = tmpsize - i;
737 }
738 else
739 {
740 /* We can only save the memory of the limbs which
741 are zero. The non-zero parts occupy the same
742 number of limbs. */
743
744 (void) __mpn_rshift (frac, tmp + (i - 1),
745 tmpsize - (i - 1),
746 BITS_PER_MP_LIMB - 1 - cnt_h);
747 fracsize = tmpsize - (i - 1);
748 }
749 }
750 used_limbs = fracsize - 1;
751 }
752 }
753 --explog;
754 }
755 while (powers != &_fpioconst_pow10[1] && exponent > 0);
756 /* All factors but 10^-1 are tested now. */
757 if (exponent > 0)
758 {
759 int cnt_l;
760
761 cy = __mpn_mul_1 (tmp, frac, fracsize, 10);
762 tmpsize = fracsize;
763 assert (cy == 0 || tmp[tmpsize - 1] < 20);
764
765 count_trailing_zeros (cnt_l, tmp[0]);
766 if (cnt_l < MIN (4, exponent))
767 {
768 cy = __mpn_lshift (frac, tmp, tmpsize,
769 BITS_PER_MP_LIMB - MIN (4, exponent));
770 if (cy != 0)
771 frac[tmpsize++] = cy;
772 }
773 else
774 (void) __mpn_rshift (frac, tmp, tmpsize, MIN (4, exponent));
775 fracsize = tmpsize;
776 exp10 |= 1;
777 assert (frac[fracsize - 1] < 10);
778 }
779 exponent = exp10;
780 }
781 else
782 {
783 /* This is a special case. We don't need a factor because the
784 numbers are in the range of 0.0 <= fp < 8.0. We simply
785 shift it to the right place and divide it by 1.0 to get the
786 leading digit. (Of course this division is not really made.) */
787 assert (0 <= exponent && exponent < 3 &&
788 exponent + to_shift < BITS_PER_MP_LIMB);
789
790 /* Now shift the input value to its right place. */
791 cy = __mpn_lshift (frac, fp_input, fracsize, (exponent + to_shift));
792 frac[fracsize++] = cy;
793 exponent = 0;
794 }
795
796 {
797 int width = info->width;
798 wchar_t *wbuffer, *wstartp, *wcp;
799 int buffer_malloced;
800 int chars_needed;
801 int expscale;
802 int intdig_max, intdig_no = 0;
803 int fracdig_min, fracdig_max, fracdig_no = 0;
804 int dig_max;
805 int significant;
806 int ngroups = 0;
807
808 if (_tolower (info->spec) == 'e')
809 {
810 type = info->spec;
811 intdig_max = 1;
812 fracdig_min = fracdig_max = info->prec < 0 ? 6 : info->prec;
813 chars_needed = 1 + 1 + fracdig_max + 1 + 1 + 4;
814 /* d . ddd e +- ddd */
815 dig_max = INT_MAX; /* Unlimited. */
816 significant = 1; /* Does not matter here. */
817 }
818 else if (_tolower (info->spec) == 'f')
819 {
820 type = info->spec;
821 fracdig_min = fracdig_max = info->prec < 0 ? 6 : info->prec;
822 if (expsign == 0)
823 {
824 intdig_max = exponent + 1;
825 /* This can be really big! */ /* XXX Maybe malloc if too big? */
826 chars_needed = exponent + 1 + 1 + fracdig_max;
827 }
828 else
829 {
830 intdig_max = 1;
831 chars_needed = 1 + 1 + fracdig_max;
832 }
833 dig_max = INT_MAX; /* Unlimited. */
834 significant = 1; /* Does not matter here. */
835 }
836 else
837 {
838 dig_max = info->prec < 0 ? 6 : (info->prec == 0 ? 1 : info->prec);
839 if ((expsign == 0 && exponent >= dig_max)
840 || (expsign != 0 && exponent > 4))
841 {
842 if ('g' - 'G' == 'e' - 'E')
843 type = 'E' + (info->spec - 'G');
844 else
845 type = isupper (info->spec) ? 'E' : 'e';
846 fracdig_max = dig_max - 1;
847 intdig_max = 1;
848 chars_needed = 1 + 1 + fracdig_max + 1 + 1 + 4;
849 }
850 else
851 {
852 type = 'f';
853 intdig_max = expsign == 0 ? exponent + 1 : 0;
854 fracdig_max = dig_max - intdig_max;
855 /* We need space for the significant digits and perhaps
856 for leading zeros when < 1.0. The number of leading
857 zeros can be as many as would be required for
858 exponential notation with a negative two-digit
859 exponent, which is 4. */
860 chars_needed = dig_max + 1 + 4;
861 }
862 fracdig_min = info->alt ? fracdig_max : 0;
863 significant = 0; /* We count significant digits. */
864 }
865
866 if (grouping)
867 {
868 /* Guess the number of groups we will make, and thus how
869 many spaces we need for separator characters. */
870 ngroups = __guess_grouping (intdig_max, grouping);
871 chars_needed += ngroups;
872 }
873
874 /* Allocate buffer for output. We need two more because while rounding
875 it is possible that we need two more characters in front of all the
876 other output. If the amount of memory we have to allocate is too
877 large use `malloc' instead of `alloca'. */
878 buffer_malloced = chars_needed > 5000;
879 if (buffer_malloced)
880 {
881 wbuffer = (wchar_t *) malloc ((2 + chars_needed) * sizeof (wchar_t));
882 if (wbuffer == NULL)
883 /* Signal an error to the caller. */
884 return -1;
885 }
886 else
887 wbuffer = (wchar_t *) alloca ((2 + chars_needed) * sizeof (wchar_t));
888 wcp = wstartp = wbuffer + 2; /* Let room for rounding. */
889
890 /* Do the real work: put digits in allocated buffer. */
891 if (expsign == 0 || _tolower (type) != 'f')
892 {
893 assert (expsign == 0 || intdig_max == 1);
894 while (intdig_no < intdig_max)
895 {
896 ++intdig_no;
897 hack_digit_callee = 1;
898 goto hack_digit;
899 hack_digit_callee1:
900 *wcp++ = hack_digit_ret;
901 }
902 significant = 1;
903 if (info->alt
904 || fracdig_min > 0
905 || (fracdig_max > 0 && (fracsize > 1 || frac[0] != 0)))
906 *wcp++ = decimalwc;
907 }
908 else
909 {
910 /* |fp| < 1.0 and the selected type is 'f', so put "0."
911 in the buffer. */
912 *wcp++ = L'0';
913 --exponent;
914 *wcp++ = decimalwc;
915 }
916
917 /* Generate the needed number of fractional digits. */
918 while (fracdig_no < fracdig_min
919 || (fracdig_no < fracdig_max && (fracsize > 1 || frac[0] != 0)))
920 {
921 ++fracdig_no;
922 hack_digit_callee = 2;
923 goto hack_digit;
924 hack_digit_callee2:
925 *wcp = hack_digit_ret;
926 if (*wcp != L'0')
927 significant = 1;
928 else if (significant == 0)
929 {
930 ++fracdig_max;
931 if (fracdig_min > 0)
932 ++fracdig_min;
933 }
934 ++wcp;
935 }
936
937 /* Do rounding. */
938 hack_digit_callee = 3;
939 goto hack_digit;
940 hack_digit_callee3:
941 digit = hack_digit_ret;
942 if (digit > L'4')
943 {
944 wchar_t *wtp = wcp;
945
946 if (digit == L'5'
947 && ((*(wcp - 1) != decimalwc && (*(wcp - 1) & 1) == 0)
948 || ((*(wcp - 1) == decimalwc && (*(wcp - 2) & 1) == 0))))
949 {
950 /* This is the critical case. */
951 if (fracsize == 1 && frac[0] == 0)
952 /* Rest of the number is zero -> round to even.
953 (IEEE 754-1985 4.1 says this is the default rounding.) */
954 goto do_expo;
955 else if (scalesize == 0)
956 {
957 /* Here we have to see whether all limbs are zero since no
958 normalization happened. */
959 size_t lcnt = fracsize;
960 while (lcnt >= 1 && frac[lcnt - 1] == 0)
961 --lcnt;
962 if (lcnt == 0)
963 /* Rest of the number is zero -> round to even.
964 (IEEE 754-1985 4.1 says this is the default rounding.) */
965 goto do_expo;
966 }
967 }
968
969 if (fracdig_no > 0)
970 {
971 /* Process fractional digits. Terminate if not rounded or
972 radix character is reached. */
973 while (*--wtp != decimalwc && *wtp == L'9')
974 *wtp = '0';
975 if (*wtp != decimalwc)
976 /* Round up. */
977 (*wtp)++;
978 }
979
980 if (fracdig_no == 0 || *wtp == decimalwc)
981 {
982 /* Round the integer digits. */
983 if (*(wtp - 1) == decimalwc)
984 --wtp;
985
986 while (--wtp >= wstartp && *wtp == L'9')
987 *wtp = L'0';
988
989 if (wtp >= wstartp)
990 /* Round up. */
991 (*wtp)++;
992 else
993 /* It is more critical. All digits were 9's. */
994 {
995 if (_tolower (type) != 'f')
996 {
997 *wstartp = '1';
998 exponent += expsign == 0 ? 1 : -1;
999 }
1000 else if (intdig_no == dig_max)
1001 {
1002 /* This is the case where for type %g the number fits
1003 really in the range for %f output but after rounding
1004 the number of digits is too big. */
1005 *--wstartp = decimalwc;
1006 *--wstartp = L'1';
1007
1008 if (info->alt || fracdig_no > 0)
1009 {
1010 /* Overwrite the old radix character. */
1011 wstartp[intdig_no + 2] = L'0';
1012 ++fracdig_no;
1013 }
1014
1015 fracdig_no += intdig_no;
1016 intdig_no = 1;
1017 fracdig_max = intdig_max - intdig_no;
1018 ++exponent;
1019 /* Now we must print the exponent. */
1020 type = isupper (info->spec) ? 'E' : 'e';
1021 }
1022 else
1023 {
1024 /* We can simply add another another digit before the
1025 radix. */
1026 *--wstartp = L'1';
1027 ++intdig_no;
1028 }
1029
1030 /* While rounding the number of digits can change.
1031 If the number now exceeds the limits remove some
1032 fractional digits. */
1033 if (intdig_no + fracdig_no > dig_max)
1034 {
1035 wcp -= intdig_no + fracdig_no - dig_max;
1036 fracdig_no -= intdig_no + fracdig_no - dig_max;
1037 }
1038 }
1039 }
1040 }
1041
1042 do_expo:
1043 /* Now remove unnecessary '0' at the end of the string. */
1044 while (fracdig_no > fracdig_min && *(wcp - 1) == L'0')
1045 {
1046 --wcp;
1047 --fracdig_no;
1048 }
1049 /* If we eliminate all fractional digits we perhaps also can remove
1050 the radix character. */
1051 if (fracdig_no == 0 && !info->alt && *(wcp - 1) == decimalwc)
1052 --wcp;
1053
1054 if (grouping)
1055 /* Add in separator characters, overwriting the same buffer. */
1056 wcp = group_number (wstartp, wcp, intdig_no, grouping, thousands_sepwc,
1057 ngroups);
1058
1059 /* Write the exponent if it is needed. */
1060 if (_tolower (type) != 'f')
1061 {
1062 *wcp++ = (wchar_t) type;
1063 *wcp++ = expsign ? L'-' : L'+';
1064
1065 /* Find the magnitude of the exponent. */
1066 expscale = 10;
1067 while (expscale <= exponent)
1068 expscale *= 10;
1069
1070 if (exponent < 10)
1071 /* Exponent always has at least two digits. */
1072 *wcp++ = L'0';
1073 else
1074 do
1075 {
1076 expscale /= 10;
1077 *wcp++ = L'0' + (exponent / expscale);
1078 exponent %= expscale;
1079 }
1080 while (expscale > 10);
1081 *wcp++ = L'0' + exponent;
1082 }
1083
1084 /* Compute number of characters which must be filled with the padding
1085 character. */
1086 if (is_neg || info->showsign || info->space)
1087 --width;
1088 width -= wcp - wstartp;
1089
1090 if (!info->left && info->pad != '0' && width > 0)
1091 PADN (info->pad, width);
1092
1093 if (is_neg)
1094 outchar ('-');
1095 else if (info->showsign)
1096 outchar ('+');
1097 else if (info->space)
1098 outchar (' ');
1099
1100 if (!info->left && info->pad == '0' && width > 0)
1101 PADN ('0', width);
1102
1103 {
1104 char *buffer = NULL;
1105 char *cp = NULL;
1106 char *tmpptr;
1107
1108 if (! wide)
1109 {
1110 /* Create the single byte string. */
1111 size_t decimal_len;
1112 size_t thousands_sep_len;
1113 wchar_t *copywc;
1114
1115 decimal_len = strlen (decimal);
1116
1117 if (thousands_sep == NULL)
1118 thousands_sep_len = 0;
1119 else
1120 thousands_sep_len = strlen (thousands_sep);
1121
1122 if (buffer_malloced)
1123 {
1124 buffer = (char *) malloc (2 + chars_needed + decimal_len
1125 + ngroups * thousands_sep_len);
1126 if (buffer == NULL)
1127 /* Signal an error to the caller. */
1128 return -1;
1129 }
1130 else
1131 buffer = (char *) alloca (2 + chars_needed + decimal_len
1132 + ngroups * thousands_sep_len);
1133
1134 /* Now copy the wide character string. Since the character
1135 (except for the decimal point and thousands separator) must
1136 be coming from the ASCII range we can esily convert the
1137 string without mapping tables. */
1138 for (cp = buffer, copywc = wstartp; copywc < wcp; ++copywc)
1139 if (*copywc == decimalwc)
1140 cp = (char *) __mempcpy (cp, decimal, decimal_len);
1141 else if (*copywc == thousands_sepwc)
1142 cp = (char *) __mempcpy (cp, thousands_sep, thousands_sep_len);
1143 else
1144 *cp++ = (char) *copywc;
1145 }
1146
1147 tmpptr = buffer;
1148 PRINT (tmpptr, wstartp, wide ? wcp - wstartp : cp - tmpptr);
1149
1150 /* Free the memory if necessary. */
1151 if (buffer_malloced)
1152 {
1153 free (buffer);
1154 free (wbuffer);
1155 }
1156 }
1157
1158 if (info->left && width > 0)
1159 PADN (info->pad, width);
1160 }
1161 return done;
1162 }
1163
1164 /* Return the number of extra grouping characters that will be inserted
1165 into a number with INTDIG_MAX integer digits. */
1166
1167 unsigned int
__guess_grouping(unsigned int intdig_max,const char * grouping)1168 __guess_grouping (unsigned int intdig_max, const char *grouping)
1169 {
1170 unsigned int groups;
1171
1172 /* We treat all negative values like CHAR_MAX. */
1173
1174 if (*grouping == CHAR_MAX || *grouping <= 0)
1175 /* No grouping should be done. */
1176 return 0;
1177
1178 groups = 0;
1179 while (intdig_max > (unsigned int) *grouping)
1180 {
1181 ++groups;
1182 intdig_max -= *grouping++;
1183
1184 if (*grouping == CHAR_MAX
1185 #if CHAR_MIN < 0
1186 || *grouping < 0
1187 #endif
1188 )
1189 /* No more grouping should be done. */
1190 break;
1191 else if (*grouping == 0)
1192 {
1193 /* Same grouping repeats. */
1194 groups += (intdig_max - 1) / grouping[-1];
1195 break;
1196 }
1197 }
1198
1199 return groups;
1200 }
1201
1202 /* Group the INTDIG_NO integer digits of the number in [BUF,BUFEND).
1203 There is guaranteed enough space past BUFEND to extend it.
1204 Return the new end of buffer. */
1205
1206 static wchar_t *
1207 internal_function
group_number(wchar_t * buf,wchar_t * bufend,unsigned int intdig_no,const char * grouping,wchar_t thousands_sep,int ngroups)1208 group_number (wchar_t *buf, wchar_t *bufend, unsigned int intdig_no,
1209 const char *grouping, wchar_t thousands_sep, int ngroups)
1210 {
1211 wchar_t *p;
1212
1213 if (ngroups == 0)
1214 return bufend;
1215
1216 /* Move the fractional part down. */
1217 __wmemmove (buf + intdig_no + ngroups, buf + intdig_no,
1218 bufend - (buf + intdig_no));
1219
1220 p = buf + intdig_no + ngroups - 1;
1221 do
1222 {
1223 unsigned int len = *grouping++;
1224 do
1225 *p-- = buf[--intdig_no];
1226 while (--len > 0);
1227 *p-- = thousands_sep;
1228
1229 if (*grouping == CHAR_MAX
1230 #if CHAR_MIN < 0
1231 || *grouping < 0
1232 #endif
1233 )
1234 /* No more grouping should be done. */
1235 break;
1236 else if (*grouping == 0)
1237 /* Same grouping repeats. */
1238 --grouping;
1239 } while (intdig_no > (unsigned int) *grouping);
1240
1241 /* Copy the remaining ungrouped digits. */
1242 do
1243 *p-- = buf[--intdig_no];
1244 while (p > buf);
1245
1246 return bufend + ngroups;
1247 }
1248