stdio-common/printf_fp.c

   1 /* Floating point output for `printf'.
   2    Copyright (C) 1995-2015 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 \f
  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 \f
 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
 136 extern mp_size_t __mpn_extract_double (mp_ptr res_ptr, mp_size_t size,
 137                                        int *expt, int *is_neg,
 138                                        double value);
 139 extern 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);
 142 extern unsigned int __guess_grouping (unsigned int intdig_max,
 143                                       const char *grouping);
 144
 145
 146 static 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
 151 struct 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
 167 static wchar_t
 168 hack_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
 211 int
 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       int res;
 336       if (__isnanl (fpnum.ldbl))
 337         {
 338           is_neg = signbit (fpnum.ldbl);
 339           if (isupper (info->spec))
 340             {
 341               special = "NAN";
 342               wspecial = L"NAN";
 343             }
 344             else
 345               {
 346                 special = "nan";
 347                 wspecial = L"nan";
 348               }
 349         }
 350       else if ((res = __isinfl (fpnum.ldbl)))
 351         {
 352           is_neg = res < 0;
 353           if (isupper (info->spec))
 354             {
 355               special = "INF";
 356               wspecial = L"INF";
 357             }
 358           else
 359             {
 360               special = "inf";
 361               wspecial = L"inf";
 362             }
 363         }
 364       else
 365         {
 366           p.fracsize = __mpn_extract_long_double (fp_input,
 367                                                 (sizeof (fp_input) /
 368                                                  sizeof (fp_input[0])),
 369                                                 &p.exponent, &is_neg,
 370                                                 fpnum.ldbl);
 371           to_shift = 1 + p.fracsize * BITS_PER_MP_LIMB - LDBL_MANT_DIG;
 372         }
 373     }
 374   else
 375 #endif  /* no long double */
 376     {
 377       fpnum.dbl = *(const double *) args[0];
 378
 379       /* Check for special values: not a number or infinity.  */
 380       int res;
 381       if (__isnan (fpnum.dbl))
 382         {
 383           union ieee754_double u = { .d = fpnum.dbl };
 384           is_neg = u.ieee.negative != 0;
 385           if (isupper (info->spec))
 386             {
 387               special = "NAN";
 388               wspecial = L"NAN";
 389             }
 390           else
 391             {
 392               special = "nan";
 393               wspecial = L"nan";
 394             }
 395         }
 396       else if ((res = __isinf (fpnum.dbl)))
 397         {
 398           is_neg = res < 0;
 399           if (isupper (info->spec))
 400             {
 401               special = "INF";
 402               wspecial = L"INF";
 403             }
 404           else
 405             {
 406               special = "inf";
 407               wspecial = L"inf";
 408             }
 409         }
 410       else
 411         {
 412           p.fracsize = __mpn_extract_double (fp_input,
 413                                            (sizeof (fp_input)
 414                                             / sizeof (fp_input[0])),
 415                                            &p.exponent, &is_neg, fpnum.dbl);
 416           to_shift = 1 + p.fracsize * BITS_PER_MP_LIMB - DBL_MANT_DIG;
 417         }
 418     }
 419
 420   if (special)
 421     {
 422       int width = info->width;
 423
 424       if (is_neg || info->showsign || info->space)
 425         --width;
 426       width -= 3;
 427
 428       if (!info->left && width > 0)
 429         PADN (' ', width);
 430
 431       if (is_neg)
 432         outchar ('-');
 433       else if (info->showsign)
 434         outchar ('+');
 435       else if (info->space)
 436         outchar (' ');
 437
 438       PRINT (special, wspecial, 3);
 439
 440       if (info->left && width > 0)
 441         PADN (' ', width);
 442
 443       return done;
 444     }
 445
 446
 447   /* We need three multiprecision variables.  Now that we have the p.exponent
 448      of the number we can allocate the needed memory.  It would be more
 449      efficient to use variables of the fixed maximum size but because this
 450      would be really big it could lead to memory problems.  */
 451   {
 452     mp_size_t bignum_size = ((ABS (p.exponent) + BITS_PER_MP_LIMB - 1)
 453                              / BITS_PER_MP_LIMB
 454                              + (LDBL_MANT_DIG / BITS_PER_MP_LIMB > 2 ? 8 : 4))
 455                             * sizeof (mp_limb_t);
 456     p.frac = (mp_limb_t *) alloca (bignum_size);
 457     p.tmp = (mp_limb_t *) alloca (bignum_size);
 458     p.scale = (mp_limb_t *) alloca (bignum_size);
 459   }
 460
 461   /* We now have to distinguish between numbers with positive and negative
 462      exponents because the method used for the one is not applicable/efficient
 463      for the other.  */
 464   p.scalesize = 0;
 465   if (p.exponent > 2)
 466     {
 467       /* |FP| >= 8.0.  */
 468       int scaleexpo = 0;
 469       int explog = LDBL_MAX_10_EXP_LOG;
 470       int exp10 = 0;
 471       const struct mp_power *powers = &_fpioconst_pow10[explog + 1];
 472       int cnt_h, cnt_l, i;
 473
 474       if ((p.exponent + to_shift) % BITS_PER_MP_LIMB == 0)
 475         {
 476           MPN_COPY_DECR (p.frac + (p.exponent + to_shift) / BITS_PER_MP_LIMB,
 477                          fp_input, p.fracsize);
 478           p.fracsize += (p.exponent + to_shift) / BITS_PER_MP_LIMB;
 479         }
 480       else
 481         {
 482           cy = __mpn_lshift (p.frac +
 483                              (p.exponent + to_shift) / BITS_PER_MP_LIMB,
 484                              fp_input, p.fracsize,
 485                              (p.exponent + to_shift) % BITS_PER_MP_LIMB);
 486           p.fracsize += (p.exponent + to_shift) / BITS_PER_MP_LIMB;
 487           if (cy)
 488             p.frac[p.fracsize++] = cy;
 489         }
 490       MPN_ZERO (p.frac, (p.exponent + to_shift) / BITS_PER_MP_LIMB);
 491
 492       assert (powers > &_fpioconst_pow10[0]);
 493       do
 494         {
 495           --powers;
 496
 497           /* The number of the product of two binary numbers with n and m
 498              bits respectively has m+n or m+n-1 bits.   */
 499           if (p.exponent >= scaleexpo + powers->p_expo - 1)
 500             {
 501               if (p.scalesize == 0)
 502                 {
 503 #ifndef __NO_LONG_DOUBLE_MATH
 504                   if (LDBL_MANT_DIG > _FPIO_CONST_OFFSET * BITS_PER_MP_LIMB
 505                       && info->is_long_double)
 506                     {
 507 #define _FPIO_CONST_SHIFT \
 508   (((LDBL_MANT_DIG + BITS_PER_MP_LIMB - 1) / BITS_PER_MP_LIMB) \
 509    - _FPIO_CONST_OFFSET)
 510                       /* 64bit const offset is not enough for
 511                          IEEE quad long double.  */
 512                       p.tmpsize = powers->arraysize + _FPIO_CONST_SHIFT;
 513                       memcpy (p.tmp + _FPIO_CONST_SHIFT,
 514                               &__tens[powers->arrayoff],
 515                               p.tmpsize * sizeof (mp_limb_t));
 516                       MPN_ZERO (p.tmp, _FPIO_CONST_SHIFT);
 517                       /* Adjust p.exponent, as scaleexpo will be this much
 518                          bigger too.  */
 519                       p.exponent += _FPIO_CONST_SHIFT * BITS_PER_MP_LIMB;
 520                     }
 521                   else
 522 #endif
 523                     {
 524                       p.tmpsize = powers->arraysize;
 525                       memcpy (p.tmp, &__tens[powers->arrayoff],
 526                               p.tmpsize * sizeof (mp_limb_t));
 527                     }
 528                 }
 529               else
 530                 {
 531                   cy = __mpn_mul (p.tmp, p.scale, p.scalesize,
 532                                   &__tens[powers->arrayoff
 533                                          + _FPIO_CONST_OFFSET],
 534                                   powers->arraysize - _FPIO_CONST_OFFSET);
 535                   p.tmpsize = p.scalesize +
 536                     powers->arraysize - _FPIO_CONST_OFFSET;
 537                   if (cy == 0)
 538                     --p.tmpsize;
 539                 }
 540
 541               if (MPN_GE (p.frac, p.tmp))
 542                 {
 543                   int cnt;
 544                   MPN_ASSIGN (p.scale, p.tmp);
 545                   count_leading_zeros (cnt, p.scale[p.scalesize - 1]);
 546                   scaleexpo = (p.scalesize - 2) * BITS_PER_MP_LIMB - cnt - 1;
 547                   exp10 |= 1 << explog;
 548                 }
 549             }
 550           --explog;
 551         }
 552       while (powers > &_fpioconst_pow10[0]);
 553       p.exponent = exp10;
 554
 555       /* Optimize number representations.  We want to represent the numbers
 556          with the lowest number of bytes possible without losing any
 557          bytes. Also the highest bit in the scaling factor has to be set
 558          (this is a requirement of the MPN division routines).  */
 559       if (p.scalesize > 0)
 560         {
 561           /* Determine minimum number of zero bits at the end of
 562              both numbers.  */
 563           for (i = 0; p.scale[i] == 0 && p.frac[i] == 0; i++)
 564             ;
 565
 566           /* Determine number of bits the scaling factor is misplaced.  */
 567           count_leading_zeros (cnt_h, p.scale[p.scalesize - 1]);
 568
 569           if (cnt_h == 0)
 570             {
 571               /* The highest bit of the scaling factor is already set.  So
 572                  we only have to remove the trailing empty limbs.  */
 573               if (i > 0)
 574                 {
 575                   MPN_COPY_INCR (p.scale, p.scale + i, p.scalesize - i);
 576                   p.scalesize -= i;
 577                   MPN_COPY_INCR (p.frac, p.frac + i, p.fracsize - i);
 578                   p.fracsize -= i;
 579                 }
 580             }
 581           else
 582             {
 583               if (p.scale[i] != 0)
 584                 {
 585                   count_trailing_zeros (cnt_l, p.scale[i]);
 586                   if (p.frac[i] != 0)
 587                     {
 588                       int cnt_l2;
 589                       count_trailing_zeros (cnt_l2, p.frac[i]);
 590                       if (cnt_l2 < cnt_l)
 591                         cnt_l = cnt_l2;
 592                     }
 593                 }
 594               else
 595                 count_trailing_zeros (cnt_l, p.frac[i]);
 596
 597               /* Now shift the numbers to their optimal position.  */
 598               if (i == 0 && BITS_PER_MP_LIMB - cnt_h > cnt_l)
 599                 {
 600                   /* We cannot save any memory.  So just roll both numbers
 601                      so that the scaling factor has its highest bit set.  */
 602
 603                   (void) __mpn_lshift (p.scale, p.scale, p.scalesize, cnt_h);
 604                   cy = __mpn_lshift (p.frac, p.frac, p.fracsize, cnt_h);
 605                   if (cy != 0)
 606                     p.frac[p.fracsize++] = cy;
 607                 }
 608               else if (BITS_PER_MP_LIMB - cnt_h <= cnt_l)
 609                 {
 610                   /* We can save memory by removing the trailing zero limbs
 611                      and by packing the non-zero limbs which gain another
 612                      free one. */
 613
 614                   (void) __mpn_rshift (p.scale, p.scale + i, p.scalesize - i,
 615                                        BITS_PER_MP_LIMB - cnt_h);
 616                   p.scalesize -= i + 1;
 617                   (void) __mpn_rshift (p.frac, p.frac + i, p.fracsize - i,
 618                                        BITS_PER_MP_LIMB - cnt_h);
 619                   p.fracsize -= p.frac[p.fracsize - i - 1] == 0 ? i + 1 : i;
 620                 }
 621               else
 622                 {
 623                   /* We can only save the memory of the limbs which are zero.
 624                      The non-zero parts occupy the same number of limbs.  */
 625
 626                   (void) __mpn_rshift (p.scale, p.scale + (i - 1),
 627                                        p.scalesize - (i - 1),
 628                                        BITS_PER_MP_LIMB - cnt_h);
 629                   p.scalesize -= i;
 630                   (void) __mpn_rshift (p.frac, p.frac + (i - 1),
 631                                        p.fracsize - (i - 1),
 632                                        BITS_PER_MP_LIMB - cnt_h);
 633                   p.fracsize -=
 634                     p.frac[p.fracsize - (i - 1) - 1] == 0 ? i : i - 1;
 635                 }
 636             }
 637         }
 638     }
 639   else if (p.exponent < 0)
 640     {
 641       /* |FP| < 1.0.  */
 642       int exp10 = 0;
 643       int explog = LDBL_MAX_10_EXP_LOG;
 644       const struct mp_power *powers = &_fpioconst_pow10[explog + 1];
 645
 646       /* Now shift the input value to its right place.  */
 647       cy = __mpn_lshift (p.frac, fp_input, p.fracsize, to_shift);
 648       p.frac[p.fracsize++] = cy;
 649       assert (cy == 1 || (p.frac[p.fracsize - 2] == 0 && p.frac[0] == 0));
 650
 651       p.expsign = 1;
 652       p.exponent = -p.exponent;
 653
 654       assert (powers != &_fpioconst_pow10[0]);
 655       do
 656         {
 657           --powers;
 658
 659           if (p.exponent >= powers->m_expo)
 660             {
 661               int i, incr, cnt_h, cnt_l;
 662               mp_limb_t topval[2];
 663
 664               /* The __mpn_mul function expects the first argument to be
 665                  bigger than the second.  */
 666               if (p.fracsize < powers->arraysize - _FPIO_CONST_OFFSET)
 667                 cy = __mpn_mul (p.tmp, &__tens[powers->arrayoff
 668                                             + _FPIO_CONST_OFFSET],
 669                                 powers->arraysize - _FPIO_CONST_OFFSET,
 670                                 p.frac, p.fracsize);
 671               else
 672                 cy = __mpn_mul (p.tmp, p.frac, p.fracsize,
 673                                 &__tens[powers->arrayoff + _FPIO_CONST_OFFSET],
 674                                 powers->arraysize - _FPIO_CONST_OFFSET);
 675               p.tmpsize = p.fracsize + powers->arraysize - _FPIO_CONST_OFFSET;
 676               if (cy == 0)
 677                 --p.tmpsize;
 678
 679               count_leading_zeros (cnt_h, p.tmp[p.tmpsize - 1]);
 680               incr = (p.tmpsize - p.fracsize) * BITS_PER_MP_LIMB
 681                      + BITS_PER_MP_LIMB - 1 - cnt_h;
 682
 683               assert (incr <= powers->p_expo);
 684
 685               /* If we increased the p.exponent by exactly 3 we have to test
 686                  for overflow.  This is done by comparing with 10 shifted
 687                  to the right position.  */
 688               if (incr == p.exponent + 3)
 689                 {
 690                   if (cnt_h <= BITS_PER_MP_LIMB - 4)
 691                     {
 692                       topval[0] = 0;
 693                       topval[1]
 694                         = ((mp_limb_t) 10) << (BITS_PER_MP_LIMB - 4 - cnt_h);
 695                     }
 696                   else
 697                     {
 698                       topval[0] = ((mp_limb_t) 10) << (BITS_PER_MP_LIMB - 4);
 699                       topval[1] = 0;
 700                       (void) __mpn_lshift (topval, topval, 2,
 701                                            BITS_PER_MP_LIMB - cnt_h);
 702                     }
 703                 }
 704
 705               /* We have to be careful when multiplying the last factor.
 706                  If the result is greater than 1.0 be have to test it
 707                  against 10.0.  If it is greater or equal to 10.0 the
 708                  multiplication was not valid.  This is because we cannot
 709                  determine the number of bits in the result in advance.  */
 710               if (incr < p.exponent + 3
 711                   || (incr == p.exponent + 3 &&
 712                       (p.tmp[p.tmpsize - 1] < topval[1]
 713                        || (p.tmp[p.tmpsize - 1] == topval[1]
 714                            && p.tmp[p.tmpsize - 2] < topval[0]))))
 715                 {
 716                   /* The factor is right.  Adapt binary and decimal
 717                      exponents.  */
 718                   p.exponent -= incr;
 719                   exp10 |= 1 << explog;
 720
 721                   /* If this factor yields a number greater or equal to
 722                      1.0, we must not shift the non-fractional digits down. */
 723                   if (p.exponent < 0)
 724                     cnt_h += -p.exponent;
 725
 726                   /* Now we optimize the number representation.  */
 727                   for (i = 0; p.tmp[i] == 0; ++i);
 728                   if (cnt_h == BITS_PER_MP_LIMB - 1)
 729                     {
 730                       MPN_COPY (p.frac, p.tmp + i, p.tmpsize - i);
 731                       p.fracsize = p.tmpsize - i;
 732                     }
 733                   else
 734                     {
 735                       count_trailing_zeros (cnt_l, p.tmp[i]);
 736
 737                       /* Now shift the numbers to their optimal position.  */
 738                       if (i == 0 && BITS_PER_MP_LIMB - 1 - cnt_h > cnt_l)
 739                         {
 740                           /* We cannot save any memory.  Just roll the
 741                              number so that the leading digit is in a
 742                              separate limb.  */
 743
 744                           cy = __mpn_lshift (p.frac, p.tmp, p.tmpsize,
 745                             cnt_h + 1);
 746                           p.fracsize = p.tmpsize + 1;
 747                           p.frac[p.fracsize - 1] = cy;
 748                         }
 749                       else if (BITS_PER_MP_LIMB - 1 - cnt_h <= cnt_l)
 750                         {
 751                           (void) __mpn_rshift (p.frac, p.tmp + i, p.tmpsize - i,
 752                                                BITS_PER_MP_LIMB - 1 - cnt_h);
 753                           p.fracsize = p.tmpsize - i;
 754                         }
 755                       else
 756                         {
 757                           /* We can only save the memory of the limbs which
 758                              are zero.  The non-zero parts occupy the same
 759                              number of limbs.  */
 760
 761                           (void) __mpn_rshift (p.frac, p.tmp + (i - 1),
 762                                                p.tmpsize - (i - 1),
 763                                                BITS_PER_MP_LIMB - 1 - cnt_h);
 764                           p.fracsize = p.tmpsize - (i - 1);
 765                         }
 766                     }
 767                 }
 768             }
 769           --explog;
 770         }
 771       while (powers != &_fpioconst_pow10[1] && p.exponent > 0);
 772       /* All factors but 10^-1 are tested now.  */
 773       if (p.exponent > 0)
 774         {
 775           int cnt_l;
 776
 777           cy = __mpn_mul_1 (p.tmp, p.frac, p.fracsize, 10);
 778           p.tmpsize = p.fracsize;
 779           assert (cy == 0 || p.tmp[p.tmpsize - 1] < 20);
 780
 781           count_trailing_zeros (cnt_l, p.tmp[0]);
 782           if (cnt_l < MIN (4, p.exponent))
 783             {
 784               cy = __mpn_lshift (p.frac, p.tmp, p.tmpsize,
 785                                  BITS_PER_MP_LIMB - MIN (4, p.exponent));
 786               if (cy != 0)
 787                 p.frac[p.tmpsize++] = cy;
 788             }
 789           else
 790             (void) __mpn_rshift (p.frac, p.tmp, p.tmpsize, MIN (4, p.exponent));
 791           p.fracsize = p.tmpsize;
 792           exp10 |= 1;
 793           assert (p.frac[p.fracsize - 1] < 10);
 794         }
 795       p.exponent = exp10;
 796     }
 797   else
 798     {
 799       /* This is a special case.  We don't need a factor because the
 800          numbers are in the range of 1.0 <= |fp| < 8.0.  We simply
 801          shift it to the right place and divide it by 1.0 to get the
 802          leading digit.  (Of course this division is not really made.)  */
 803       assert (0 <= p.exponent && p.exponent < 3 &&
 804               p.exponent + to_shift < BITS_PER_MP_LIMB);
 805
 806       /* Now shift the input value to its right place.  */
 807       cy = __mpn_lshift (p.frac, fp_input, p.fracsize, (p.exponent + to_shift));
 808       p.frac[p.fracsize++] = cy;
 809       p.exponent = 0;
 810     }
 811
 812   {
 813     int width = info->width;
 814     wchar_t *wstartp, *wcp;
 815     size_t chars_needed;
 816     int expscale;
 817     int intdig_max, intdig_no = 0;
 818     int fracdig_min;
 819     int fracdig_max;
 820     int dig_max;
 821     int significant;
 822     int ngroups = 0;
 823     char spec = _tolower (info->spec);
 824
 825     if (spec == 'e')
 826       {
 827         p.type = info->spec;
 828         intdig_max = 1;
 829         fracdig_min = fracdig_max = info->prec < 0 ? 6 : info->prec;
 830         chars_needed = 1 + 1 + (size_t) fracdig_max + 1 + 1 + 4;
 831         /*             d   .     ddd         e   +-  ddd  */
 832         dig_max = INT_MAX;              /* Unlimited.  */
 833         significant = 1;                /* Does not matter here.  */
 834       }
 835     else if (spec == 'f')
 836       {
 837         p.type = 'f';
 838         fracdig_min = fracdig_max = info->prec < 0 ? 6 : info->prec;
 839         dig_max = INT_MAX;              /* Unlimited.  */
 840         significant = 1;                /* Does not matter here.  */
 841         if (p.expsign == 0)
 842           {
 843             intdig_max = p.exponent + 1;
 844             /* This can be really big!  */  /* XXX Maybe malloc if too big? */
 845             chars_needed = (size_t) p.exponent + 1 + 1 + (size_t) fracdig_max;
 846           }
 847         else
 848           {
 849             intdig_max = 1;
 850             chars_needed = 1 + 1 + (size_t) fracdig_max;
 851           }
 852       }
 853     else
 854       {
 855         dig_max = info->prec < 0 ? 6 : (info->prec == 0 ? 1 : info->prec);
 856         if ((p.expsign == 0 && p.exponent >= dig_max)
 857             || (p.expsign != 0 && p.exponent > 4))
 858           {
 859             if ('g' - 'G' == 'e' - 'E')
 860               p.type = 'E' + (info->spec - 'G');
 861             else
 862               p.type = isupper (info->spec) ? 'E' : 'e';
 863             fracdig_max = dig_max - 1;
 864             intdig_max = 1;
 865             chars_needed = 1 + 1 + (size_t) fracdig_max + 1 + 1 + 4;
 866           }
 867         else
 868           {
 869             p.type = 'f';
 870             intdig_max = p.expsign == 0 ? p.exponent + 1 : 0;
 871             fracdig_max = dig_max - intdig_max;
 872             /* We need space for the significant digits and perhaps
 873                for leading zeros when < 1.0.  The number of leading
 874                zeros can be as many as would be required for
 875                exponential notation with a negative two-digit
 876                p.exponent, which is 4.  */
 877             chars_needed = (size_t) dig_max + 1 + 4;
 878           }
 879         fracdig_min = info->alt ? fracdig_max : 0;
 880         significant = 0;                /* We count significant digits.  */
 881       }
 882
 883     if (grouping)
 884       {
 885         /* Guess the number of groups we will make, and thus how
 886            many spaces we need for separator characters.  */
 887         ngroups = __guess_grouping (intdig_max, grouping);
 888         /* Allocate one more character in case rounding increases the
 889            number of groups.  */
 890         chars_needed += ngroups + 1;
 891       }
 892
 893     /* Allocate buffer for output.  We need two more because while rounding
 894        it is possible that we need two more characters in front of all the
 895        other output.  If the amount of memory we have to allocate is too
 896        large use `malloc' instead of `alloca'.  */
 897     if (__builtin_expect (chars_needed >= (size_t) -1 / sizeof (wchar_t) - 2
 898                           || chars_needed < fracdig_max, 0))
 899       {
 900         /* Some overflow occurred.  */
 901         __set_errno (ERANGE);
 902         return -1;
 903       }
 904     size_t wbuffer_to_alloc = (2 + chars_needed) * sizeof (wchar_t);
 905     buffer_malloced = ! __libc_use_alloca (wbuffer_to_alloc);
 906     if (__builtin_expect (buffer_malloced, 0))
 907       {
 908         wbuffer = (wchar_t *) malloc (wbuffer_to_alloc);
 909         if (wbuffer == NULL)
 910           /* Signal an error to the caller.  */
 911           return -1;
 912       }
 913     else
 914       wbuffer = (wchar_t *) alloca (wbuffer_to_alloc);
 915     wcp = wstartp = wbuffer + 2;        /* Let room for rounding.  */
 916
 917     /* Do the real work: put digits in allocated buffer.  */
 918     if (p.expsign == 0 || p.type != 'f')
 919       {
 920         assert (p.expsign == 0 || intdig_max == 1);
 921         while (intdig_no < intdig_max)
 922           {
 923             ++intdig_no;
 924             *wcp++ = hack_digit (&p);
 925           }
 926         significant = 1;
 927         if (info->alt
 928             || fracdig_min > 0
 929             || (fracdig_max > 0 && (p.fracsize > 1 || p.frac[0] != 0)))
 930           *wcp++ = decimalwc;
 931       }
 932     else
 933       {
 934         /* |fp| < 1.0 and the selected p.type is 'f', so put "0."
 935            in the buffer.  */
 936         *wcp++ = L'0';
 937         --p.exponent;
 938         *wcp++ = decimalwc;
 939       }
 940
 941     /* Generate the needed number of fractional digits.  */
 942     int fracdig_no = 0;
 943     int added_zeros = 0;
 944     while (fracdig_no < fracdig_min + added_zeros
 945            || (fracdig_no < fracdig_max && (p.fracsize > 1 || p.frac[0] != 0)))
 946       {
 947         ++fracdig_no;
 948         *wcp = hack_digit (&p);
 949         if (*wcp++ != L'0')
 950           significant = 1;
 951         else if (significant == 0)
 952           {
 953             ++fracdig_max;
 954             if (fracdig_min > 0)
 955               ++added_zeros;
 956           }
 957       }
 958
 959     /* Do rounding.  */
 960     wchar_t last_digit = wcp[-1] != decimalwc ? wcp[-1] : wcp[-2];
 961     wchar_t next_digit = hack_digit (&p);
 962     bool more_bits;
 963     if (next_digit != L'0' && next_digit != L'5')
 964       more_bits = true;
 965     else if (p.fracsize == 1 && p.frac[0] == 0)
 966       /* Rest of the number is zero.  */
 967       more_bits = false;
 968     else if (p.scalesize == 0)
 969       {
 970         /* Here we have to see whether all limbs are zero since no
 971            normalization happened.  */
 972         size_t lcnt = p.fracsize;
 973         while (lcnt >= 1 && p.frac[lcnt - 1] == 0)
 974           --lcnt;
 975         more_bits = lcnt > 0;
 976       }
 977     else
 978       more_bits = true;
 979     int rounding_mode = get_rounding_mode ();
 980     if (round_away (is_neg, (last_digit - L'0') & 1, next_digit >= L'5',
 981                     more_bits, rounding_mode))
 982       {
 983         wchar_t *wtp = wcp;
 984
 985         if (fracdig_no > 0)
 986           {
 987             /* Process fractional digits.  Terminate if not rounded or
 988                radix character is reached.  */
 989             int removed = 0;
 990             while (*--wtp != decimalwc && *wtp == L'9')
 991               {
 992                 *wtp = L'0';
 993                 ++removed;
 994               }
 995             if (removed == fracdig_min && added_zeros > 0)
 996               --added_zeros;
 997             if (*wtp != decimalwc)
 998               /* Round up.  */
 999               (*wtp)++;
1000             else if (__builtin_expect (spec == 'g' && p.type == 'f' && info->alt
1001                                        && wtp == wstartp + 1
1002                                        && wstartp[0] == L'0',
1003                                        0))
1004               /* This is a special case: the rounded number is 1.0,
1005                  the format is 'g' or 'G', and the alternative format
1006                  is selected.  This means the result must be "1.".  */
1007               --added_zeros;
1008           }
1009
1010         if (fracdig_no == 0 || *wtp == decimalwc)
1011           {
1012             /* Round the integer digits.  */
1013             if (*(wtp - 1) == decimalwc)
1014               --wtp;
1015
1016             while (--wtp >= wstartp && *wtp == L'9')
1017               *wtp = L'0';
1018
1019             if (wtp >= wstartp)
1020               /* Round up.  */
1021               (*wtp)++;
1022             else
1023               /* It is more critical.  All digits were 9's.  */
1024               {
1025                 if (p.type != 'f')
1026                   {
1027                     *wstartp = '1';
1028                     p.exponent += p.expsign == 0 ? 1 : -1;
1029
1030                     /* The above p.exponent adjustment could lead to 1.0e-00,
1031                        e.g. for 0.999999999.  Make sure p.exponent 0 always
1032                        uses + sign.  */
1033                     if (p.exponent == 0)
1034                       p.expsign = 0;
1035                   }
1036                 else if (intdig_no == dig_max)
1037                   {
1038                     /* This is the case where for p.type %g the number fits
1039                        really in the range for %f output but after rounding
1040                        the number of digits is too big.  */
1041                     *--wstartp = decimalwc;
1042                     *--wstartp = L'1';
1043
1044                     if (info->alt || fracdig_no > 0)
1045                       {
1046                         /* Overwrite the old radix character.  */
1047                         wstartp[intdig_no + 2] = L'0';
1048                         ++fracdig_no;
1049                       }
1050
1051                     fracdig_no += intdig_no;
1052                     intdig_no = 1;
1053                     fracdig_max = intdig_max - intdig_no;
1054                     ++p.exponent;
1055                     /* Now we must print the p.exponent.        */
1056                     p.type = isupper (info->spec) ? 'E' : 'e';
1057                   }
1058                 else
1059                   {
1060                     /* We can simply add another another digit before the
1061                        radix.  */
1062                     *--wstartp = L'1';
1063                     ++intdig_no;
1064                   }
1065
1066                 /* While rounding the number of digits can change.
1067                    If the number now exceeds the limits remove some
1068                    fractional digits.  */
1069                 if (intdig_no + fracdig_no > dig_max)
1070                   {
1071                     wcp -= intdig_no + fracdig_no - dig_max;
1072                     fracdig_no -= intdig_no + fracdig_no - dig_max;
1073                   }
1074               }
1075           }
1076       }
1077
1078     /* Now remove unnecessary '0' at the end of the string.  */
1079     while (fracdig_no > fracdig_min + added_zeros && *(wcp - 1) == L'0')
1080       {
1081         --wcp;
1082         --fracdig_no;
1083       }
1084     /* If we eliminate all fractional digits we perhaps also can remove
1085        the radix character.  */
1086     if (fracdig_no == 0 && !info->alt && *(wcp - 1) == decimalwc)
1087       --wcp;
1088
1089     if (grouping)
1090       {
1091         /* Rounding might have changed the number of groups.  We allocated
1092            enough memory but we need here the correct number of groups.  */
1093         if (intdig_no != intdig_max)
1094           ngroups = __guess_grouping (intdig_no, grouping);
1095
1096         /* Add in separator characters, overwriting the same buffer.  */
1097         wcp = group_number (wstartp, wcp, intdig_no, grouping, thousands_sepwc,
1098                             ngroups);
1099       }
1100
1101     /* Write the p.exponent if it is needed.  */
1102     if (p.type != 'f')
1103       {
1104         if (__glibc_unlikely (p.expsign != 0 && p.exponent == 4 && spec == 'g'))
1105           {
1106             /* This is another special case.  The p.exponent of the number is
1107                really smaller than -4, which requires the 'e'/'E' format.
1108                But after rounding the number has an p.exponent of -4.  */
1109             assert (wcp >= wstartp + 1);
1110             assert (wstartp[0] == L'1');
1111             __wmemcpy (wstartp, L"0.0001", 6);
1112             wstartp[1] = decimalwc;
1113             if (wcp >= wstartp + 2)
1114               {
1115                 __wmemset (wstartp + 6, L'0', wcp - (wstartp + 2));
1116                 wcp += 4;
1117               }
1118             else
1119               wcp += 5;
1120           }
1121         else
1122           {
1123             *wcp++ = (wchar_t) p.type;
1124             *wcp++ = p.expsign ? L'-' : L'+';
1125
1126             /* Find the magnitude of the p.exponent.    */
1127             expscale = 10;
1128             while (expscale <= p.exponent)
1129               expscale *= 10;
1130
1131             if (p.exponent < 10)
1132               /* Exponent always has at least two digits.  */
1133               *wcp++ = L'0';
1134             else
1135               do
1136                 {
1137                   expscale /= 10;
1138                   *wcp++ = L'0' + (p.exponent / expscale);
1139                   p.exponent %= expscale;
1140                 }
1141               while (expscale > 10);
1142             *wcp++ = L'0' + p.exponent;
1143           }
1144       }
1145
1146     /* Compute number of characters which must be filled with the padding
1147        character.  */
1148     if (is_neg || info->showsign || info->space)
1149       --width;
1150     width -= wcp - wstartp;
1151
1152     if (!info->left && info->pad != '0' && width > 0)
1153       PADN (info->pad, width);
1154
1155     if (is_neg)
1156       outchar ('-');
1157     else if (info->showsign)
1158       outchar ('+');
1159     else if (info->space)
1160       outchar (' ');
1161
1162     if (!info->left && info->pad == '0' && width > 0)
1163       PADN ('0', width);
1164
1165     {
1166       char *buffer = NULL;
1167       char *buffer_end = NULL;
1168       char *cp = NULL;
1169       char *tmpptr;
1170
1171       if (! wide)
1172         {
1173           /* Create the single byte string.  */
1174           size_t decimal_len;
1175           size_t thousands_sep_len;
1176           wchar_t *copywc;
1177           size_t factor = (info->i18n
1178                            ? _NL_CURRENT_WORD (LC_CTYPE, _NL_CTYPE_MB_CUR_MAX)
1179                            : 1);
1180
1181           decimal_len = strlen (decimal);
1182
1183           if (thousands_sep == NULL)
1184             thousands_sep_len = 0;
1185           else
1186             thousands_sep_len = strlen (thousands_sep);
1187
1188           size_t nbuffer = (2 + chars_needed * factor + decimal_len
1189                             + ngroups * thousands_sep_len);
1190           if (__glibc_unlikely (buffer_malloced))
1191             {
1192               buffer = (char *) malloc (nbuffer);
1193               if (buffer == NULL)
1194                 {
1195                   /* Signal an error to the caller.  */
1196                   free (wbuffer);
1197                   return -1;
1198                 }
1199             }
1200           else
1201             buffer = (char *) alloca (nbuffer);
1202           buffer_end = buffer + nbuffer;
1203
1204           /* Now copy the wide character string.  Since the character
1205              (except for the decimal point and thousands separator) must
1206              be coming from the ASCII range we can esily convert the
1207              string without mapping tables.  */
1208           for (cp = buffer, copywc = wstartp; copywc < wcp; ++copywc)
1209             if (*copywc == decimalwc)
1210               cp = (char *) __mempcpy (cp, decimal, decimal_len);
1211             else if (*copywc == thousands_sepwc)
1212               cp = (char *) __mempcpy (cp, thousands_sep, thousands_sep_len);
1213             else
1214               *cp++ = (char) *copywc;
1215         }
1216
1217       tmpptr = buffer;
1218       if (__glibc_unlikely (info->i18n))
1219         {
1220 #ifdef COMPILE_WPRINTF
1221           wstartp = _i18n_number_rewrite (wstartp, wcp,
1222                                           wbuffer + wbuffer_to_alloc);
1223           wcp = wbuffer + wbuffer_to_alloc;
1224           assert ((uintptr_t) wbuffer <= (uintptr_t) wstartp);
1225           assert ((uintptr_t) wstartp
1226                   < (uintptr_t) wbuffer + wbuffer_to_alloc);
1227 #else
1228           tmpptr = _i18n_number_rewrite (tmpptr, cp, buffer_end);
1229           cp = buffer_end;
1230           assert ((uintptr_t) buffer <= (uintptr_t) tmpptr);
1231           assert ((uintptr_t) tmpptr < (uintptr_t) buffer_end);
1232 #endif
1233         }
1234
1235       PRINT (tmpptr, wstartp, wide ? wcp - wstartp : cp - tmpptr);
1236
1237       /* Free the memory if necessary.  */
1238       if (__glibc_unlikely (buffer_malloced))
1239         {
1240           free (buffer);
1241           free (wbuffer);
1242         }
1243     }
1244
1245     if (info->left && width > 0)
1246       PADN (info->pad, width);
1247   }
1248   return done;
1249 }
1250 ldbl_hidden_def (___printf_fp, __printf_fp)
1251 ldbl_strong_alias (___printf_fp, __printf_fp)
1252 \f
1253 /* Return the number of extra grouping characters that will be inserted
1254    into a number with INTDIG_MAX integer digits.  */
1255
1256 unsigned int
1257 __guess_grouping (unsigned int intdig_max, const char *grouping)
1258 {
1259   unsigned int groups;
1260
1261   /* We treat all negative values like CHAR_MAX.  */
1262
1263   if (*grouping == CHAR_MAX || *grouping <= 0)
1264     /* No grouping should be done.  */
1265     return 0;
1266
1267   groups = 0;
1268   while (intdig_max > (unsigned int) *grouping)
1269     {
1270       ++groups;
1271       intdig_max -= *grouping++;
1272
1273       if (*grouping == CHAR_MAX
1274 #if CHAR_MIN < 0
1275           || *grouping < 0
1276 #endif
1277           )
1278         /* No more grouping should be done.  */
1279         break;
1280       else if (*grouping == 0)
1281         {
1282           /* Same grouping repeats.  */
1283           groups += (intdig_max - 1) / grouping[-1];
1284           break;
1285         }
1286     }
1287
1288   return groups;
1289 }
1290
1291 /* Group the INTDIG_NO integer digits of the number in [BUF,BUFEND).
1292    There is guaranteed enough space past BUFEND to extend it.
1293    Return the new end of buffer.  */
1294
1295 static wchar_t *
1296 internal_function
1297 group_number (wchar_t *buf, wchar_t *bufend, unsigned int intdig_no,
1298               const char *grouping, wchar_t thousands_sep, int ngroups)
1299 {
1300   wchar_t *p;
1301
1302   if (ngroups == 0)
1303     return bufend;
1304
1305   /* Move the fractional part down.  */
1306   __wmemmove (buf + intdig_no + ngroups, buf + intdig_no,
1307               bufend - (buf + intdig_no));
1308
1309   p = buf + intdig_no + ngroups - 1;
1310   do
1311     {
1312       unsigned int len = *grouping++;
1313       do
1314         *p-- = buf[--intdig_no];
1315       while (--len > 0);
1316       *p-- = thousands_sep;
1317
1318       if (*grouping == CHAR_MAX
1319 #if CHAR_MIN < 0
1320           || *grouping < 0
1321 #endif
1322           )
1323         /* No more grouping should be done.  */
1324         break;
1325       else if (*grouping == 0)
1326         /* Same grouping repeats.  */
1327         --grouping;
1328     } while (intdig_no > (unsigned int) *grouping);
1329
1330   /* Copy the remaining ungrouped digits.  */
1331   do
1332     *p-- = buf[--intdig_no];
1333   while (p > buf);
1334
1335   return bufend + ngroups;
1336 }