Daily bump.
[official-gcc.git] / libquadmath / printf / printf_fp.c
blob8effcee88faaee1519d9fca05fb08102ea91b25d
1 /* Floating point output for `printf'.
2 Copyright (C) 1995-2012 Free Software Foundation, Inc.
4 This file is part of the GNU C Library.
5 Written by Ulrich Drepper <drepper@gnu.ai.mit.edu>, 1995.
7 The GNU C Library is free software; you can redistribute it and/or
8 modify it under the terms of the GNU Lesser General Public
9 License as published by the Free Software Foundation; either
10 version 2.1 of the License, or (at your option) any later version.
12 The GNU C Library is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
15 Lesser General Public License for more details.
17 You should have received a copy of the GNU Lesser General Public
18 License along with the GNU C Library; if not, see
19 <http://www.gnu.org/licenses/>. */
21 #include <config.h>
22 #include <float.h>
23 #include <limits.h>
24 #include <math.h>
25 #include <string.h>
26 #include <unistd.h>
27 #include <stdlib.h>
28 #include <stdbool.h>
29 #define NDEBUG
30 #include <assert.h>
31 #ifdef HAVE_ERRNO_H
32 #include <errno.h>
33 #endif
34 #include <stdio.h>
35 #include <stdarg.h>
36 #ifdef HAVE_FENV_H
37 #include "quadmath-rounding-mode.h"
38 #endif
39 #include "quadmath-printf.h"
40 #include "fpioconst.h"
42 #ifdef USE_I18N_NUMBER_H
43 #include "_i18n_number.h"
44 #endif
47 /* Macros for doing the actual output. */
49 #define outchar(ch) \
50 do \
51 { \
52 register const int outc = (ch); \
53 if (PUTC (outc, fp) == EOF) \
54 { \
55 if (buffer_malloced) \
56 free (wbuffer); \
57 return -1; \
58 } \
59 ++done; \
60 } while (0)
62 #define PRINT(ptr, wptr, len) \
63 do \
64 { \
65 register size_t outlen = (len); \
66 if (len > 20) \
67 { \
68 if (PUT (fp, wide ? (const char *) wptr : ptr, outlen) != outlen) \
69 { \
70 if (buffer_malloced) \
71 free (wbuffer); \
72 return -1; \
73 } \
74 ptr += outlen; \
75 done += outlen; \
76 } \
77 else \
78 { \
79 if (wide) \
80 while (outlen-- > 0) \
81 outchar (*wptr++); \
82 else \
83 while (outlen-- > 0) \
84 outchar (*ptr++); \
85 } \
86 } while (0)
88 #define PADN(ch, len) \
89 do \
90 { \
91 if (PAD (fp, ch, len) != len) \
92 { \
93 if (buffer_malloced) \
94 free (wbuffer); \
95 return -1; \
96 } \
97 done += len; \
98 } \
99 while (0)
102 /* We use the GNU MP library to handle large numbers.
104 An MP variable occupies a varying number of entries in its array. We keep
105 track of this number for efficiency reasons. Otherwise we would always
106 have to process the whole array. */
107 #define MPN_VAR(name) mp_limb_t *name; mp_size_t name##size
109 #define MPN_ASSIGN(dst,src) \
110 memcpy (dst, src, (dst##size = src##size) * sizeof (mp_limb_t))
111 #define MPN_GE(u,v) \
112 (u##size > v##size || (u##size == v##size && mpn_cmp (u, v, u##size) >= 0))
114 extern mp_size_t mpn_extract_flt128 (mp_ptr res_ptr, mp_size_t size,
115 int *expt, int *is_neg,
116 __float128 value) attribute_hidden;
117 static unsigned int guess_grouping (unsigned int intdig_max,
118 const char *grouping);
121 static wchar_t *group_number (wchar_t *buf, wchar_t *bufend,
122 unsigned int intdig_no, const char *grouping,
123 wchar_t thousands_sep, int ngroups);
127 __quadmath_printf_fp (struct __quadmath_printf_file *fp,
128 const struct printf_info *info,
129 const void *const *args)
131 /* The floating-point value to output. */
132 __float128 fpnum;
134 /* Locale-dependent representation of decimal point. */
135 const char *decimal;
136 wchar_t decimalwc;
138 /* Locale-dependent thousands separator and grouping specification. */
139 const char *thousands_sep = NULL;
140 wchar_t thousands_sepwc = L_('\0');
141 const char *grouping;
143 /* "NaN" or "Inf" for the special cases. */
144 const char *special = NULL;
145 const wchar_t *wspecial = NULL;
147 /* We need just a few limbs for the input before shifting to the right
148 position. */
149 mp_limb_t fp_input[(FLT128_MANT_DIG + BITS_PER_MP_LIMB - 1) / BITS_PER_MP_LIMB];
150 /* We need to shift the contents of fp_input by this amount of bits. */
151 int to_shift = 0;
153 /* The fraction of the floting-point value in question */
154 MPN_VAR(frac);
155 /* and the exponent. */
156 int exponent;
157 /* Sign of the exponent. */
158 int expsign = 0;
159 /* Sign of float number. */
160 int is_neg = 0;
162 /* Scaling factor. */
163 MPN_VAR(scale);
165 /* Temporary bignum value. */
166 MPN_VAR(tmp);
168 /* Digit which is result of last hack_digit() call. */
169 wchar_t last_digit, next_digit;
170 bool more_bits;
172 /* The type of output format that will be used: 'e'/'E' or 'f'. */
173 int type;
175 /* Counter for number of written characters. */
176 int done = 0;
178 /* General helper (carry limb). */
179 mp_limb_t cy;
181 /* Nonzero if this is output on a wide character stream. */
182 int wide = info->wide;
184 /* Buffer in which we produce the output. */
185 wchar_t *wbuffer = NULL;
186 /* Flag whether wbuffer is malloc'ed or not. */
187 int buffer_malloced = 0;
189 auto wchar_t hack_digit (void);
191 wchar_t hack_digit (void)
193 mp_limb_t hi;
195 if (expsign != 0 && type == 'f' && exponent-- > 0)
196 hi = 0;
197 else if (scalesize == 0)
199 hi = frac[fracsize - 1];
200 frac[fracsize - 1] = mpn_mul_1 (frac, frac, fracsize - 1, 10);
202 else
204 if (fracsize < scalesize)
205 hi = 0;
206 else
208 hi = mpn_divmod (tmp, frac, fracsize, scale, scalesize);
209 tmp[fracsize - scalesize] = hi;
210 hi = tmp[0];
212 fracsize = scalesize;
213 while (fracsize != 0 && frac[fracsize - 1] == 0)
214 --fracsize;
215 if (fracsize == 0)
217 /* We're not prepared for an mpn variable with zero
218 limbs. */
219 fracsize = 1;
220 return L_('0') + hi;
224 mp_limb_t _cy = mpn_mul_1 (frac, frac, fracsize, 10);
225 if (_cy != 0)
226 frac[fracsize++] = _cy;
229 return L_('0') + hi;
232 /* Figure out the decimal point character. */
233 #ifdef USE_NL_LANGINFO
234 if (info->extra == 0)
235 decimal = nl_langinfo (DECIMAL_POINT);
236 else
238 decimal = nl_langinfo (MON_DECIMAL_POINT);
239 if (*decimal == '\0')
240 decimal = nl_langinfo (DECIMAL_POINT);
242 /* The decimal point character must never be zero. */
243 assert (*decimal != '\0');
244 #elif defined USE_LOCALECONV
245 const struct lconv *lc = localeconv ();
246 if (info->extra == 0)
247 decimal = lc->decimal_point;
248 else
250 decimal = lc->mon_decimal_point;
251 if (decimal == NULL || *decimal == '\0')
252 decimal = lc->decimal_point;
254 if (decimal == NULL || *decimal == '\0')
255 decimal = ".";
256 #else
257 decimal = ".";
258 #endif
259 #ifdef USE_NL_LANGINFO_WC
260 if (info->extra == 0)
261 decimalwc = nl_langinfo_wc (_NL_NUMERIC_DECIMAL_POINT_WC);
262 else
264 decimalwc = nl_langinfo_wc (_NL_MONETARY_DECIMAL_POINT_WC);
265 if (decimalwc == L_('\0'))
266 decimalwc = nl_langinfo_wc (_NL_NUMERIC_DECIMAL_POINT_WC);
268 /* The decimal point character must never be zero. */
269 assert (decimalwc != L_('\0'));
270 #else
271 decimalwc = L_('.');
272 #endif
274 #if defined USE_NL_LANGINFO && defined USE_NL_LANGINFO_WC
275 if (info->group)
277 if (info->extra == 0)
278 grouping = nl_langinfo (GROUPING);
279 else
280 grouping = nl_langinfo (MON_GROUPING);
282 if (*grouping <= 0 || *grouping == CHAR_MAX)
283 grouping = NULL;
284 else
286 /* Figure out the thousands separator character. */
287 if (wide)
289 if (info->extra == 0)
290 thousands_sepwc = nl_langinfo_wc (_NL_NUMERIC_THOUSANDS_SEP_WC);
291 else
292 thousands_sepwc = nl_langinfo_wc (_NL_MONETARY_THOUSANDS_SEP_WC);
294 if (thousands_sepwc == L_('\0'))
295 grouping = NULL;
297 else
299 if (info->extra == 0)
300 thousands_sep = nl_langinfo (THOUSANDS_SEP);
301 else
302 thousands_sep = nl_langinfo (MON_THOUSANDS_SEP);
303 if (*thousands_sep == '\0')
304 grouping = NULL;
308 else
309 #elif defined USE_NL_LANGINFO
310 if (info->group && !wide)
312 if (info->extra == 0)
313 grouping = nl_langinfo (GROUPING);
314 else
315 grouping = nl_langinfo (MON_GROUPING);
317 if (*grouping <= 0 || *grouping == CHAR_MAX)
318 grouping = NULL;
319 else
321 /* Figure out the thousands separator character. */
322 if (info->extra == 0)
323 thousands_sep = nl_langinfo (THOUSANDS_SEP);
324 else
325 thousands_sep = nl_langinfo (MON_THOUSANDS_SEP);
327 if (*thousands_sep == '\0')
328 grouping = NULL;
331 else
332 #elif defined USE_LOCALECONV
333 if (info->group && !wide)
335 if (info->extra == 0)
336 grouping = lc->grouping;
337 else
338 grouping = lc->mon_grouping;
340 if (grouping == NULL || *grouping <= 0 || *grouping == CHAR_MAX)
341 grouping = NULL;
342 else
344 /* Figure out the thousands separator character. */
345 if (info->extra == 0)
346 thousands_sep = lc->thousands_sep;
347 else
348 thousands_sep = lc->mon_thousands_sep;
350 if (thousands_sep == NULL || *thousands_sep == '\0')
351 grouping = NULL;
354 else
355 #endif
356 grouping = NULL;
357 if (grouping != NULL && !wide)
358 /* If we are printing multibyte characters and there is a
359 multibyte representation for the thousands separator,
360 we must ensure the wide character thousands separator
361 is available, even if it is fake. */
362 thousands_sepwc = (wchar_t) 0xfffffffe;
364 /* Fetch the argument value. */
366 fpnum = **(const __float128 **) args[0];
368 /* Check for special values: not a number or infinity. */
369 if (isnanq (fpnum))
371 ieee854_float128 u = { .value = fpnum };
372 is_neg = u.ieee.negative != 0;
373 if (isupper (info->spec))
375 special = "NAN";
376 wspecial = L_("NAN");
378 else
380 special = "nan";
381 wspecial = L_("nan");
384 else if (isinfq (fpnum))
386 is_neg = fpnum < 0;
387 if (isupper (info->spec))
389 special = "INF";
390 wspecial = L_("INF");
392 else
394 special = "inf";
395 wspecial = L_("inf");
398 else
400 fracsize = mpn_extract_flt128 (fp_input,
401 (sizeof (fp_input) /
402 sizeof (fp_input[0])),
403 &exponent, &is_neg, fpnum);
404 to_shift = 1 + fracsize * BITS_PER_MP_LIMB - FLT128_MANT_DIG;
408 if (special)
410 int width = info->width;
412 if (is_neg || info->showsign || info->space)
413 --width;
414 width -= 3;
416 if (!info->left && width > 0)
417 PADN (' ', width);
419 if (is_neg)
420 outchar ('-');
421 else if (info->showsign)
422 outchar ('+');
423 else if (info->space)
424 outchar (' ');
426 PRINT (special, wspecial, 3);
428 if (info->left && width > 0)
429 PADN (' ', width);
431 return done;
435 /* We need three multiprecision variables. Now that we have the exponent
436 of the number we can allocate the needed memory. It would be more
437 efficient to use variables of the fixed maximum size but because this
438 would be really big it could lead to memory problems. */
440 mp_size_t bignum_size = ((ABS (exponent) + BITS_PER_MP_LIMB - 1)
441 / BITS_PER_MP_LIMB
442 + (FLT128_MANT_DIG / BITS_PER_MP_LIMB > 2 ? 8 : 4))
443 * sizeof (mp_limb_t);
444 frac = (mp_limb_t *) alloca (bignum_size);
445 tmp = (mp_limb_t *) alloca (bignum_size);
446 scale = (mp_limb_t *) alloca (bignum_size);
449 /* We now have to distinguish between numbers with positive and negative
450 exponents because the method used for the one is not applicable/efficient
451 for the other. */
452 scalesize = 0;
453 if (exponent > 2)
455 /* |FP| >= 8.0. */
456 int scaleexpo = 0;
457 int explog = FLT128_MAX_10_EXP_LOG;
458 int exp10 = 0;
459 const struct mp_power *powers = &_fpioconst_pow10[explog + 1];
460 int cnt_h, cnt_l, i;
462 if ((exponent + to_shift) % BITS_PER_MP_LIMB == 0)
464 MPN_COPY_DECR (frac + (exponent + to_shift) / BITS_PER_MP_LIMB,
465 fp_input, fracsize);
466 fracsize += (exponent + to_shift) / BITS_PER_MP_LIMB;
468 else
470 cy = mpn_lshift (frac + (exponent + to_shift) / BITS_PER_MP_LIMB,
471 fp_input, fracsize,
472 (exponent + to_shift) % BITS_PER_MP_LIMB);
473 fracsize += (exponent + to_shift) / BITS_PER_MP_LIMB;
474 if (cy)
475 frac[fracsize++] = cy;
477 MPN_ZERO (frac, (exponent + to_shift) / BITS_PER_MP_LIMB);
479 assert (powers > &_fpioconst_pow10[0]);
482 --powers;
484 /* The number of the product of two binary numbers with n and m
485 bits respectively has m+n or m+n-1 bits. */
486 if (exponent >= scaleexpo + powers->p_expo - 1)
488 if (scalesize == 0)
490 if (FLT128_MANT_DIG > _FPIO_CONST_OFFSET * BITS_PER_MP_LIMB)
492 #define _FPIO_CONST_SHIFT \
493 (((FLT128_MANT_DIG + BITS_PER_MP_LIMB - 1) / BITS_PER_MP_LIMB) \
494 - _FPIO_CONST_OFFSET)
495 /* 64bit const offset is not enough for
496 IEEE quad long double. */
497 tmpsize = powers->arraysize + _FPIO_CONST_SHIFT;
498 memcpy (tmp + _FPIO_CONST_SHIFT,
499 &__tens[powers->arrayoff],
500 tmpsize * sizeof (mp_limb_t));
501 MPN_ZERO (tmp, _FPIO_CONST_SHIFT);
502 /* Adjust exponent, as scaleexpo will be this much
503 bigger too. */
504 exponent += _FPIO_CONST_SHIFT * BITS_PER_MP_LIMB;
506 else
508 tmpsize = powers->arraysize;
509 memcpy (tmp, &__tens[powers->arrayoff],
510 tmpsize * sizeof (mp_limb_t));
513 else
515 cy = mpn_mul (tmp, scale, scalesize,
516 &__tens[powers->arrayoff
517 + _FPIO_CONST_OFFSET],
518 powers->arraysize - _FPIO_CONST_OFFSET);
519 tmpsize = scalesize + powers->arraysize - _FPIO_CONST_OFFSET;
520 if (cy == 0)
521 --tmpsize;
524 if (MPN_GE (frac, tmp))
526 int cnt;
527 MPN_ASSIGN (scale, tmp);
528 count_leading_zeros (cnt, scale[scalesize - 1]);
529 scaleexpo = (scalesize - 2) * BITS_PER_MP_LIMB - cnt - 1;
530 exp10 |= 1 << explog;
533 --explog;
535 while (powers > &_fpioconst_pow10[0]);
536 exponent = exp10;
538 /* Optimize number representations. We want to represent the numbers
539 with the lowest number of bytes possible without losing any
540 bytes. Also the highest bit in the scaling factor has to be set
541 (this is a requirement of the MPN division routines). */
542 if (scalesize > 0)
544 /* Determine minimum number of zero bits at the end of
545 both numbers. */
546 for (i = 0; scale[i] == 0 && frac[i] == 0; i++)
549 /* Determine number of bits the scaling factor is misplaced. */
550 count_leading_zeros (cnt_h, scale[scalesize - 1]);
552 if (cnt_h == 0)
554 /* The highest bit of the scaling factor is already set. So
555 we only have to remove the trailing empty limbs. */
556 if (i > 0)
558 MPN_COPY_INCR (scale, scale + i, scalesize - i);
559 scalesize -= i;
560 MPN_COPY_INCR (frac, frac + i, fracsize - i);
561 fracsize -= i;
564 else
566 if (scale[i] != 0)
568 count_trailing_zeros (cnt_l, scale[i]);
569 if (frac[i] != 0)
571 int cnt_l2;
572 count_trailing_zeros (cnt_l2, frac[i]);
573 if (cnt_l2 < cnt_l)
574 cnt_l = cnt_l2;
577 else
578 count_trailing_zeros (cnt_l, frac[i]);
580 /* Now shift the numbers to their optimal position. */
581 if (i == 0 && BITS_PER_MP_LIMB - cnt_h > cnt_l)
583 /* We cannot save any memory. So just roll both numbers
584 so that the scaling factor has its highest bit set. */
586 (void) mpn_lshift (scale, scale, scalesize, cnt_h);
587 cy = mpn_lshift (frac, frac, fracsize, cnt_h);
588 if (cy != 0)
589 frac[fracsize++] = cy;
591 else if (BITS_PER_MP_LIMB - cnt_h <= cnt_l)
593 /* We can save memory by removing the trailing zero limbs
594 and by packing the non-zero limbs which gain another
595 free one. */
597 (void) mpn_rshift (scale, scale + i, scalesize - i,
598 BITS_PER_MP_LIMB - cnt_h);
599 scalesize -= i + 1;
600 (void) mpn_rshift (frac, frac + i, fracsize - i,
601 BITS_PER_MP_LIMB - cnt_h);
602 fracsize -= frac[fracsize - i - 1] == 0 ? i + 1 : i;
604 else
606 /* We can only save the memory of the limbs which are zero.
607 The non-zero parts occupy the same number of limbs. */
609 (void) mpn_rshift (scale, scale + (i - 1),
610 scalesize - (i - 1),
611 BITS_PER_MP_LIMB - cnt_h);
612 scalesize -= i;
613 (void) mpn_rshift (frac, frac + (i - 1),
614 fracsize - (i - 1),
615 BITS_PER_MP_LIMB - cnt_h);
616 fracsize -= frac[fracsize - (i - 1) - 1] == 0 ? i : i - 1;
621 else if (exponent < 0)
623 /* |FP| < 1.0. */
624 int exp10 = 0;
625 int explog = FLT128_MAX_10_EXP_LOG;
626 const struct mp_power *powers = &_fpioconst_pow10[explog + 1];
628 /* Now shift the input value to its right place. */
629 cy = mpn_lshift (frac, fp_input, fracsize, to_shift);
630 frac[fracsize++] = cy;
631 assert (cy == 1 || (frac[fracsize - 2] == 0 && frac[0] == 0));
633 expsign = 1;
634 exponent = -exponent;
636 assert (powers != &_fpioconst_pow10[0]);
639 --powers;
641 if (exponent >= powers->m_expo)
643 int i, incr, cnt_h, cnt_l;
644 mp_limb_t topval[2];
646 /* The mpn_mul function expects the first argument to be
647 bigger than the second. */
648 if (fracsize < powers->arraysize - _FPIO_CONST_OFFSET)
649 cy = mpn_mul (tmp, &__tens[powers->arrayoff
650 + _FPIO_CONST_OFFSET],
651 powers->arraysize - _FPIO_CONST_OFFSET,
652 frac, fracsize);
653 else
654 cy = mpn_mul (tmp, frac, fracsize,
655 &__tens[powers->arrayoff + _FPIO_CONST_OFFSET],
656 powers->arraysize - _FPIO_CONST_OFFSET);
657 tmpsize = fracsize + powers->arraysize - _FPIO_CONST_OFFSET;
658 if (cy == 0)
659 --tmpsize;
661 count_leading_zeros (cnt_h, tmp[tmpsize - 1]);
662 incr = (tmpsize - fracsize) * BITS_PER_MP_LIMB
663 + BITS_PER_MP_LIMB - 1 - cnt_h;
665 assert (incr <= powers->p_expo);
667 /* If we increased the exponent by exactly 3 we have to test
668 for overflow. This is done by comparing with 10 shifted
669 to the right position. */
670 if (incr == exponent + 3)
672 if (cnt_h <= BITS_PER_MP_LIMB - 4)
674 topval[0] = 0;
675 topval[1]
676 = ((mp_limb_t) 10) << (BITS_PER_MP_LIMB - 4 - cnt_h);
678 else
680 topval[0] = ((mp_limb_t) 10) << (BITS_PER_MP_LIMB - 4);
681 topval[1] = 0;
682 (void) mpn_lshift (topval, topval, 2,
683 BITS_PER_MP_LIMB - cnt_h);
687 /* We have to be careful when multiplying the last factor.
688 If the result is greater than 1.0 be have to test it
689 against 10.0. If it is greater or equal to 10.0 the
690 multiplication was not valid. This is because we cannot
691 determine the number of bits in the result in advance. */
692 if (incr < exponent + 3
693 || (incr == exponent + 3 &&
694 (tmp[tmpsize - 1] < topval[1]
695 || (tmp[tmpsize - 1] == topval[1]
696 && tmp[tmpsize - 2] < topval[0]))))
698 /* The factor is right. Adapt binary and decimal
699 exponents. */
700 exponent -= incr;
701 exp10 |= 1 << explog;
703 /* If this factor yields a number greater or equal to
704 1.0, we must not shift the non-fractional digits down. */
705 if (exponent < 0)
706 cnt_h += -exponent;
708 /* Now we optimize the number representation. */
709 for (i = 0; tmp[i] == 0; ++i);
710 if (cnt_h == BITS_PER_MP_LIMB - 1)
712 MPN_COPY (frac, tmp + i, tmpsize - i);
713 fracsize = tmpsize - i;
715 else
717 count_trailing_zeros (cnt_l, tmp[i]);
719 /* Now shift the numbers to their optimal position. */
720 if (i == 0 && BITS_PER_MP_LIMB - 1 - cnt_h > cnt_l)
722 /* We cannot save any memory. Just roll the
723 number so that the leading digit is in a
724 separate limb. */
726 cy = mpn_lshift (frac, tmp, tmpsize, cnt_h + 1);
727 fracsize = tmpsize + 1;
728 frac[fracsize - 1] = cy;
730 else if (BITS_PER_MP_LIMB - 1 - cnt_h <= cnt_l)
732 (void) mpn_rshift (frac, tmp + i, tmpsize - i,
733 BITS_PER_MP_LIMB - 1 - cnt_h);
734 fracsize = tmpsize - i;
736 else
738 /* We can only save the memory of the limbs which
739 are zero. The non-zero parts occupy the same
740 number of limbs. */
742 (void) mpn_rshift (frac, tmp + (i - 1),
743 tmpsize - (i - 1),
744 BITS_PER_MP_LIMB - 1 - cnt_h);
745 fracsize = tmpsize - (i - 1);
750 --explog;
752 while (powers != &_fpioconst_pow10[1] && exponent > 0);
753 /* All factors but 10^-1 are tested now. */
754 if (exponent > 0)
756 int cnt_l;
758 cy = mpn_mul_1 (tmp, frac, fracsize, 10);
759 tmpsize = fracsize;
760 assert (cy == 0 || tmp[tmpsize - 1] < 20);
762 count_trailing_zeros (cnt_l, tmp[0]);
763 if (cnt_l < MIN (4, exponent))
765 cy = mpn_lshift (frac, tmp, tmpsize,
766 BITS_PER_MP_LIMB - MIN (4, exponent));
767 if (cy != 0)
768 frac[tmpsize++] = cy;
770 else
771 (void) mpn_rshift (frac, tmp, tmpsize, MIN (4, exponent));
772 fracsize = tmpsize;
773 exp10 |= 1;
774 assert (frac[fracsize - 1] < 10);
776 exponent = exp10;
778 else
780 /* This is a special case. We don't need a factor because the
781 numbers are in the range of 1.0 <= |fp| < 8.0. We simply
782 shift it to the right place and divide it by 1.0 to get the
783 leading digit. (Of course this division is not really made.) */
784 assert (0 <= exponent && exponent < 3 &&
785 exponent + to_shift < BITS_PER_MP_LIMB);
787 /* Now shift the input value to its right place. */
788 cy = mpn_lshift (frac, fp_input, fracsize, (exponent + to_shift));
789 frac[fracsize++] = cy;
790 exponent = 0;
794 int width = info->width;
795 wchar_t *wstartp, *wcp;
796 size_t chars_needed;
797 int expscale;
798 int intdig_max, intdig_no = 0;
799 int fracdig_min;
800 int fracdig_max;
801 int dig_max;
802 int significant;
803 int ngroups = 0;
804 char spec = tolower (info->spec);
806 if (spec == 'e')
808 type = info->spec;
809 intdig_max = 1;
810 fracdig_min = fracdig_max = info->prec < 0 ? 6 : info->prec;
811 chars_needed = 1 + 1 + (size_t) fracdig_max + 1 + 1 + 4;
812 /* d . ddd e +- ddd */
813 dig_max = __INT_MAX__; /* Unlimited. */
814 significant = 1; /* Does not matter here. */
816 else if (spec == 'f')
818 type = 'f';
819 fracdig_min = fracdig_max = info->prec < 0 ? 6 : info->prec;
820 dig_max = __INT_MAX__; /* Unlimited. */
821 significant = 1; /* Does not matter here. */
822 if (expsign == 0)
824 intdig_max = exponent + 1;
825 /* This can be really big! */ /* XXX Maybe malloc if too big? */
826 chars_needed = (size_t) exponent + 1 + 1 + (size_t) fracdig_max;
828 else
830 intdig_max = 1;
831 chars_needed = 1 + 1 + (size_t) fracdig_max;
834 else
836 dig_max = info->prec < 0 ? 6 : (info->prec == 0 ? 1 : info->prec);
837 if ((expsign == 0 && exponent >= dig_max)
838 || (expsign != 0 && exponent > 4))
840 if ('g' - 'G' == 'e' - 'E')
841 type = 'E' + (info->spec - 'G');
842 else
843 type = isupper (info->spec) ? 'E' : 'e';
844 fracdig_max = dig_max - 1;
845 intdig_max = 1;
846 chars_needed = 1 + 1 + (size_t) fracdig_max + 1 + 1 + 4;
848 else
850 type = 'f';
851 intdig_max = expsign == 0 ? exponent + 1 : 0;
852 fracdig_max = dig_max - intdig_max;
853 /* We need space for the significant digits and perhaps
854 for leading zeros when < 1.0. The number of leading
855 zeros can be as many as would be required for
856 exponential notation with a negative two-digit
857 exponent, which is 4. */
858 chars_needed = (size_t) dig_max + 1 + 4;
860 fracdig_min = info->alt ? fracdig_max : 0;
861 significant = 0; /* We count significant digits. */
864 if (grouping)
866 /* Guess the number of groups we will make, and thus how
867 many spaces we need for separator characters. */
868 ngroups = guess_grouping (intdig_max, grouping);
869 /* Allocate one more character in case rounding increases the
870 number of groups. */
871 chars_needed += ngroups + 1;
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 if (__builtin_expect (chars_needed >= (size_t) -1 / sizeof (wchar_t) - 2
879 || chars_needed < fracdig_max, 0))
881 /* Some overflow occurred. */
882 #if defined HAVE_ERRNO_H && defined ERANGE
883 errno = ERANGE;
884 #endif
885 return -1;
887 size_t wbuffer_to_alloc = (2 + chars_needed) * sizeof (wchar_t);
888 buffer_malloced = wbuffer_to_alloc >= 4096;
889 if (__builtin_expect (buffer_malloced, 0))
891 wbuffer = (wchar_t *) malloc (wbuffer_to_alloc);
892 if (wbuffer == NULL)
893 /* Signal an error to the caller. */
894 return -1;
896 else
897 wbuffer = (wchar_t *) alloca (wbuffer_to_alloc);
898 wcp = wstartp = wbuffer + 2; /* Let room for rounding. */
900 /* Do the real work: put digits in allocated buffer. */
901 if (expsign == 0 || type != 'f')
903 assert (expsign == 0 || intdig_max == 1);
904 while (intdig_no < intdig_max)
906 ++intdig_no;
907 *wcp++ = hack_digit ();
909 significant = 1;
910 if (info->alt
911 || fracdig_min > 0
912 || (fracdig_max > 0 && (fracsize > 1 || frac[0] != 0)))
913 *wcp++ = decimalwc;
915 else
917 /* |fp| < 1.0 and the selected type is 'f', so put "0."
918 in the buffer. */
919 *wcp++ = L_('0');
920 --exponent;
921 *wcp++ = decimalwc;
924 /* Generate the needed number of fractional digits. */
925 int fracdig_no = 0;
926 int added_zeros = 0;
927 while (fracdig_no < fracdig_min + added_zeros
928 || (fracdig_no < fracdig_max && (fracsize > 1 || frac[0] != 0)))
930 ++fracdig_no;
931 *wcp = hack_digit ();
932 if (*wcp++ != L_('0'))
933 significant = 1;
934 else if (significant == 0)
936 ++fracdig_max;
937 if (fracdig_min > 0)
938 ++added_zeros;
942 /* Do rounding. */
943 last_digit = wcp[-1] != decimalwc ? wcp[-1] : wcp[-2];
944 next_digit =hack_digit ();
946 if (next_digit != L_('0') && next_digit != L_('5'))
947 more_bits = true;
948 else if (fracsize == 1 && frac[0] == 0)
949 /* Rest of the number is zero. */
950 more_bits = false;
951 else if (scalesize == 0)
953 /* Here we have to see whether all limbs are zero since no
954 normalization happened. */
955 size_t lcnt = fracsize;
956 while (lcnt >= 1 && frac[lcnt - 1] == 0)
957 --lcnt;
958 more_bits = lcnt > 0;
960 else
961 more_bits = true;
963 #ifdef HAVE_FENV_H
964 int rounding_mode = get_rounding_mode ();
965 if (round_away (is_neg, (last_digit - L_('0')) & 1, next_digit >= L_('5'),
966 more_bits, rounding_mode))
968 wchar_t *wtp = wcp;
970 if (fracdig_no > 0)
972 /* Process fractional digits. Terminate if not rounded or
973 radix character is reached. */
974 int removed = 0;
975 while (*--wtp != decimalwc && *wtp == L_('9'))
977 *wtp = L_('0');
978 ++removed;
980 if (removed == fracdig_min && added_zeros > 0)
981 --added_zeros;
982 if (*wtp != decimalwc)
983 /* Round up. */
984 (*wtp)++;
985 else if (__builtin_expect (spec == 'g' && type == 'f' && info->alt
986 && wtp == wstartp + 1
987 && wstartp[0] == L_('0'),
989 /* This is a special case: the rounded number is 1.0,
990 the format is 'g' or 'G', and the alternative format
991 is selected. This means the result must be "1.". */
992 --added_zeros;
995 if (fracdig_no == 0 || *wtp == decimalwc)
997 /* Round the integer digits. */
998 if (*(wtp - 1) == decimalwc)
999 --wtp;
1001 while (--wtp >= wstartp && *wtp == L_('9'))
1002 *wtp = L_('0');
1004 if (wtp >= wstartp)
1005 /* Round up. */
1006 (*wtp)++;
1007 else
1008 /* It is more critical. All digits were 9's. */
1010 if (type != 'f')
1012 *wstartp = '1';
1013 exponent += expsign == 0 ? 1 : -1;
1015 /* The above exponent adjustment could lead to 1.0e-00,
1016 e.g. for 0.999999999. Make sure exponent 0 always
1017 uses + sign. */
1018 if (exponent == 0)
1019 expsign = 0;
1021 else if (intdig_no == dig_max)
1023 /* This is the case where for type %g the number fits
1024 really in the range for %f output but after rounding
1025 the number of digits is too big. */
1026 *--wstartp = decimalwc;
1027 *--wstartp = L_('1');
1029 if (info->alt || fracdig_no > 0)
1031 /* Overwrite the old radix character. */
1032 wstartp[intdig_no + 2] = L_('0');
1033 ++fracdig_no;
1036 fracdig_no += intdig_no;
1037 intdig_no = 1;
1038 fracdig_max = intdig_max - intdig_no;
1039 ++exponent;
1040 /* Now we must print the exponent. */
1041 type = isupper (info->spec) ? 'E' : 'e';
1043 else
1045 /* We can simply add another another digit before the
1046 radix. */
1047 *--wstartp = L_('1');
1048 ++intdig_no;
1051 /* While rounding the number of digits can change.
1052 If the number now exceeds the limits remove some
1053 fractional digits. */
1054 if (intdig_no + fracdig_no > dig_max)
1056 wcp -= intdig_no + fracdig_no - dig_max;
1057 fracdig_no -= intdig_no + fracdig_no - dig_max;
1062 #endif
1064 /* Now remove unnecessary '0' at the end of the string. */
1065 while (fracdig_no > fracdig_min + added_zeros && *(wcp - 1) == L_('0'))
1067 --wcp;
1068 --fracdig_no;
1070 /* If we eliminate all fractional digits we perhaps also can remove
1071 the radix character. */
1072 if (fracdig_no == 0 && !info->alt && *(wcp - 1) == decimalwc)
1073 --wcp;
1075 if (grouping)
1077 /* Rounding might have changed the number of groups. We allocated
1078 enough memory but we need here the correct number of groups. */
1079 if (intdig_no != intdig_max)
1080 ngroups = guess_grouping (intdig_no, grouping);
1082 /* Add in separator characters, overwriting the same buffer. */
1083 wcp = group_number (wstartp, wcp, intdig_no, grouping, thousands_sepwc,
1084 ngroups);
1087 /* Write the exponent if it is needed. */
1088 if (type != 'f')
1090 if (__builtin_expect (expsign != 0 && exponent == 4 && spec == 'g', 0))
1092 /* This is another special case. The exponent of the number is
1093 really smaller than -4, which requires the 'e'/'E' format.
1094 But after rounding the number has an exponent of -4. */
1095 assert (wcp >= wstartp + 1);
1096 assert (wstartp[0] == L_('1'));
1097 memcpy (wstartp, L_("0.0001"), 6 * sizeof (wchar_t));
1098 wstartp[1] = decimalwc;
1099 if (wcp >= wstartp + 2)
1101 size_t cnt;
1102 for (cnt = 0; cnt < wcp - (wstartp + 2); cnt++)
1103 wstartp[6 + cnt] = L_('0');
1104 wcp += 4;
1106 else
1107 wcp += 5;
1109 else
1111 *wcp++ = (wchar_t) type;
1112 *wcp++ = expsign ? L_('-') : L_('+');
1114 /* Find the magnitude of the exponent. */
1115 expscale = 10;
1116 while (expscale <= exponent)
1117 expscale *= 10;
1119 if (exponent < 10)
1120 /* Exponent always has at least two digits. */
1121 *wcp++ = L_('0');
1122 else
1125 expscale /= 10;
1126 *wcp++ = L_('0') + (exponent / expscale);
1127 exponent %= expscale;
1129 while (expscale > 10);
1130 *wcp++ = L_('0') + exponent;
1134 /* Compute number of characters which must be filled with the padding
1135 character. */
1136 if (is_neg || info->showsign || info->space)
1137 --width;
1138 width -= wcp - wstartp;
1140 if (!info->left && info->pad != '0' && width > 0)
1141 PADN (info->pad, width);
1143 if (is_neg)
1144 outchar ('-');
1145 else if (info->showsign)
1146 outchar ('+');
1147 else if (info->space)
1148 outchar (' ');
1150 if (!info->left && info->pad == '0' && width > 0)
1151 PADN ('0', width);
1154 char *buffer = NULL;
1155 char *buffer_end __attribute__((__unused__)) = NULL;
1156 char *cp = NULL;
1157 char *tmpptr;
1159 if (! wide)
1161 /* Create the single byte string. */
1162 size_t decimal_len;
1163 size_t thousands_sep_len;
1164 wchar_t *copywc;
1165 #ifdef USE_I18N_NUMBER_H
1166 size_t factor = (info->i18n
1167 ? nl_langinfo_wc (_NL_CTYPE_MB_CUR_MAX)
1168 : 1);
1169 #else
1170 size_t factor = 1;
1171 #endif
1173 decimal_len = strlen (decimal);
1175 if (thousands_sep == NULL)
1176 thousands_sep_len = 0;
1177 else
1178 thousands_sep_len = strlen (thousands_sep);
1180 size_t nbuffer = (2 + chars_needed * factor + decimal_len
1181 + ngroups * thousands_sep_len);
1182 if (__builtin_expect (buffer_malloced, 0))
1184 buffer = (char *) malloc (nbuffer);
1185 if (buffer == NULL)
1187 /* Signal an error to the caller. */
1188 free (wbuffer);
1189 return -1;
1192 else
1193 buffer = (char *) alloca (nbuffer);
1194 buffer_end = buffer + nbuffer;
1196 /* Now copy the wide character string. Since the character
1197 (except for the decimal point and thousands separator) must
1198 be coming from the ASCII range we can esily convert the
1199 string without mapping tables. */
1200 for (cp = buffer, copywc = wstartp; copywc < wcp; ++copywc)
1201 if (*copywc == decimalwc)
1202 memcpy (cp, decimal, decimal_len), cp += decimal_len;
1203 else if (*copywc == thousands_sepwc)
1204 memcpy (cp, thousands_sep, thousands_sep_len), cp += thousands_sep_len;
1205 else
1206 *cp++ = (char) *copywc;
1209 tmpptr = buffer;
1210 #if USE_I18N_NUMBER_H
1211 if (__builtin_expect (info->i18n, 0))
1213 tmpptr = _i18n_number_rewrite (tmpptr, cp, buffer_end);
1214 cp = buffer_end;
1215 assert ((uintptr_t) buffer <= (uintptr_t) tmpptr);
1216 assert ((uintptr_t) tmpptr < (uintptr_t) buffer_end);
1218 #endif
1220 PRINT (tmpptr, wstartp, wide ? wcp - wstartp : cp - tmpptr);
1222 /* Free the memory if necessary. */
1223 if (__builtin_expect (buffer_malloced, 0))
1225 free (buffer);
1226 free (wbuffer);
1230 if (info->left && width > 0)
1231 PADN (info->pad, width);
1233 return done;
1236 /* Return the number of extra grouping characters that will be inserted
1237 into a number with INTDIG_MAX integer digits. */
1239 static unsigned int
1240 guess_grouping (unsigned int intdig_max, const char *grouping)
1242 unsigned int groups;
1244 /* We treat all negative values like CHAR_MAX. */
1246 if (*grouping == CHAR_MAX || *grouping <= 0)
1247 /* No grouping should be done. */
1248 return 0;
1250 groups = 0;
1251 while (intdig_max > (unsigned int) *grouping)
1253 ++groups;
1254 intdig_max -= *grouping++;
1256 if (*grouping == CHAR_MAX
1257 #if CHAR_MIN < 0
1258 || *grouping < 0
1259 #endif
1261 /* No more grouping should be done. */
1262 break;
1263 else if (*grouping == 0)
1265 /* Same grouping repeats. */
1266 groups += (intdig_max - 1) / grouping[-1];
1267 break;
1271 return groups;
1274 /* Group the INTDIG_NO integer digits of the number in [BUF,BUFEND).
1275 There is guaranteed enough space past BUFEND to extend it.
1276 Return the new end of buffer. */
1278 static wchar_t *
1279 group_number (wchar_t *buf, wchar_t *bufend, unsigned int intdig_no,
1280 const char *grouping, wchar_t thousands_sep, int ngroups)
1282 wchar_t *p;
1284 if (ngroups == 0)
1285 return bufend;
1287 /* Move the fractional part down. */
1288 memmove (buf + intdig_no + ngroups, buf + intdig_no,
1289 (bufend - (buf + intdig_no)) * sizeof (wchar_t));
1291 p = buf + intdig_no + ngroups - 1;
1294 unsigned int len = *grouping++;
1296 *p-- = buf[--intdig_no];
1297 while (--len > 0);
1298 *p-- = thousands_sep;
1300 if (*grouping == CHAR_MAX
1301 #if CHAR_MIN < 0
1302 || *grouping < 0
1303 #endif
1305 /* No more grouping should be done. */
1306 break;
1307 else if (*grouping == 0)
1308 /* Same grouping repeats. */
1309 --grouping;
1310 } while (intdig_no > (unsigned int) *grouping);
1312 /* Copy the remaining ungrouped digits. */
1314 *p-- = buf[--intdig_no];
1315 while (p > buf);
1317 return bufend + ngroups;