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