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