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