1 /* Convert string representing a number to float value, using given locale.
2 Copyright (C) 1997-2017 Free Software Foundation, Inc.
3 This file is part of the GNU C Library.
4 Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
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.
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.
16 You should have received a copy of the GNU Lesser General Public
17 License along with the GNU C Library; if not, see
18 <http://www.gnu.org/licenses/>. */
22 extern double ____strtod_l_internal (const char *, char **, int, __locale_t
);
24 /* Configuration part. These macros are defined by `strtold.c',
25 `strtof.c', `wcstod.c', `wcstold.c', and `wcstof.c' to produce the
26 `long double' and `float' versions of the reader. */
28 # include <math_ldbl_opt.h>
32 # define STRTOF wcstod_l
33 # define __STRTOF __wcstod_l
34 # define STRTOF_NAN __wcstod_nan
36 # define STRTOF strtod_l
37 # define __STRTOF __strtod_l
38 # define STRTOF_NAN __strtod_nan
40 # define MPN2FLOAT __mpn_construct_double
41 # define FLOAT_HUGE_VAL HUGE_VAL
43 /* End of configuration part. */
48 #include "../locale/localeinfo.h"
51 #include <math_private.h>
55 #include <rounding-mode.h>
58 /* The gmp headers need some configuration frobs. */
61 /* Include gmp-mparam.h first, such that definitions of _SHORT_LIMB
62 and _LONG_LONG_LIMB in it can take effect into gmp.h. */
63 #include <gmp-mparam.h>
67 #include "fpioconst.h"
72 /* We use this code for the extended locale handling where the
73 function gets as an additional argument the locale which has to be
74 used. To access the values we have to redefine the _NL_CURRENT and
75 _NL_CURRENT_WORD macros. */
77 #define _NL_CURRENT(category, item) \
78 (current->values[_NL_ITEM_INDEX (item)].string)
79 #undef _NL_CURRENT_WORD
80 #define _NL_CURRENT_WORD(category, item) \
81 ((uint32_t) current->values[_NL_ITEM_INDEX (item)].word)
83 #if defined _LIBC || defined HAVE_WCHAR_H
89 # define STRING_TYPE wchar_t
90 # define CHAR_TYPE wint_t
92 # define ISSPACE(Ch) __iswspace_l ((Ch), loc)
93 # define ISDIGIT(Ch) __iswdigit_l ((Ch), loc)
94 # define ISXDIGIT(Ch) __iswxdigit_l ((Ch), loc)
95 # define TOLOWER(Ch) __towlower_l ((Ch), loc)
96 # define TOLOWER_C(Ch) __towlower_l ((Ch), _nl_C_locobj_ptr)
97 # define STRNCASECMP(S1, S2, N) \
98 __wcsncasecmp_l ((S1), (S2), (N), _nl_C_locobj_ptr)
100 # define STRING_TYPE char
101 # define CHAR_TYPE char
103 # define ISSPACE(Ch) __isspace_l ((Ch), loc)
104 # define ISDIGIT(Ch) __isdigit_l ((Ch), loc)
105 # define ISXDIGIT(Ch) __isxdigit_l ((Ch), loc)
106 # define TOLOWER(Ch) __tolower_l ((Ch), loc)
107 # define TOLOWER_C(Ch) __tolower_l ((Ch), _nl_C_locobj_ptr)
108 # define STRNCASECMP(S1, S2, N) \
109 __strncasecmp_l ((S1), (S2), (N), _nl_C_locobj_ptr)
113 /* Constants we need from float.h; select the set for the FLOAT precision. */
114 #define MANT_DIG PASTE(FLT,_MANT_DIG)
115 #define DIG PASTE(FLT,_DIG)
116 #define MAX_EXP PASTE(FLT,_MAX_EXP)
117 #define MIN_EXP PASTE(FLT,_MIN_EXP)
118 #define MAX_10_EXP PASTE(FLT,_MAX_10_EXP)
119 #define MIN_10_EXP PASTE(FLT,_MIN_10_EXP)
120 #define MAX_VALUE PASTE(FLT,_MAX)
121 #define MIN_VALUE PASTE(FLT,_MIN)
123 /* Extra macros required to get FLT expanded before the pasting. */
124 #define PASTE(a,b) PASTE1(a,b)
125 #define PASTE1(a,b) a##b
127 /* Function to construct a floating point number from an MP integer
128 containing the fraction bits, a base 2 exponent, and a sign flag. */
129 extern FLOAT
MPN2FLOAT (mp_srcptr mpn
, int exponent
, int negative
);
131 /* Definitions according to limb size used. */
132 #if BITS_PER_MP_LIMB == 32
133 # define MAX_DIG_PER_LIMB 9
134 # define MAX_FAC_PER_LIMB 1000000000UL
135 #elif BITS_PER_MP_LIMB == 64
136 # define MAX_DIG_PER_LIMB 19
137 # define MAX_FAC_PER_LIMB 10000000000000000000ULL
139 # error "mp_limb_t size " BITS_PER_MP_LIMB "not accounted for"
142 extern const mp_limb_t _tens_in_limb
[MAX_DIG_PER_LIMB
+ 1];
145 #define howmany(x,y) (((x)+((y)-1))/(y))
147 #define SWAP(x, y) ({ typeof(x) _tmp = x; x = y; y = _tmp; })
149 #define RETURN_LIMB_SIZE howmany (MANT_DIG, BITS_PER_MP_LIMB)
151 #define RETURN(val,end) \
152 do { if (endptr != NULL) *endptr = (STRING_TYPE *) (end); \
153 return val; } while (0)
155 /* Maximum size necessary for mpn integers to hold floating point
156 numbers. The largest number we need to hold is 10^n where 2^-n is
157 1/4 ulp of the smallest representable value (that is, n = MANT_DIG
158 - MIN_EXP + 2). Approximate using 10^3 < 2^10. */
159 #define MPNSIZE (howmany (1 + ((MANT_DIG - MIN_EXP + 2) * 10) / 3, \
160 BITS_PER_MP_LIMB) + 2)
161 /* Declare an mpn integer variable that big. */
162 #define MPN_VAR(name) mp_limb_t name[MPNSIZE]; mp_size_t name##size
163 /* Copy an mpn integer value. */
164 #define MPN_ASSIGN(dst, src) \
165 memcpy (dst, src, (dst##size = src##size) * sizeof (mp_limb_t))
168 /* Set errno and return an overflowing value with sign specified by
171 overflow_value (int negative
)
173 __set_errno (ERANGE
);
174 FLOAT result
= math_narrow_eval ((negative
? -MAX_VALUE
: MAX_VALUE
)
180 /* Set errno and return an underflowing value with sign specified by
183 underflow_value (int negative
)
185 __set_errno (ERANGE
);
186 FLOAT result
= math_narrow_eval ((negative
? -MIN_VALUE
: MIN_VALUE
)
192 /* Return a floating point number of the needed type according to the given
193 multi-precision number after possible rounding. */
195 round_and_return (mp_limb_t
*retval
, intmax_t exponent
, int negative
,
196 mp_limb_t round_limb
, mp_size_t round_bit
, int more_bits
)
198 int mode
= get_rounding_mode ();
200 if (exponent
< MIN_EXP
- 1)
202 if (exponent
< MIN_EXP
- 1 - MANT_DIG
)
203 return underflow_value (negative
);
205 mp_size_t shift
= MIN_EXP
- 1 - exponent
;
208 more_bits
|= (round_limb
& ((((mp_limb_t
) 1) << round_bit
) - 1)) != 0;
209 if (shift
== MANT_DIG
)
210 /* This is a special case to handle the very seldom case where
211 the mantissa will be empty after the shift. */
215 round_limb
= retval
[RETURN_LIMB_SIZE
- 1];
216 round_bit
= (MANT_DIG
- 1) % BITS_PER_MP_LIMB
;
217 for (i
= 0; i
< RETURN_LIMB_SIZE
- 1; ++i
)
218 more_bits
|= retval
[i
] != 0;
219 MPN_ZERO (retval
, RETURN_LIMB_SIZE
);
221 else if (shift
>= BITS_PER_MP_LIMB
)
225 round_limb
= retval
[(shift
- 1) / BITS_PER_MP_LIMB
];
226 round_bit
= (shift
- 1) % BITS_PER_MP_LIMB
;
227 for (i
= 0; i
< (shift
- 1) / BITS_PER_MP_LIMB
; ++i
)
228 more_bits
|= retval
[i
] != 0;
229 more_bits
|= ((round_limb
& ((((mp_limb_t
) 1) << round_bit
) - 1))
232 /* __mpn_rshift requires 0 < shift < BITS_PER_MP_LIMB. */
233 if ((shift
% BITS_PER_MP_LIMB
) != 0)
234 (void) __mpn_rshift (retval
, &retval
[shift
/ BITS_PER_MP_LIMB
],
235 RETURN_LIMB_SIZE
- (shift
/ BITS_PER_MP_LIMB
),
236 shift
% BITS_PER_MP_LIMB
);
238 for (i
= 0; i
< RETURN_LIMB_SIZE
- (shift
/ BITS_PER_MP_LIMB
); i
++)
239 retval
[i
] = retval
[i
+ (shift
/ BITS_PER_MP_LIMB
)];
240 MPN_ZERO (&retval
[RETURN_LIMB_SIZE
- (shift
/ BITS_PER_MP_LIMB
)],
241 shift
/ BITS_PER_MP_LIMB
);
245 if (TININESS_AFTER_ROUNDING
&& shift
== 1)
247 /* Whether the result counts as tiny depends on whether,
248 after rounding to the normal precision, it still has
249 a subnormal exponent. */
250 mp_limb_t retval_normal
[RETURN_LIMB_SIZE
];
251 if (round_away (negative
,
252 (retval
[0] & 1) != 0,
254 & (((mp_limb_t
) 1) << round_bit
)) != 0,
257 & ((((mp_limb_t
) 1) << round_bit
) - 1))
261 mp_limb_t cy
= __mpn_add_1 (retval_normal
, retval
,
262 RETURN_LIMB_SIZE
, 1);
264 if (((MANT_DIG
% BITS_PER_MP_LIMB
) == 0 && cy
) ||
265 ((MANT_DIG
% BITS_PER_MP_LIMB
) != 0 &&
266 ((retval_normal
[RETURN_LIMB_SIZE
- 1]
267 & (((mp_limb_t
) 1) << (MANT_DIG
% BITS_PER_MP_LIMB
)))
272 round_limb
= retval
[0];
273 round_bit
= shift
- 1;
274 (void) __mpn_rshift (retval
, retval
, RETURN_LIMB_SIZE
, shift
);
276 /* This is a hook for the m68k long double format, where the
277 exponent bias is the same for normalized and denormalized
280 # define DENORM_EXP (MIN_EXP - 2)
282 exponent
= DENORM_EXP
;
284 && ((round_limb
& (((mp_limb_t
) 1) << round_bit
)) != 0
286 || (round_limb
& ((((mp_limb_t
) 1) << round_bit
) - 1)) != 0))
288 __set_errno (ERANGE
);
289 FLOAT force_underflow
= MIN_VALUE
* MIN_VALUE
;
290 math_force_eval (force_underflow
);
294 if (exponent
> MAX_EXP
)
297 bool half_bit
= (round_limb
& (((mp_limb_t
) 1) << round_bit
)) != 0;
298 bool more_bits_nonzero
300 || (round_limb
& ((((mp_limb_t
) 1) << round_bit
) - 1)) != 0);
301 if (round_away (negative
,
302 (retval
[0] & 1) != 0,
307 mp_limb_t cy
= __mpn_add_1 (retval
, retval
, RETURN_LIMB_SIZE
, 1);
309 if (((MANT_DIG
% BITS_PER_MP_LIMB
) == 0 && cy
) ||
310 ((MANT_DIG
% BITS_PER_MP_LIMB
) != 0 &&
311 (retval
[RETURN_LIMB_SIZE
- 1]
312 & (((mp_limb_t
) 1) << (MANT_DIG
% BITS_PER_MP_LIMB
))) != 0))
315 (void) __mpn_rshift (retval
, retval
, RETURN_LIMB_SIZE
, 1);
316 retval
[RETURN_LIMB_SIZE
- 1]
317 |= ((mp_limb_t
) 1) << ((MANT_DIG
- 1) % BITS_PER_MP_LIMB
);
319 else if (exponent
== DENORM_EXP
320 && (retval
[RETURN_LIMB_SIZE
- 1]
321 & (((mp_limb_t
) 1) << ((MANT_DIG
- 1) % BITS_PER_MP_LIMB
)))
323 /* The number was denormalized but now normalized. */
324 exponent
= MIN_EXP
- 1;
327 if (exponent
> MAX_EXP
)
329 return overflow_value (negative
);
331 if (half_bit
|| more_bits_nonzero
)
333 FLOAT force_inexact
= (FLOAT
) 1 + MIN_VALUE
;
334 math_force_eval (force_inexact
);
336 return MPN2FLOAT (retval
, exponent
, negative
);
340 /* Read a multi-precision integer starting at STR with exactly DIGCNT digits
341 into N. Return the size of the number limbs in NSIZE at the first
342 character od the string that is not part of the integer as the function
343 value. If the EXPONENT is small enough to be taken as an additional
344 factor for the resulting number (see code) multiply by it. */
345 static const STRING_TYPE
*
346 str_to_mpn (const STRING_TYPE
*str
, int digcnt
, mp_limb_t
*n
, mp_size_t
*nsize
,
348 #ifndef USE_WIDE_CHAR
349 , const char *decimal
, size_t decimal_len
, const char *thousands
354 /* Number of digits for actual limb. */
363 if (cnt
== MAX_DIG_PER_LIMB
)
373 cy
= __mpn_mul_1 (n
, n
, *nsize
, MAX_FAC_PER_LIMB
);
374 cy
+= __mpn_add_1 (n
, n
, *nsize
, low
);
377 assert (*nsize
< MPNSIZE
);
386 /* There might be thousands separators or radix characters in
387 the string. But these all can be ignored because we know the
388 format of the number is correct and we have an exact number
389 of characters to read. */
391 if (*str
< L
'0' || *str
> L
'9')
394 if (*str
< '0' || *str
> '9')
397 if (thousands
!= NULL
&& *str
== *thousands
398 && ({ for (inner
= 1; thousands
[inner
] != '\0'; ++inner
)
399 if (thousands
[inner
] != str
[inner
])
401 thousands
[inner
] == '\0'; }))
407 low
= low
* 10 + *str
++ - L_('0');
410 while (--digcnt
> 0);
412 if (*exponent
> 0 && *exponent
<= MAX_DIG_PER_LIMB
- cnt
)
414 low
*= _tens_in_limb
[*exponent
];
415 start
= _tens_in_limb
[cnt
+ *exponent
];
419 start
= _tens_in_limb
[cnt
];
429 cy
= __mpn_mul_1 (n
, n
, *nsize
, start
);
430 cy
+= __mpn_add_1 (n
, n
, *nsize
, low
);
433 assert (*nsize
< MPNSIZE
);
442 /* Shift {PTR, SIZE} COUNT bits to the left, and fill the vacated bits
443 with the COUNT most significant bits of LIMB.
445 Implemented as a macro, so that __builtin_constant_p works even at -O0.
447 Tege doesn't like this macro so I have to write it here myself. :)
449 #define __mpn_lshift_1(ptr, size, count, limb) \
452 mp_limb_t *__ptr = (ptr); \
453 if (__builtin_constant_p (count) && count == BITS_PER_MP_LIMB) \
456 for (i = (size) - 1; i > 0; --i) \
457 __ptr[i] = __ptr[i - 1]; \
462 /* We assume count > 0 && count < BITS_PER_MP_LIMB here. */ \
463 unsigned int __count = (count); \
464 (void) __mpn_lshift (__ptr, __ptr, size, __count); \
465 __ptr[0] |= (limb) >> (BITS_PER_MP_LIMB - __count); \
471 #define INTERNAL(x) INTERNAL1(x)
472 #define INTERNAL1(x) __##x##_internal
473 #ifndef ____STRTOF_INTERNAL
474 # define ____STRTOF_INTERNAL INTERNAL (__STRTOF)
477 /* This file defines a function to check for correct grouping. */
478 #include "grouping.h"
481 /* Return a floating point number with the value of the given string NPTR.
482 Set *ENDPTR to the character after the last used one. If the number is
483 smaller than the smallest representable number, set `errno' to ERANGE and
484 return 0.0. If the number is too big to be represented, set `errno' to
485 ERANGE and return HUGE_VAL with the appropriate sign. */
487 ____STRTOF_INTERNAL (const STRING_TYPE
*nptr
, STRING_TYPE
**endptr
, int group
,
490 int negative
; /* The sign of the number. */
491 MPN_VAR (num
); /* MP representation of the number. */
492 intmax_t exponent
; /* Exponent of the number. */
494 /* Numbers starting `0X' or `0x' have to be processed with base 16. */
497 /* When we have to compute fractional digits we form a fraction with a
498 second multi-precision number (and we sometimes need a second for
499 temporary results). */
502 /* Representation for the return value. */
503 mp_limb_t retval
[RETURN_LIMB_SIZE
];
504 /* Number of bits currently in result value. */
507 /* Running pointer after the last character processed in the string. */
508 const STRING_TYPE
*cp
, *tp
;
509 /* Start of significant part of the number. */
510 const STRING_TYPE
*startp
, *start_of_digits
;
511 /* Points at the character following the integer and fractional digits. */
512 const STRING_TYPE
*expp
;
513 /* Total number of digit and number of digits in integer part. */
514 size_t dig_no
, int_no
, lead_zero
;
515 /* Contains the last character read. */
518 /* We should get wint_t from <stddef.h>, but not all GCC versions define it
519 there. So define it ourselves if it remains undefined. */
521 typedef unsigned int wint_t;
523 /* The radix character of the current locale. */
530 /* The thousands character of the current locale. */
532 wchar_t thousands
= L
'\0';
534 const char *thousands
= NULL
;
536 /* The numeric grouping specification of the current locale,
537 in the format described in <locale.h>. */
538 const char *grouping
;
539 /* Used in several places. */
542 struct __locale_data
*current
= loc
->__locales
[LC_NUMERIC
];
544 if (__glibc_unlikely (group
))
546 grouping
= _NL_CURRENT (LC_NUMERIC
, GROUPING
);
547 if (*grouping
<= 0 || *grouping
== CHAR_MAX
)
551 /* Figure out the thousands separator character. */
553 thousands
= _NL_CURRENT_WORD (LC_NUMERIC
,
554 _NL_NUMERIC_THOUSANDS_SEP_WC
);
555 if (thousands
== L
'\0')
558 thousands
= _NL_CURRENT (LC_NUMERIC
, THOUSANDS_SEP
);
559 if (*thousands
== '\0')
570 /* Find the locale's decimal point character. */
572 decimal
= _NL_CURRENT_WORD (LC_NUMERIC
, _NL_NUMERIC_DECIMAL_POINT_WC
);
573 assert (decimal
!= L
'\0');
574 # define decimal_len 1
576 decimal
= _NL_CURRENT (LC_NUMERIC
, DECIMAL_POINT
);
577 decimal_len
= strlen (decimal
);
578 assert (decimal_len
> 0);
581 /* Prepare number representation. */
586 /* Parse string to get maximal legal prefix. We need the number of
587 characters of the integer part, the fractional part and the exponent. */
589 /* Ignore leading white space. */
594 /* Get sign of the result. */
600 else if (c
== L_('+'))
603 /* Return 0.0 if no legal string is found.
604 No character is used even if a sign was found. */
606 if (c
== (wint_t) decimal
607 && (wint_t) cp
[1] >= L
'0' && (wint_t) cp
[1] <= L
'9')
609 /* We accept it. This funny construct is here only to indent
610 the code correctly. */
613 for (cnt
= 0; decimal
[cnt
] != '\0'; ++cnt
)
614 if (cp
[cnt
] != decimal
[cnt
])
616 if (decimal
[cnt
] == '\0' && cp
[cnt
] >= '0' && cp
[cnt
] <= '9')
618 /* We accept it. This funny construct is here only to indent
619 the code correctly. */
622 else if (c
< L_('0') || c
> L_('9'))
624 /* Check for `INF' or `INFINITY'. */
625 CHAR_TYPE lowc
= TOLOWER_C (c
);
627 if (lowc
== L_('i') && STRNCASECMP (cp
, L_("inf"), 3) == 0)
629 /* Return +/- infinity. */
631 *endptr
= (STRING_TYPE
*)
632 (cp
+ (STRNCASECMP (cp
+ 3, L_("inity"), 5) == 0
635 return negative
? -FLOAT_HUGE_VAL
: FLOAT_HUGE_VAL
;
638 if (lowc
== L_('n') && STRNCASECMP (cp
, L_("nan"), 3) == 0)
645 /* Match `(n-char-sequence-digit)'. */
648 const STRING_TYPE
*startp
= cp
;
650 retval
= STRTOF_NAN (cp
+ 1, &endp
, L_(')'));
651 if (*endp
== L_(')'))
652 /* Consume the closing parenthesis. */
655 /* Only match the NAN part. */
660 *endptr
= (STRING_TYPE
*) cp
;
665 /* It is really a text we do not recognize. */
669 /* First look whether we are faced with a hexadecimal number. */
670 if (c
== L_('0') && TOLOWER (cp
[1]) == L_('x'))
672 /* Okay, it is a hexa-decimal number. Remember this and skip
673 the characters. BTW: hexadecimal numbers must not be
681 /* Record the start of the digits, in case we will check their grouping. */
682 start_of_digits
= startp
= cp
;
684 /* Ignore leading zeroes. This helps us to avoid useless computations. */
686 while (c
== L
'0' || ((wint_t) thousands
!= L
'\0' && c
== (wint_t) thousands
))
689 if (__glibc_likely (thousands
== NULL
))
694 /* We also have the multibyte thousands string. */
699 for (cnt
= 0; thousands
[cnt
] != '\0'; ++cnt
)
700 if (thousands
[cnt
] != cp
[cnt
])
702 if (thousands
[cnt
] != '\0')
711 /* If no other digit but a '0' is found the result is 0.0.
712 Return current read pointer. */
713 CHAR_TYPE lowc
= TOLOWER (c
);
714 if (!((c
>= L_('0') && c
<= L_('9'))
715 || (base
== 16 && lowc
>= L_('a') && lowc
<= L_('f'))
718 c
== (wint_t) decimal
720 ({ for (cnt
= 0; decimal
[cnt
] != '\0'; ++cnt
)
721 if (decimal
[cnt
] != cp
[cnt
])
723 decimal
[cnt
] == '\0'; })
725 /* '0x.' alone is not a valid hexadecimal number.
726 '.' alone is not valid either, but that has been checked
729 || cp
!= start_of_digits
730 || (cp
[decimal_len
] >= L_('0') && cp
[decimal_len
] <= L_('9'))
731 || ({ CHAR_TYPE lo
= TOLOWER (cp
[decimal_len
]);
732 lo
>= L_('a') && lo
<= L_('f'); })))
733 || (base
== 16 && (cp
!= start_of_digits
735 || (base
!= 16 && lowc
== L_('e'))))
738 tp
= __correctly_grouped_prefixwc (start_of_digits
, cp
, thousands
,
741 tp
= __correctly_grouped_prefixmb (start_of_digits
, cp
, thousands
,
744 /* If TP is at the start of the digits, there was no correctly
745 grouped prefix of the string; so no number found. */
746 RETURN (negative
? -0.0 : 0.0,
747 tp
== start_of_digits
? (base
== 16 ? cp
- 1 : nptr
) : tp
);
750 /* Remember first significant digit and read following characters until the
751 decimal point, exponent character or any non-FP number character. */
756 if ((c
>= L_('0') && c
<= L_('9'))
758 && ({ CHAR_TYPE lo
= TOLOWER (c
);
759 lo
>= L_('a') && lo
<= L_('f'); })))
764 if (__builtin_expect ((wint_t) thousands
== L
'\0', 1)
765 || c
!= (wint_t) thousands
)
766 /* Not a digit or separator: end of the integer part. */
769 if (__glibc_likely (thousands
== NULL
))
773 for (cnt
= 0; thousands
[cnt
] != '\0'; ++cnt
)
774 if (thousands
[cnt
] != cp
[cnt
])
776 if (thousands
[cnt
] != '\0')
785 if (__builtin_expect (grouping
!= NULL
, 0) && cp
> start_of_digits
)
787 /* Check the grouping of the digits. */
789 tp
= __correctly_grouped_prefixwc (start_of_digits
, cp
, thousands
,
792 tp
= __correctly_grouped_prefixmb (start_of_digits
, cp
, thousands
,
797 /* Less than the entire string was correctly grouped. */
799 if (tp
== start_of_digits
)
800 /* No valid group of numbers at all: no valid number. */
804 /* The number is validly grouped, but consists
805 only of zeroes. The whole value is zero. */
806 RETURN (negative
? -0.0 : 0.0, tp
);
808 /* Recompute DIG_NO so we won't read more digits than
809 are properly grouped. */
812 for (tp
= startp
; tp
< cp
; ++tp
)
813 if (*tp
>= L_('0') && *tp
<= L_('9'))
823 /* We have the number of digits in the integer part. Whether these
824 are all or any is really a fractional digit will be decided
827 lead_zero
= int_no
== 0 ? (size_t) -1 : 0;
829 /* Read the fractional digits. A special case are the 'american
830 style' numbers like `16.' i.e. with decimal point but without
834 c
== (wint_t) decimal
836 ({ for (cnt
= 0; decimal
[cnt
] != '\0'; ++cnt
)
837 if (decimal
[cnt
] != cp
[cnt
])
839 decimal
[cnt
] == '\0'; })
845 while ((c
>= L_('0') && c
<= L_('9')) ||
846 (base
== 16 && ({ CHAR_TYPE lo
= TOLOWER (c
);
847 lo
>= L_('a') && lo
<= L_('f'); })))
849 if (c
!= L_('0') && lead_zero
== (size_t) -1)
850 lead_zero
= dig_no
- int_no
;
855 assert (dig_no
<= (uintmax_t) INTMAX_MAX
);
857 /* Remember start of exponent (if any). */
862 if ((base
== 16 && lowc
== L_('p'))
863 || (base
!= 16 && lowc
== L_('e')))
865 int exp_negative
= 0;
873 else if (c
== L_('+'))
876 if (c
>= L_('0') && c
<= L_('9'))
880 /* Get the exponent limit. */
885 assert (int_no
<= (uintmax_t) (INTMAX_MAX
886 + MIN_EXP
- MANT_DIG
) / 4);
887 exp_limit
= -MIN_EXP
+ MANT_DIG
+ 4 * (intmax_t) int_no
;
893 assert (lead_zero
== 0
894 && int_no
<= (uintmax_t) INTMAX_MAX
/ 4);
895 exp_limit
= MAX_EXP
- 4 * (intmax_t) int_no
+ 3;
897 else if (lead_zero
== (size_t) -1)
899 /* The number is zero and this limit is
901 exp_limit
= MAX_EXP
+ 3;
906 <= (uintmax_t) (INTMAX_MAX
- MAX_EXP
- 3) / 4);
908 + 4 * (intmax_t) lead_zero
918 <= (uintmax_t) (INTMAX_MAX
+ MIN_10_EXP
- MANT_DIG
));
919 exp_limit
= -MIN_10_EXP
+ MANT_DIG
+ (intmax_t) int_no
;
925 assert (lead_zero
== 0
926 && int_no
<= (uintmax_t) INTMAX_MAX
);
927 exp_limit
= MAX_10_EXP
- (intmax_t) int_no
+ 1;
929 else if (lead_zero
== (size_t) -1)
931 /* The number is zero and this limit is
933 exp_limit
= MAX_10_EXP
+ 1;
938 <= (uintmax_t) (INTMAX_MAX
- MAX_10_EXP
- 1));
939 exp_limit
= MAX_10_EXP
+ (intmax_t) lead_zero
+ 1;
949 if (__builtin_expect ((exponent
> exp_limit
/ 10
950 || (exponent
== exp_limit
/ 10
951 && c
- L_('0') > exp_limit
% 10)), 0))
952 /* The exponent is too large/small to represent a valid
957 /* We have to take care for special situation: a joker
958 might have written "0.0e100000" which is in fact
960 if (lead_zero
== (size_t) -1)
961 result
= negative
? -0.0 : 0.0;
964 /* Overflow or underflow. */
965 result
= (exp_negative
966 ? underflow_value (negative
)
967 : overflow_value (negative
));
970 /* Accept all following digits as part of the exponent. */
973 while (*cp
>= L_('0') && *cp
<= L_('9'));
980 exponent
+= c
- L_('0');
984 while (c
>= L_('0') && c
<= L_('9'));
987 exponent
= -exponent
;
993 /* We don't want to have to work with trailing zeroes after the radix. */
996 while (expp
[-1] == L_('0'))
1001 assert (dig_no
>= int_no
);
1004 if (dig_no
== int_no
&& dig_no
> 0 && exponent
< 0)
1007 while (! (base
== 16 ? ISXDIGIT (expp
[-1]) : ISDIGIT (expp
[-1])))
1010 if (expp
[-1] != L_('0'))
1016 exponent
+= base
== 16 ? 4 : 1;
1018 while (dig_no
> 0 && exponent
< 0);
1022 /* The whole string is parsed. Store the address of the next character. */
1024 *endptr
= (STRING_TYPE
*) cp
;
1027 return negative
? -0.0 : 0.0;
1031 /* Find the decimal point */
1032 #ifdef USE_WIDE_CHAR
1033 while (*startp
!= decimal
)
1038 if (*startp
== decimal
[0])
1040 for (cnt
= 1; decimal
[cnt
] != '\0'; ++cnt
)
1041 if (decimal
[cnt
] != startp
[cnt
])
1043 if (decimal
[cnt
] == '\0')
1049 startp
+= lead_zero
+ decimal_len
;
1050 assert (lead_zero
<= (base
== 16
1051 ? (uintmax_t) INTMAX_MAX
/ 4
1052 : (uintmax_t) INTMAX_MAX
));
1053 assert (lead_zero
<= (base
== 16
1054 ? ((uintmax_t) exponent
1055 - (uintmax_t) INTMAX_MIN
) / 4
1056 : ((uintmax_t) exponent
- (uintmax_t) INTMAX_MIN
)));
1057 exponent
-= base
== 16 ? 4 * (intmax_t) lead_zero
: (intmax_t) lead_zero
;
1058 dig_no
-= lead_zero
;
1061 /* If the BASE is 16 we can use a simpler algorithm. */
1064 static const int nbits
[16] = { 0, 1, 2, 2, 3, 3, 3, 3,
1065 4, 4, 4, 4, 4, 4, 4, 4 };
1066 int idx
= (MANT_DIG
- 1) / BITS_PER_MP_LIMB
;
1067 int pos
= (MANT_DIG
- 1) % BITS_PER_MP_LIMB
;
1070 while (!ISXDIGIT (*startp
))
1072 while (*startp
== L_('0'))
1074 if (ISDIGIT (*startp
))
1075 val
= *startp
++ - L_('0');
1077 val
= 10 + TOLOWER (*startp
++) - L_('a');
1079 /* We cannot have a leading zero. */
1082 if (pos
+ 1 >= 4 || pos
+ 1 >= bits
)
1084 /* We don't have to care for wrapping. This is the normal
1085 case so we add the first clause in the `if' expression as
1086 an optimization. It is a compile-time constant and so does
1087 not cost anything. */
1088 retval
[idx
] = val
<< (pos
- bits
+ 1);
1093 retval
[idx
--] = val
>> (bits
- pos
- 1);
1094 retval
[idx
] = val
<< (BITS_PER_MP_LIMB
- (bits
- pos
- 1));
1095 pos
= BITS_PER_MP_LIMB
- 1 - (bits
- pos
- 1);
1098 /* Adjust the exponent for the bits we are shifting in. */
1099 assert (int_no
<= (uintmax_t) (exponent
< 0
1100 ? (INTMAX_MAX
- bits
+ 1) / 4
1101 : (INTMAX_MAX
- exponent
- bits
+ 1) / 4));
1102 exponent
+= bits
- 1 + ((intmax_t) int_no
- 1) * 4;
1104 while (--dig_no
> 0 && idx
>= 0)
1106 if (!ISXDIGIT (*startp
))
1107 startp
+= decimal_len
;
1108 if (ISDIGIT (*startp
))
1109 val
= *startp
++ - L_('0');
1111 val
= 10 + TOLOWER (*startp
++) - L_('a');
1115 retval
[idx
] |= val
<< (pos
- 4 + 1);
1120 retval
[idx
--] |= val
>> (4 - pos
- 1);
1121 val
<<= BITS_PER_MP_LIMB
- (4 - pos
- 1);
1124 int rest_nonzero
= 0;
1125 while (--dig_no
> 0)
1127 if (*startp
!= L_('0'))
1134 return round_and_return (retval
, exponent
, negative
, val
,
1135 BITS_PER_MP_LIMB
- 1, rest_nonzero
);
1139 pos
= BITS_PER_MP_LIMB
- 1 - (4 - pos
- 1);
1143 /* We ran out of digits. */
1144 MPN_ZERO (retval
, idx
);
1146 return round_and_return (retval
, exponent
, negative
, 0, 0, 0);
1149 /* Now we have the number of digits in total and the integer digits as well
1150 as the exponent and its sign. We can decide whether the read digits are
1151 really integer digits or belong to the fractional part; i.e. we normalize
1154 intmax_t incr
= (exponent
< 0
1155 ? MAX (-(intmax_t) int_no
, exponent
)
1156 : MIN ((intmax_t) dig_no
- (intmax_t) int_no
, exponent
));
1161 if (__glibc_unlikely (exponent
> MAX_10_EXP
+ 1 - (intmax_t) int_no
))
1162 return overflow_value (negative
);
1164 /* 10^(MIN_10_EXP-1) is not normal. Thus, 10^(MIN_10_EXP-1) /
1165 2^MANT_DIG is below half the least subnormal, so anything with a
1166 base-10 exponent less than the base-10 exponent (which is
1167 MIN_10_EXP - 1 - ceil(MANT_DIG*log10(2))) of that value
1168 underflows. DIG is floor((MANT_DIG-1)log10(2)), so an exponent
1169 below MIN_10_EXP - (DIG + 3) underflows. But EXPONENT is
1170 actually an exponent multiplied only by a fractional part, not an
1171 integer part, so an exponent below MIN_10_EXP - (DIG + 2)
1173 if (__glibc_unlikely (exponent
< MIN_10_EXP
- (DIG
+ 2)))
1174 return underflow_value (negative
);
1178 /* Read the integer part as a multi-precision number to NUM. */
1179 startp
= str_to_mpn (startp
, int_no
, num
, &numsize
, &exponent
1180 #ifndef USE_WIDE_CHAR
1181 , decimal
, decimal_len
, thousands
1187 /* We now multiply the gained number by the given power of ten. */
1188 mp_limb_t
*psrc
= num
;
1189 mp_limb_t
*pdest
= den
;
1191 const struct mp_power
*ttab
= &_fpioconst_pow10
[0];
1195 if ((exponent
& expbit
) != 0)
1197 size_t size
= ttab
->arraysize
- _FPIO_CONST_OFFSET
;
1201 /* FIXME: not the whole multiplication has to be
1202 done. If we have the needed number of bits we
1203 only need the information whether more non-zero
1205 if (numsize
>= ttab
->arraysize
- _FPIO_CONST_OFFSET
)
1206 cy
= __mpn_mul (pdest
, psrc
, numsize
,
1207 &__tens
[ttab
->arrayoff
1208 + _FPIO_CONST_OFFSET
],
1211 cy
= __mpn_mul (pdest
, &__tens
[ttab
->arrayoff
1212 + _FPIO_CONST_OFFSET
],
1213 size
, psrc
, numsize
);
1217 (void) SWAP (psrc
, pdest
);
1222 while (exponent
!= 0);
1225 memcpy (num
, den
, numsize
* sizeof (mp_limb_t
));
1228 /* Determine how many bits of the result we already have. */
1229 count_leading_zeros (bits
, num
[numsize
- 1]);
1230 bits
= numsize
* BITS_PER_MP_LIMB
- bits
;
1232 /* Now we know the exponent of the number in base two.
1233 Check it against the maximum possible exponent. */
1234 if (__glibc_unlikely (bits
> MAX_EXP
))
1235 return overflow_value (negative
);
1237 /* We have already the first BITS bits of the result. Together with
1238 the information whether more non-zero bits follow this is enough
1239 to determine the result. */
1240 if (bits
> MANT_DIG
)
1243 const mp_size_t least_idx
= (bits
- MANT_DIG
) / BITS_PER_MP_LIMB
;
1244 const mp_size_t least_bit
= (bits
- MANT_DIG
) % BITS_PER_MP_LIMB
;
1245 const mp_size_t round_idx
= least_bit
== 0 ? least_idx
- 1
1247 const mp_size_t round_bit
= least_bit
== 0 ? BITS_PER_MP_LIMB
- 1
1251 memcpy (retval
, &num
[least_idx
],
1252 RETURN_LIMB_SIZE
* sizeof (mp_limb_t
));
1255 for (i
= least_idx
; i
< numsize
- 1; ++i
)
1256 retval
[i
- least_idx
] = (num
[i
] >> least_bit
)
1258 << (BITS_PER_MP_LIMB
- least_bit
));
1259 if (i
- least_idx
< RETURN_LIMB_SIZE
)
1260 retval
[RETURN_LIMB_SIZE
- 1] = num
[i
] >> least_bit
;
1263 /* Check whether any limb beside the ones in RETVAL are non-zero. */
1264 for (i
= 0; num
[i
] == 0; ++i
)
1267 return round_and_return (retval
, bits
- 1, negative
,
1268 num
[round_idx
], round_bit
,
1269 int_no
< dig_no
|| i
< round_idx
);
1272 else if (dig_no
== int_no
)
1274 const mp_size_t target_bit
= (MANT_DIG
- 1) % BITS_PER_MP_LIMB
;
1275 const mp_size_t is_bit
= (bits
- 1) % BITS_PER_MP_LIMB
;
1277 if (target_bit
== is_bit
)
1279 memcpy (&retval
[RETURN_LIMB_SIZE
- numsize
], num
,
1280 numsize
* sizeof (mp_limb_t
));
1281 /* FIXME: the following loop can be avoided if we assume a
1282 maximal MANT_DIG value. */
1283 MPN_ZERO (retval
, RETURN_LIMB_SIZE
- numsize
);
1285 else if (target_bit
> is_bit
)
1287 (void) __mpn_lshift (&retval
[RETURN_LIMB_SIZE
- numsize
],
1288 num
, numsize
, target_bit
- is_bit
);
1289 /* FIXME: the following loop can be avoided if we assume a
1290 maximal MANT_DIG value. */
1291 MPN_ZERO (retval
, RETURN_LIMB_SIZE
- numsize
);
1296 assert (numsize
< RETURN_LIMB_SIZE
);
1298 cy
= __mpn_rshift (&retval
[RETURN_LIMB_SIZE
- numsize
],
1299 num
, numsize
, is_bit
- target_bit
);
1300 retval
[RETURN_LIMB_SIZE
- numsize
- 1] = cy
;
1301 /* FIXME: the following loop can be avoided if we assume a
1302 maximal MANT_DIG value. */
1303 MPN_ZERO (retval
, RETURN_LIMB_SIZE
- numsize
- 1);
1306 return round_and_return (retval
, bits
- 1, negative
, 0, 0, 0);
1310 /* Store the bits we already have. */
1311 memcpy (retval
, num
, numsize
* sizeof (mp_limb_t
));
1312 #if RETURN_LIMB_SIZE > 1
1313 if (numsize
< RETURN_LIMB_SIZE
)
1314 # if RETURN_LIMB_SIZE == 2
1315 retval
[numsize
] = 0;
1317 MPN_ZERO (retval
+ numsize
, RETURN_LIMB_SIZE
- numsize
);
1322 /* We have to compute at least some of the fractional digits. */
1324 /* We construct a fraction and the result of the division gives us
1325 the needed digits. The denominator is 1.0 multiplied by the
1326 exponent of the lowest digit; i.e. 0.123 gives 123 / 1000 and
1327 123e-6 gives 123 / 1000000. */
1332 int need_frac_digits
;
1334 mp_limb_t
*psrc
= den
;
1335 mp_limb_t
*pdest
= num
;
1336 const struct mp_power
*ttab
= &_fpioconst_pow10
[0];
1338 assert (dig_no
> int_no
1340 && exponent
>= MIN_10_EXP
- (DIG
+ 2));
1342 /* We need to compute MANT_DIG - BITS fractional bits that lie
1343 within the mantissa of the result, the following bit for
1344 rounding, and to know whether any subsequent bit is 0.
1345 Computing a bit with value 2^-n means looking at n digits after
1346 the decimal point. */
1349 /* The bits required are those immediately after the point. */
1350 assert (int_no
> 0 && exponent
== 0);
1351 need_frac_digits
= 1 + MANT_DIG
- bits
;
1355 /* The number is in the form .123eEXPONENT. */
1356 assert (int_no
== 0 && *startp
!= L_('0'));
1357 /* The number is at least 10^(EXPONENT-1), and 10^3 <
1359 int neg_exp_2
= ((1 - exponent
) * 10) / 3 + 1;
1360 /* The number is at least 2^-NEG_EXP_2. We need up to
1361 MANT_DIG bits following that bit. */
1362 need_frac_digits
= neg_exp_2
+ MANT_DIG
;
1363 /* However, we never need bits beyond 1/4 ulp of the smallest
1364 representable value. (That 1/4 ulp bit is only needed to
1365 determine tinyness on machines where tinyness is determined
1367 if (need_frac_digits
> MANT_DIG
- MIN_EXP
+ 2)
1368 need_frac_digits
= MANT_DIG
- MIN_EXP
+ 2;
1369 /* At this point, NEED_FRAC_DIGITS is the total number of
1370 digits needed after the point, but some of those may be
1372 need_frac_digits
+= exponent
;
1373 /* Any cases underflowing enough that none of the fractional
1374 digits are needed should have been caught earlier (such
1375 cases are on the order of 10^-n or smaller where 2^-n is
1376 the least subnormal). */
1377 assert (need_frac_digits
> 0);
1380 if (need_frac_digits
> (intmax_t) dig_no
- (intmax_t) int_no
)
1381 need_frac_digits
= (intmax_t) dig_no
- (intmax_t) int_no
;
1383 if ((intmax_t) dig_no
> (intmax_t) int_no
+ need_frac_digits
)
1385 dig_no
= int_no
+ need_frac_digits
;
1391 neg_exp
= (intmax_t) dig_no
- (intmax_t) int_no
- exponent
;
1393 /* Construct the denominator. */
1398 if ((neg_exp
& expbit
) != 0)
1405 densize
= ttab
->arraysize
- _FPIO_CONST_OFFSET
;
1406 memcpy (psrc
, &__tens
[ttab
->arrayoff
+ _FPIO_CONST_OFFSET
],
1407 densize
* sizeof (mp_limb_t
));
1411 cy
= __mpn_mul (pdest
, &__tens
[ttab
->arrayoff
1412 + _FPIO_CONST_OFFSET
],
1413 ttab
->arraysize
- _FPIO_CONST_OFFSET
,
1415 densize
+= ttab
->arraysize
- _FPIO_CONST_OFFSET
;
1418 (void) SWAP (psrc
, pdest
);
1424 while (neg_exp
!= 0);
1427 memcpy (den
, num
, densize
* sizeof (mp_limb_t
));
1429 /* Read the fractional digits from the string. */
1430 (void) str_to_mpn (startp
, dig_no
- int_no
, num
, &numsize
, &exponent
1431 #ifndef USE_WIDE_CHAR
1432 , decimal
, decimal_len
, thousands
1436 /* We now have to shift both numbers so that the highest bit in the
1437 denominator is set. In the same process we copy the numerator to
1438 a high place in the array so that the division constructs the wanted
1439 digits. This is done by a "quasi fix point" number representation.
1441 num: ddddddddddd . 0000000000000000000000
1443 den: ddddddddddd n >= m
1447 count_leading_zeros (cnt
, den
[densize
- 1]);
1451 /* Don't call `mpn_shift' with a count of zero since the specification
1452 does not allow this. */
1453 (void) __mpn_lshift (den
, den
, densize
, cnt
);
1454 cy
= __mpn_lshift (num
, num
, numsize
, cnt
);
1456 num
[numsize
++] = cy
;
1459 /* Now we are ready for the division. But it is not necessary to
1460 do a full multi-precision division because we only need a small
1461 number of bits for the result. So we do not use __mpn_divmod
1462 here but instead do the division here by hand and stop whenever
1463 the needed number of bits is reached. The code itself comes
1464 from the GNU MP Library by Torbj\"orn Granlund. */
1472 mp_limb_t d
, n
, quot
;
1477 assert (numsize
== 1 && n
< d
);
1481 udiv_qrnnd (quot
, n
, n
, 0, d
);
1488 cnt = BITS_PER_MP_LIMB; \
1490 count_leading_zeros (cnt, quot); \
1492 if (BITS_PER_MP_LIMB - cnt > MANT_DIG) \
1494 used = MANT_DIG + cnt; \
1495 retval[0] = quot >> (BITS_PER_MP_LIMB - used); \
1496 bits = MANT_DIG + 1; \
1500 /* Note that we only clear the second element. */ \
1501 /* The conditional is determined at compile time. */ \
1502 if (RETURN_LIMB_SIZE > 1) \
1508 else if (bits + BITS_PER_MP_LIMB <= MANT_DIG) \
1509 __mpn_lshift_1 (retval, RETURN_LIMB_SIZE, BITS_PER_MP_LIMB, \
1513 used = MANT_DIG - bits; \
1515 __mpn_lshift_1 (retval, RETURN_LIMB_SIZE, used, quot); \
1517 bits += BITS_PER_MP_LIMB
1521 while (bits
<= MANT_DIG
);
1523 return round_and_return (retval
, exponent
- 1, negative
,
1524 quot
, BITS_PER_MP_LIMB
- 1 - used
,
1525 more_bits
|| n
!= 0);
1529 mp_limb_t d0
, d1
, n0
, n1
;
1536 if (numsize
< densize
)
1540 /* The numerator of the number occupies fewer bits than
1541 the denominator but the one limb is bigger than the
1542 high limb of the numerator. */
1549 exponent
-= BITS_PER_MP_LIMB
;
1552 if (bits
+ BITS_PER_MP_LIMB
<= MANT_DIG
)
1553 __mpn_lshift_1 (retval
, RETURN_LIMB_SIZE
,
1554 BITS_PER_MP_LIMB
, 0);
1557 used
= MANT_DIG
- bits
;
1559 __mpn_lshift_1 (retval
, RETURN_LIMB_SIZE
, used
, 0);
1561 bits
+= BITS_PER_MP_LIMB
;
1573 while (bits
<= MANT_DIG
)
1579 /* QUOT should be either 111..111 or 111..110. We need
1580 special treatment of this rare case as normal division
1581 would give overflow. */
1582 quot
= ~(mp_limb_t
) 0;
1585 if (r
< d1
) /* Carry in the addition? */
1587 add_ssaaaa (n1
, n0
, r
- d0
, 0, 0, d0
);
1590 n1
= d0
- (d0
!= 0);
1595 udiv_qrnnd (quot
, r
, n1
, n0
, d1
);
1596 umul_ppmm (n1
, n0
, d0
, quot
);
1600 if (n1
> r
|| (n1
== r
&& n0
> 0))
1602 /* The estimated QUOT was too large. */
1605 sub_ddmmss (n1
, n0
, n1
, n0
, 0, d0
);
1607 if (r
>= d1
) /* If not carry, test QUOT again. */
1610 sub_ddmmss (n1
, n0
, r
, 0, n1
, n0
);
1616 return round_and_return (retval
, exponent
- 1, negative
,
1617 quot
, BITS_PER_MP_LIMB
- 1 - used
,
1618 more_bits
|| n1
!= 0 || n0
!= 0);
1623 mp_limb_t cy
, dX
, d1
, n0
, n1
;
1627 dX
= den
[densize
- 1];
1628 d1
= den
[densize
- 2];
1630 /* The division does not work if the upper limb of the two-limb
1631 numerator is greater than the denominator. */
1632 if (__mpn_cmp (num
, &den
[densize
- numsize
], numsize
) > 0)
1635 if (numsize
< densize
)
1637 mp_size_t empty
= densize
- numsize
;
1641 exponent
-= empty
* BITS_PER_MP_LIMB
;
1644 if (bits
+ empty
* BITS_PER_MP_LIMB
<= MANT_DIG
)
1646 /* We make a difference here because the compiler
1647 cannot optimize the `else' case that good and
1648 this reflects all currently used FLOAT types
1649 and GMP implementations. */
1650 #if RETURN_LIMB_SIZE <= 2
1651 assert (empty
== 1);
1652 __mpn_lshift_1 (retval
, RETURN_LIMB_SIZE
,
1653 BITS_PER_MP_LIMB
, 0);
1655 for (i
= RETURN_LIMB_SIZE
- 1; i
>= empty
; --i
)
1656 retval
[i
] = retval
[i
- empty
];
1663 used
= MANT_DIG
- bits
;
1664 if (used
>= BITS_PER_MP_LIMB
)
1667 (void) __mpn_lshift (&retval
[used
1668 / BITS_PER_MP_LIMB
],
1671 - used
/ BITS_PER_MP_LIMB
),
1672 used
% BITS_PER_MP_LIMB
);
1673 for (i
= used
/ BITS_PER_MP_LIMB
- 1; i
>= 0; --i
)
1677 __mpn_lshift_1 (retval
, RETURN_LIMB_SIZE
, used
, 0);
1679 bits
+= empty
* BITS_PER_MP_LIMB
;
1681 for (i
= numsize
; i
> 0; --i
)
1682 num
[i
+ empty
] = num
[i
- 1];
1683 MPN_ZERO (num
, empty
+ 1);
1688 assert (numsize
== densize
);
1689 for (i
= numsize
; i
> 0; --i
)
1690 num
[i
] = num
[i
- 1];
1697 while (bits
<= MANT_DIG
)
1700 /* This might over-estimate QUOT, but it's probably not
1701 worth the extra code here to find out. */
1702 quot
= ~(mp_limb_t
) 0;
1707 udiv_qrnnd (quot
, r
, n0
, num
[densize
- 1], dX
);
1708 umul_ppmm (n1
, n0
, d1
, quot
);
1710 while (n1
> r
|| (n1
== r
&& n0
> num
[densize
- 2]))
1714 if (r
< dX
) /* I.e. "carry in previous addition?" */
1721 /* Possible optimization: We already have (q * n0) and (1 * n1)
1722 after the calculation of QUOT. Taking advantage of this, we
1723 could make this loop make two iterations less. */
1725 cy
= __mpn_submul_1 (num
, den
, densize
+ 1, quot
);
1727 if (num
[densize
] != cy
)
1729 cy
= __mpn_add_n (num
, num
, den
, densize
);
1733 n0
= num
[densize
] = num
[densize
- 1];
1734 for (i
= densize
- 1; i
> 0; --i
)
1735 num
[i
] = num
[i
- 1];
1741 for (i
= densize
; i
>= 0 && num
[i
] == 0; --i
)
1743 return round_and_return (retval
, exponent
- 1, negative
,
1744 quot
, BITS_PER_MP_LIMB
- 1 - used
,
1745 more_bits
|| i
>= 0);
1752 #if defined _LIBC && !defined USE_WIDE_CHAR
1753 libc_hidden_def (____STRTOF_INTERNAL
)
1756 /* External user entry point. */
1759 #ifdef weak_function
1762 __STRTOF (const STRING_TYPE
*nptr
, STRING_TYPE
**endptr
, __locale_t loc
)
1764 return ____STRTOF_INTERNAL (nptr
, endptr
, 0, loc
);
1767 libc_hidden_def (__STRTOF
)
1768 libc_hidden_ver (__STRTOF
, STRTOF
)
1770 weak_alias (__STRTOF
, STRTOF
)
1772 #ifdef LONG_DOUBLE_COMPAT
1773 # if LONG_DOUBLE_COMPAT(libc, GLIBC_2_1)
1774 # ifdef USE_WIDE_CHAR
1775 compat_symbol (libc
, __wcstod_l
, __wcstold_l
, GLIBC_2_1
);
1777 compat_symbol (libc
, __strtod_l
, __strtold_l
, GLIBC_2_1
);
1780 # if LONG_DOUBLE_COMPAT(libc, GLIBC_2_3)
1781 # ifdef USE_WIDE_CHAR
1782 compat_symbol (libc
, wcstod_l
, wcstold_l
, GLIBC_2_3
);
1784 compat_symbol (libc
, strtod_l
, strtold_l
, GLIBC_2_3
);