1 /* Convert string representing a number to float value, using given locale.
2 Copyright (C) 1997-2012 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
);
23 extern unsigned long long int ____strtoull_l_internal (const char *, char **,
24 int, int, __locale_t
);
26 /* Configuration part. These macros are defined by `strtold.c',
27 `strtof.c', `wcstod.c', `wcstold.c', and `wcstof.c' to produce the
28 `long double' and `float' versions of the reader. */
30 # include <math_ldbl_opt.h>
34 # define STRTOF wcstod_l
35 # define __STRTOF __wcstod_l
37 # define STRTOF strtod_l
38 # define __STRTOF __strtod_l
40 # define MPN2FLOAT __mpn_construct_double
41 # define FLOAT_HUGE_VAL HUGE_VAL
42 # define SET_MANTISSA(flt, mant) \
43 do { union ieee754_double u; \
45 if ((mant & 0xfffffffffffffULL) == 0) \
46 mant = 0x8000000000000ULL; \
47 u.ieee.mantissa0 = ((mant) >> 32) & 0xfffff; \
48 u.ieee.mantissa1 = (mant) & 0xffffffff; \
52 /* End of configuration part. */
58 #include "../locale/localeinfo.h"
64 /* The gmp headers need some configuration frobs. */
67 /* Include gmp-mparam.h first, such that definitions of _SHORT_LIMB
68 and _LONG_LONG_LIMB in it can take effect into gmp.h. */
69 #include <gmp-mparam.h>
73 #include "fpioconst.h"
79 /* We use this code for the extended locale handling where the
80 function gets as an additional argument the locale which has to be
81 used. To access the values we have to redefine the _NL_CURRENT and
82 _NL_CURRENT_WORD macros. */
84 #define _NL_CURRENT(category, item) \
85 (current->values[_NL_ITEM_INDEX (item)].string)
86 #undef _NL_CURRENT_WORD
87 #define _NL_CURRENT_WORD(category, item) \
88 ((uint32_t) current->values[_NL_ITEM_INDEX (item)].word)
90 #if defined _LIBC || defined HAVE_WCHAR_H
96 # define STRING_TYPE wchar_t
97 # define CHAR_TYPE wint_t
99 # define ISSPACE(Ch) __iswspace_l ((Ch), loc)
100 # define ISDIGIT(Ch) __iswdigit_l ((Ch), loc)
101 # define ISXDIGIT(Ch) __iswxdigit_l ((Ch), loc)
102 # define TOLOWER(Ch) __towlower_l ((Ch), loc)
103 # define TOLOWER_C(Ch) __towlower_l ((Ch), _nl_C_locobj_ptr)
104 # define STRNCASECMP(S1, S2, N) \
105 __wcsncasecmp_l ((S1), (S2), (N), _nl_C_locobj_ptr)
106 # define STRTOULL(S, E, B) ____wcstoull_l_internal ((S), (E), (B), 0, loc)
108 # define STRING_TYPE char
109 # define CHAR_TYPE char
111 # define ISSPACE(Ch) __isspace_l ((Ch), loc)
112 # define ISDIGIT(Ch) __isdigit_l ((Ch), loc)
113 # define ISXDIGIT(Ch) __isxdigit_l ((Ch), loc)
114 # define TOLOWER(Ch) __tolower_l ((Ch), loc)
115 # define TOLOWER_C(Ch) __tolower_l ((Ch), _nl_C_locobj_ptr)
116 # define STRNCASECMP(S1, S2, N) \
117 __strncasecmp_l ((S1), (S2), (N), _nl_C_locobj_ptr)
118 # define STRTOULL(S, E, B) ____strtoull_l_internal ((S), (E), (B), 0, loc)
122 /* Constants we need from float.h; select the set for the FLOAT precision. */
123 #define MANT_DIG PASTE(FLT,_MANT_DIG)
124 #define DIG PASTE(FLT,_DIG)
125 #define MAX_EXP PASTE(FLT,_MAX_EXP)
126 #define MIN_EXP PASTE(FLT,_MIN_EXP)
127 #define MAX_10_EXP PASTE(FLT,_MAX_10_EXP)
128 #define MIN_10_EXP PASTE(FLT,_MIN_10_EXP)
130 /* Extra macros required to get FLT expanded before the pasting. */
131 #define PASTE(a,b) PASTE1(a,b)
132 #define PASTE1(a,b) a##b
134 /* Function to construct a floating point number from an MP integer
135 containing the fraction bits, a base 2 exponent, and a sign flag. */
136 extern FLOAT
MPN2FLOAT (mp_srcptr mpn
, int exponent
, int negative
);
138 /* Definitions according to limb size used. */
139 #if BITS_PER_MP_LIMB == 32
140 # define MAX_DIG_PER_LIMB 9
141 # define MAX_FAC_PER_LIMB 1000000000UL
142 #elif BITS_PER_MP_LIMB == 64
143 # define MAX_DIG_PER_LIMB 19
144 # define MAX_FAC_PER_LIMB 10000000000000000000ULL
146 # error "mp_limb_t size " BITS_PER_MP_LIMB "not accounted for"
149 extern const mp_limb_t _tens_in_limb
[MAX_DIG_PER_LIMB
+ 1];
152 #define howmany(x,y) (((x)+((y)-1))/(y))
154 #define SWAP(x, y) ({ typeof(x) _tmp = x; x = y; y = _tmp; })
156 #define NDIG (MAX_10_EXP - MIN_10_EXP + 2 * MANT_DIG)
157 #define HEXNDIG ((MAX_EXP - MIN_EXP + 7) / 8 + 2 * MANT_DIG)
158 #define RETURN_LIMB_SIZE howmany (MANT_DIG, BITS_PER_MP_LIMB)
160 #define RETURN(val,end) \
161 do { if (endptr != NULL) *endptr = (STRING_TYPE *) (end); \
162 return val; } while (0)
164 /* Maximum size necessary for mpn integers to hold floating point numbers. */
165 #define MPNSIZE (howmany (MAX_EXP + 2 * MANT_DIG, BITS_PER_MP_LIMB) \
167 /* Declare an mpn integer variable that big. */
168 #define MPN_VAR(name) mp_limb_t name[MPNSIZE]; mp_size_t name##size
169 /* Copy an mpn integer value. */
170 #define MPN_ASSIGN(dst, src) \
171 memcpy (dst, src, (dst##size = src##size) * sizeof (mp_limb_t))
174 /* Return a floating point number of the needed type according to the given
175 multi-precision number after possible rounding. */
177 round_and_return (mp_limb_t
*retval
, int exponent
, int negative
,
178 mp_limb_t round_limb
, mp_size_t round_bit
, int more_bits
)
180 if (exponent
< MIN_EXP
- 1)
182 mp_size_t shift
= MIN_EXP
- 1 - exponent
;
184 if (shift
> MANT_DIG
)
186 __set_errno (ERANGE
);
190 more_bits
|= (round_limb
& ((((mp_limb_t
) 1) << round_bit
) - 1)) != 0;
191 if (shift
== MANT_DIG
)
192 /* This is a special case to handle the very seldom case where
193 the mantissa will be empty after the shift. */
197 round_limb
= retval
[RETURN_LIMB_SIZE
- 1];
198 round_bit
= (MANT_DIG
- 1) % BITS_PER_MP_LIMB
;
199 for (i
= 0; i
< RETURN_LIMB_SIZE
; ++i
)
200 more_bits
|= retval
[i
] != 0;
201 MPN_ZERO (retval
, RETURN_LIMB_SIZE
);
203 else if (shift
>= BITS_PER_MP_LIMB
)
207 round_limb
= retval
[(shift
- 1) / BITS_PER_MP_LIMB
];
208 round_bit
= (shift
- 1) % BITS_PER_MP_LIMB
;
209 for (i
= 0; i
< (shift
- 1) / BITS_PER_MP_LIMB
; ++i
)
210 more_bits
|= retval
[i
] != 0;
211 more_bits
|= ((round_limb
& ((((mp_limb_t
) 1) << round_bit
) - 1))
214 (void) __mpn_rshift (retval
, &retval
[shift
/ BITS_PER_MP_LIMB
],
215 RETURN_LIMB_SIZE
- (shift
/ BITS_PER_MP_LIMB
),
216 shift
% BITS_PER_MP_LIMB
);
217 MPN_ZERO (&retval
[RETURN_LIMB_SIZE
- (shift
/ BITS_PER_MP_LIMB
)],
218 shift
/ BITS_PER_MP_LIMB
);
222 round_limb
= retval
[0];
223 round_bit
= shift
- 1;
224 (void) __mpn_rshift (retval
, retval
, RETURN_LIMB_SIZE
, shift
);
226 /* This is a hook for the m68k long double format, where the
227 exponent bias is the same for normalized and denormalized
230 # define DENORM_EXP (MIN_EXP - 2)
232 exponent
= DENORM_EXP
;
233 __set_errno (ERANGE
);
236 if ((round_limb
& (((mp_limb_t
) 1) << round_bit
)) != 0
237 && (more_bits
|| (retval
[0] & 1) != 0
238 || (round_limb
& ((((mp_limb_t
) 1) << round_bit
) - 1)) != 0))
240 mp_limb_t cy
= __mpn_add_1 (retval
, retval
, RETURN_LIMB_SIZE
, 1);
242 if (((MANT_DIG
% BITS_PER_MP_LIMB
) == 0 && cy
) ||
243 ((MANT_DIG
% BITS_PER_MP_LIMB
) != 0 &&
244 (retval
[RETURN_LIMB_SIZE
- 1]
245 & (((mp_limb_t
) 1) << (MANT_DIG
% BITS_PER_MP_LIMB
))) != 0))
248 (void) __mpn_rshift (retval
, retval
, RETURN_LIMB_SIZE
, 1);
249 retval
[RETURN_LIMB_SIZE
- 1]
250 |= ((mp_limb_t
) 1) << ((MANT_DIG
- 1) % BITS_PER_MP_LIMB
);
252 else if (exponent
== DENORM_EXP
253 && (retval
[RETURN_LIMB_SIZE
- 1]
254 & (((mp_limb_t
) 1) << ((MANT_DIG
- 1) % BITS_PER_MP_LIMB
)))
256 /* The number was denormalized but now normalized. */
257 exponent
= MIN_EXP
- 1;
260 if (exponent
> MAX_EXP
)
261 return negative
? -FLOAT_HUGE_VAL
: FLOAT_HUGE_VAL
;
263 return MPN2FLOAT (retval
, exponent
, negative
);
267 /* Read a multi-precision integer starting at STR with exactly DIGCNT digits
268 into N. Return the size of the number limbs in NSIZE at the first
269 character od the string that is not part of the integer as the function
270 value. If the EXPONENT is small enough to be taken as an additional
271 factor for the resulting number (see code) multiply by it. */
272 static const STRING_TYPE
*
273 str_to_mpn (const STRING_TYPE
*str
, int digcnt
, mp_limb_t
*n
, mp_size_t
*nsize
,
275 #ifndef USE_WIDE_CHAR
276 , const char *decimal
, size_t decimal_len
, const char *thousands
281 /* Number of digits for actual limb. */
290 if (cnt
== MAX_DIG_PER_LIMB
)
300 cy
= __mpn_mul_1 (n
, n
, *nsize
, MAX_FAC_PER_LIMB
);
301 cy
+= __mpn_add_1 (n
, n
, *nsize
, low
);
312 /* There might be thousands separators or radix characters in
313 the string. But these all can be ignored because we know the
314 format of the number is correct and we have an exact number
315 of characters to read. */
317 if (*str
< L
'0' || *str
> L
'9')
320 if (*str
< '0' || *str
> '9')
323 if (thousands
!= NULL
&& *str
== *thousands
324 && ({ for (inner
= 1; thousands
[inner
] != '\0'; ++inner
)
325 if (thousands
[inner
] != str
[inner
])
327 thousands
[inner
] == '\0'; }))
333 low
= low
* 10 + *str
++ - L_('0');
336 while (--digcnt
> 0);
338 if (*exponent
> 0 && cnt
+ *exponent
<= MAX_DIG_PER_LIMB
)
340 low
*= _tens_in_limb
[*exponent
];
341 start
= _tens_in_limb
[cnt
+ *exponent
];
345 start
= _tens_in_limb
[cnt
];
355 cy
= __mpn_mul_1 (n
, n
, *nsize
, start
);
356 cy
+= __mpn_add_1 (n
, n
, *nsize
, low
);
365 /* Shift {PTR, SIZE} COUNT bits to the left, and fill the vacated bits
366 with the COUNT most significant bits of LIMB.
368 Tege doesn't like this function so I have to write it here myself. :)
371 __attribute ((always_inline
))
372 __mpn_lshift_1 (mp_limb_t
*ptr
, mp_size_t size
, unsigned int count
,
375 if (__builtin_constant_p (count
) && count
== BITS_PER_MP_LIMB
)
377 /* Optimize the case of shifting by exactly a word:
378 just copy words, with no actual bit-shifting. */
380 for (i
= size
- 1; i
> 0; --i
)
386 (void) __mpn_lshift (ptr
, ptr
, size
, count
);
387 ptr
[0] |= limb
>> (BITS_PER_MP_LIMB
- count
);
392 #define INTERNAL(x) INTERNAL1(x)
393 #define INTERNAL1(x) __##x##_internal
394 #ifndef ____STRTOF_INTERNAL
395 # define ____STRTOF_INTERNAL INTERNAL (__STRTOF)
398 /* This file defines a function to check for correct grouping. */
399 #include "grouping.h"
402 /* Return a floating point number with the value of the given string NPTR.
403 Set *ENDPTR to the character after the last used one. If the number is
404 smaller than the smallest representable number, set `errno' to ERANGE and
405 return 0.0. If the number is too big to be represented, set `errno' to
406 ERANGE and return HUGE_VAL with the appropriate sign. */
408 ____STRTOF_INTERNAL (nptr
, endptr
, group
, loc
)
409 const STRING_TYPE
*nptr
;
410 STRING_TYPE
**endptr
;
414 int negative
; /* The sign of the number. */
415 MPN_VAR (num
); /* MP representation of the number. */
416 int exponent
; /* Exponent of the number. */
418 /* Numbers starting `0X' or `0x' have to be processed with base 16. */
421 /* When we have to compute fractional digits we form a fraction with a
422 second multi-precision number (and we sometimes need a second for
423 temporary results). */
426 /* Representation for the return value. */
427 mp_limb_t retval
[RETURN_LIMB_SIZE
];
428 /* Number of bits currently in result value. */
431 /* Running pointer after the last character processed in the string. */
432 const STRING_TYPE
*cp
, *tp
;
433 /* Start of significant part of the number. */
434 const STRING_TYPE
*startp
, *start_of_digits
;
435 /* Points at the character following the integer and fractional digits. */
436 const STRING_TYPE
*expp
;
437 /* Total number of digit and number of digits in integer part. */
438 int dig_no
, int_no
, lead_zero
;
439 /* Contains the last character read. */
442 /* We should get wint_t from <stddef.h>, but not all GCC versions define it
443 there. So define it ourselves if it remains undefined. */
445 typedef unsigned int wint_t;
447 /* The radix character of the current locale. */
454 /* The thousands character of the current locale. */
456 wchar_t thousands
= L
'\0';
458 const char *thousands
= NULL
;
460 /* The numeric grouping specification of the current locale,
461 in the format described in <locale.h>. */
462 const char *grouping
;
463 /* Used in several places. */
466 struct __locale_data
*current
= loc
->__locales
[LC_NUMERIC
];
468 if (__builtin_expect (group
, 0))
470 grouping
= _NL_CURRENT (LC_NUMERIC
, GROUPING
);
471 if (*grouping
<= 0 || *grouping
== CHAR_MAX
)
475 /* Figure out the thousands separator character. */
477 thousands
= _NL_CURRENT_WORD (LC_NUMERIC
,
478 _NL_NUMERIC_THOUSANDS_SEP_WC
);
479 if (thousands
== L
'\0')
482 thousands
= _NL_CURRENT (LC_NUMERIC
, THOUSANDS_SEP
);
483 if (*thousands
== '\0')
494 /* Find the locale's decimal point character. */
496 decimal
= _NL_CURRENT_WORD (LC_NUMERIC
, _NL_NUMERIC_DECIMAL_POINT_WC
);
497 assert (decimal
!= L
'\0');
498 # define decimal_len 1
500 decimal
= _NL_CURRENT (LC_NUMERIC
, DECIMAL_POINT
);
501 decimal_len
= strlen (decimal
);
502 assert (decimal_len
> 0);
505 /* Prepare number representation. */
510 /* Parse string to get maximal legal prefix. We need the number of
511 characters of the integer part, the fractional part and the exponent. */
513 /* Ignore leading white space. */
518 /* Get sign of the result. */
524 else if (c
== L_('+'))
527 /* Return 0.0 if no legal string is found.
528 No character is used even if a sign was found. */
530 if (c
== (wint_t) decimal
531 && (wint_t) cp
[1] >= L
'0' && (wint_t) cp
[1] <= L
'9')
533 /* We accept it. This funny construct is here only to indent
534 the code correctly. */
537 for (cnt
= 0; decimal
[cnt
] != '\0'; ++cnt
)
538 if (cp
[cnt
] != decimal
[cnt
])
540 if (decimal
[cnt
] == '\0' && cp
[cnt
] >= '0' && cp
[cnt
] <= '9')
542 /* We accept it. This funny construct is here only to indent
543 the code correctly. */
546 else if (c
< L_('0') || c
> L_('9'))
548 /* Check for `INF' or `INFINITY'. */
549 CHAR_TYPE lowc
= TOLOWER_C (c
);
551 if (lowc
== L_('i') && STRNCASECMP (cp
, L_("inf"), 3) == 0)
553 /* Return +/- infinity. */
555 *endptr
= (STRING_TYPE
*)
556 (cp
+ (STRNCASECMP (cp
+ 3, L_("inity"), 5) == 0
559 return negative
? -FLOAT_HUGE_VAL
: FLOAT_HUGE_VAL
;
562 if (lowc
== L_('n') && STRNCASECMP (cp
, L_("nan"), 3) == 0)
569 /* Match `(n-char-sequence-digit)'. */
572 const STRING_TYPE
*startp
= cp
;
575 while ((*cp
>= L_('0') && *cp
<= L_('9'))
576 || ({ CHAR_TYPE lo
= TOLOWER (*cp
);
577 lo
>= L_('a') && lo
<= L_('z'); })
581 /* The closing brace is missing. Only match the NAN
586 /* This is a system-dependent way to specify the
587 bitmask used for the NaN. We expect it to be
588 a number which is put in the mantissa of the
591 unsigned long long int mant
;
593 mant
= STRTOULL (startp
+ 1, &endp
, 0);
595 SET_MANTISSA (retval
, mant
);
597 /* Consume the closing brace. */
603 *endptr
= (STRING_TYPE
*) cp
;
608 /* It is really a text we do not recognize. */
612 /* First look whether we are faced with a hexadecimal number. */
613 if (c
== L_('0') && TOLOWER (cp
[1]) == L_('x'))
615 /* Okay, it is a hexa-decimal number. Remember this and skip
616 the characters. BTW: hexadecimal numbers must not be
624 /* Record the start of the digits, in case we will check their grouping. */
625 start_of_digits
= startp
= cp
;
627 /* Ignore leading zeroes. This helps us to avoid useless computations. */
629 while (c
== L
'0' || ((wint_t) thousands
!= L
'\0' && c
== (wint_t) thousands
))
632 if (__builtin_expect (thousands
== NULL
, 1))
637 /* We also have the multibyte thousands string. */
642 for (cnt
= 0; thousands
[cnt
] != '\0'; ++cnt
)
643 if (thousands
[cnt
] != cp
[cnt
])
645 if (thousands
[cnt
] != '\0')
654 /* If no other digit but a '0' is found the result is 0.0.
655 Return current read pointer. */
656 CHAR_TYPE lowc
= TOLOWER (c
);
657 if (!((c
>= L_('0') && c
<= L_('9'))
658 || (base
== 16 && lowc
>= L_('a') && lowc
<= L_('f'))
661 c
== (wint_t) decimal
663 ({ for (cnt
= 0; decimal
[cnt
] != '\0'; ++cnt
)
664 if (decimal
[cnt
] != cp
[cnt
])
666 decimal
[cnt
] == '\0'; })
668 /* '0x.' alone is not a valid hexadecimal number.
669 '.' alone is not valid either, but that has been checked
672 || cp
!= start_of_digits
673 || (cp
[decimal_len
] >= L_('0') && cp
[decimal_len
] <= L_('9'))
674 || ({ CHAR_TYPE lo
= TOLOWER (cp
[decimal_len
]);
675 lo
>= L_('a') && lo
<= L_('f'); })))
676 || (base
== 16 && (cp
!= start_of_digits
678 || (base
!= 16 && lowc
== L_('e'))))
681 tp
= __correctly_grouped_prefixwc (start_of_digits
, cp
, thousands
,
684 tp
= __correctly_grouped_prefixmb (start_of_digits
, cp
, thousands
,
687 /* If TP is at the start of the digits, there was no correctly
688 grouped prefix of the string; so no number found. */
689 RETURN (negative
? -0.0 : 0.0,
690 tp
== start_of_digits
? (base
== 16 ? cp
- 1 : nptr
) : tp
);
693 /* Remember first significant digit and read following characters until the
694 decimal point, exponent character or any non-FP number character. */
699 if ((c
>= L_('0') && c
<= L_('9'))
701 && ({ CHAR_TYPE lo
= TOLOWER (c
);
702 lo
>= L_('a') && lo
<= L_('f'); })))
707 if (__builtin_expect ((wint_t) thousands
== L
'\0', 1)
708 || c
!= (wint_t) thousands
)
709 /* Not a digit or separator: end of the integer part. */
712 if (__builtin_expect (thousands
== NULL
, 1))
716 for (cnt
= 0; thousands
[cnt
] != '\0'; ++cnt
)
717 if (thousands
[cnt
] != cp
[cnt
])
719 if (thousands
[cnt
] != '\0')
728 if (__builtin_expect (grouping
!= NULL
, 0) && cp
> start_of_digits
)
730 /* Check the grouping of the digits. */
732 tp
= __correctly_grouped_prefixwc (start_of_digits
, cp
, thousands
,
735 tp
= __correctly_grouped_prefixmb (start_of_digits
, cp
, thousands
,
740 /* Less than the entire string was correctly grouped. */
742 if (tp
== start_of_digits
)
743 /* No valid group of numbers at all: no valid number. */
747 /* The number is validly grouped, but consists
748 only of zeroes. The whole value is zero. */
749 RETURN (negative
? -0.0 : 0.0, tp
);
751 /* Recompute DIG_NO so we won't read more digits than
752 are properly grouped. */
755 for (tp
= startp
; tp
< cp
; ++tp
)
756 if (*tp
>= L_('0') && *tp
<= L_('9'))
766 /* We have the number of digits in the integer part. Whether these
767 are all or any is really a fractional digit will be decided
770 lead_zero
= int_no
== 0 ? -1 : 0;
772 /* Read the fractional digits. A special case are the 'american
773 style' numbers like `16.' i.e. with decimal point but without
777 c
== (wint_t) decimal
779 ({ for (cnt
= 0; decimal
[cnt
] != '\0'; ++cnt
)
780 if (decimal
[cnt
] != cp
[cnt
])
782 decimal
[cnt
] == '\0'; })
788 while ((c
>= L_('0') && c
<= L_('9')) ||
789 (base
== 16 && ({ CHAR_TYPE lo
= TOLOWER (c
);
790 lo
>= L_('a') && lo
<= L_('f'); })))
792 if (c
!= L_('0') && lead_zero
== -1)
793 lead_zero
= dig_no
- int_no
;
799 /* Remember start of exponent (if any). */
804 if ((base
== 16 && lowc
== L_('p'))
805 || (base
!= 16 && lowc
== L_('e')))
807 int exp_negative
= 0;
815 else if (c
== L_('+'))
818 if (c
>= L_('0') && c
<= L_('9'))
822 /* Get the exponent limit. */
824 exp_limit
= (exp_negative
?
825 -MIN_EXP
+ MANT_DIG
+ 4 * int_no
:
826 MAX_EXP
- 4 * int_no
+ 4 * lead_zero
+ 3);
828 exp_limit
= (exp_negative
?
829 -MIN_10_EXP
+ MANT_DIG
+ int_no
:
830 MAX_10_EXP
- int_no
+ lead_zero
+ 1);
835 exponent
+= c
- L_('0');
837 if (__builtin_expect (exponent
> exp_limit
, 0))
838 /* The exponent is too large/small to represent a valid
843 /* We have to take care for special situation: a joker
844 might have written "0.0e100000" which is in fact
847 result
= negative
? -0.0 : 0.0;
850 /* Overflow or underflow. */
851 __set_errno (ERANGE
);
852 result
= (exp_negative
? (negative
? -0.0 : 0.0) :
853 negative
? -FLOAT_HUGE_VAL
: FLOAT_HUGE_VAL
);
856 /* Accept all following digits as part of the exponent. */
859 while (*cp
>= L_('0') && *cp
<= L_('9'));
867 while (c
>= L_('0') && c
<= L_('9'));
870 exponent
= -exponent
;
876 /* We don't want to have to work with trailing zeroes after the radix. */
879 while (expp
[-1] == L_('0'))
884 assert (dig_no
>= int_no
);
887 if (dig_no
== int_no
&& dig_no
> 0 && exponent
< 0)
890 while (! (base
== 16 ? ISXDIGIT (expp
[-1]) : ISDIGIT (expp
[-1])))
893 if (expp
[-1] != L_('0'))
899 exponent
+= base
== 16 ? 4 : 1;
901 while (dig_no
> 0 && exponent
< 0);
905 /* The whole string is parsed. Store the address of the next character. */
907 *endptr
= (STRING_TYPE
*) cp
;
910 return negative
? -0.0 : 0.0;
914 /* Find the decimal point */
916 while (*startp
!= decimal
)
921 if (*startp
== decimal
[0])
923 for (cnt
= 1; decimal
[cnt
] != '\0'; ++cnt
)
924 if (decimal
[cnt
] != startp
[cnt
])
926 if (decimal
[cnt
] == '\0')
932 startp
+= lead_zero
+ decimal_len
;
933 exponent
-= base
== 16 ? 4 * lead_zero
: lead_zero
;
937 /* If the BASE is 16 we can use a simpler algorithm. */
940 static const int nbits
[16] = { 0, 1, 2, 2, 3, 3, 3, 3,
941 4, 4, 4, 4, 4, 4, 4, 4 };
942 int idx
= (MANT_DIG
- 1) / BITS_PER_MP_LIMB
;
943 int pos
= (MANT_DIG
- 1) % BITS_PER_MP_LIMB
;
946 while (!ISXDIGIT (*startp
))
948 while (*startp
== L_('0'))
950 if (ISDIGIT (*startp
))
951 val
= *startp
++ - L_('0');
953 val
= 10 + TOLOWER (*startp
++) - L_('a');
955 /* We cannot have a leading zero. */
958 if (pos
+ 1 >= 4 || pos
+ 1 >= bits
)
960 /* We don't have to care for wrapping. This is the normal
961 case so we add the first clause in the `if' expression as
962 an optimization. It is a compile-time constant and so does
963 not cost anything. */
964 retval
[idx
] = val
<< (pos
- bits
+ 1);
969 retval
[idx
--] = val
>> (bits
- pos
- 1);
970 retval
[idx
] = val
<< (BITS_PER_MP_LIMB
- (bits
- pos
- 1));
971 pos
= BITS_PER_MP_LIMB
- 1 - (bits
- pos
- 1);
974 /* Adjust the exponent for the bits we are shifting in. */
975 exponent
+= bits
- 1 + (int_no
- 1) * 4;
977 while (--dig_no
> 0 && idx
>= 0)
979 if (!ISXDIGIT (*startp
))
980 startp
+= decimal_len
;
981 if (ISDIGIT (*startp
))
982 val
= *startp
++ - L_('0');
984 val
= 10 + TOLOWER (*startp
++) - L_('a');
988 retval
[idx
] |= val
<< (pos
- 4 + 1);
993 retval
[idx
--] |= val
>> (4 - pos
- 1);
994 val
<<= BITS_PER_MP_LIMB
- (4 - pos
- 1);
997 int rest_nonzero
= 0;
1000 if (*startp
!= L_('0'))
1007 return round_and_return (retval
, exponent
, negative
, val
,
1008 BITS_PER_MP_LIMB
- 1, rest_nonzero
);
1012 pos
= BITS_PER_MP_LIMB
- 1 - (4 - pos
- 1);
1016 /* We ran out of digits. */
1017 MPN_ZERO (retval
, idx
);
1019 return round_and_return (retval
, exponent
, negative
, 0, 0, 0);
1022 /* Now we have the number of digits in total and the integer digits as well
1023 as the exponent and its sign. We can decide whether the read digits are
1024 really integer digits or belong to the fractional part; i.e. we normalize
1027 register int incr
= (exponent
< 0 ? MAX (-int_no
, exponent
)
1028 : MIN (dig_no
- int_no
, exponent
));
1033 if (__builtin_expect (int_no
+ exponent
> MAX_10_EXP
+ 1, 0))
1035 __set_errno (ERANGE
);
1036 return negative
? -FLOAT_HUGE_VAL
: FLOAT_HUGE_VAL
;
1039 if (__builtin_expect (exponent
< MIN_10_EXP
- (DIG
+ 1), 0))
1041 __set_errno (ERANGE
);
1042 return negative
? -0.0 : 0.0;
1047 /* Read the integer part as a multi-precision number to NUM. */
1048 startp
= str_to_mpn (startp
, int_no
, num
, &numsize
, &exponent
1049 #ifndef USE_WIDE_CHAR
1050 , decimal
, decimal_len
, thousands
1056 /* We now multiply the gained number by the given power of ten. */
1057 mp_limb_t
*psrc
= num
;
1058 mp_limb_t
*pdest
= den
;
1060 const struct mp_power
*ttab
= &_fpioconst_pow10
[0];
1064 if ((exponent
& expbit
) != 0)
1066 size_t size
= ttab
->arraysize
- _FPIO_CONST_OFFSET
;
1070 /* FIXME: not the whole multiplication has to be
1071 done. If we have the needed number of bits we
1072 only need the information whether more non-zero
1074 if (numsize
>= ttab
->arraysize
- _FPIO_CONST_OFFSET
)
1075 cy
= __mpn_mul (pdest
, psrc
, numsize
,
1076 &__tens
[ttab
->arrayoff
1077 + _FPIO_CONST_OFFSET
],
1080 cy
= __mpn_mul (pdest
, &__tens
[ttab
->arrayoff
1081 + _FPIO_CONST_OFFSET
],
1082 size
, psrc
, numsize
);
1086 (void) SWAP (psrc
, pdest
);
1091 while (exponent
!= 0);
1094 memcpy (num
, den
, numsize
* sizeof (mp_limb_t
));
1097 /* Determine how many bits of the result we already have. */
1098 count_leading_zeros (bits
, num
[numsize
- 1]);
1099 bits
= numsize
* BITS_PER_MP_LIMB
- bits
;
1101 /* Now we know the exponent of the number in base two.
1102 Check it against the maximum possible exponent. */
1103 if (__builtin_expect (bits
> MAX_EXP
, 0))
1105 __set_errno (ERANGE
);
1106 return negative
? -FLOAT_HUGE_VAL
: FLOAT_HUGE_VAL
;
1109 /* We have already the first BITS bits of the result. Together with
1110 the information whether more non-zero bits follow this is enough
1111 to determine the result. */
1112 if (bits
> MANT_DIG
)
1115 const mp_size_t least_idx
= (bits
- MANT_DIG
) / BITS_PER_MP_LIMB
;
1116 const mp_size_t least_bit
= (bits
- MANT_DIG
) % BITS_PER_MP_LIMB
;
1117 const mp_size_t round_idx
= least_bit
== 0 ? least_idx
- 1
1119 const mp_size_t round_bit
= least_bit
== 0 ? BITS_PER_MP_LIMB
- 1
1123 memcpy (retval
, &num
[least_idx
],
1124 RETURN_LIMB_SIZE
* sizeof (mp_limb_t
));
1127 for (i
= least_idx
; i
< numsize
- 1; ++i
)
1128 retval
[i
- least_idx
] = (num
[i
] >> least_bit
)
1130 << (BITS_PER_MP_LIMB
- least_bit
));
1131 if (i
- least_idx
< RETURN_LIMB_SIZE
)
1132 retval
[RETURN_LIMB_SIZE
- 1] = num
[i
] >> least_bit
;
1135 /* Check whether any limb beside the ones in RETVAL are non-zero. */
1136 for (i
= 0; num
[i
] == 0; ++i
)
1139 return round_and_return (retval
, bits
- 1, negative
,
1140 num
[round_idx
], round_bit
,
1141 int_no
< dig_no
|| i
< round_idx
);
1144 else if (dig_no
== int_no
)
1146 const mp_size_t target_bit
= (MANT_DIG
- 1) % BITS_PER_MP_LIMB
;
1147 const mp_size_t is_bit
= (bits
- 1) % BITS_PER_MP_LIMB
;
1149 if (target_bit
== is_bit
)
1151 memcpy (&retval
[RETURN_LIMB_SIZE
- numsize
], num
,
1152 numsize
* sizeof (mp_limb_t
));
1153 /* FIXME: the following loop can be avoided if we assume a
1154 maximal MANT_DIG value. */
1155 MPN_ZERO (retval
, RETURN_LIMB_SIZE
- numsize
);
1157 else if (target_bit
> is_bit
)
1159 (void) __mpn_lshift (&retval
[RETURN_LIMB_SIZE
- numsize
],
1160 num
, numsize
, target_bit
- is_bit
);
1161 /* FIXME: the following loop can be avoided if we assume a
1162 maximal MANT_DIG value. */
1163 MPN_ZERO (retval
, RETURN_LIMB_SIZE
- numsize
);
1168 assert (numsize
< RETURN_LIMB_SIZE
);
1170 cy
= __mpn_rshift (&retval
[RETURN_LIMB_SIZE
- numsize
],
1171 num
, numsize
, is_bit
- target_bit
);
1172 retval
[RETURN_LIMB_SIZE
- numsize
- 1] = cy
;
1173 /* FIXME: the following loop can be avoided if we assume a
1174 maximal MANT_DIG value. */
1175 MPN_ZERO (retval
, RETURN_LIMB_SIZE
- numsize
- 1);
1178 return round_and_return (retval
, bits
- 1, negative
, 0, 0, 0);
1182 /* Store the bits we already have. */
1183 memcpy (retval
, num
, numsize
* sizeof (mp_limb_t
));
1184 #if RETURN_LIMB_SIZE > 1
1185 if (numsize
< RETURN_LIMB_SIZE
)
1186 # if RETURN_LIMB_SIZE == 2
1187 retval
[numsize
] = 0;
1189 MPN_ZERO (retval
+ numsize
, RETURN_LIMB_SIZE
- numsize
);
1194 /* We have to compute at least some of the fractional digits. */
1196 /* We construct a fraction and the result of the division gives us
1197 the needed digits. The denominator is 1.0 multiplied by the
1198 exponent of the lowest digit; i.e. 0.123 gives 123 / 1000 and
1199 123e-6 gives 123 / 1000000. */
1205 mp_limb_t
*psrc
= den
;
1206 mp_limb_t
*pdest
= num
;
1207 const struct mp_power
*ttab
= &_fpioconst_pow10
[0];
1209 assert (dig_no
> int_no
&& exponent
<= 0);
1212 /* For the fractional part we need not process too many digits. One
1213 decimal digits gives us log_2(10) ~ 3.32 bits. If we now compute
1215 digits we should have enough bits for the result. The remaining
1216 decimal digits give us the information that more bits are following.
1217 This can be used while rounding. (Two added as a safety margin.) */
1218 if (dig_no
- int_no
> (MANT_DIG
- bits
+ 2) / 3 + 2)
1220 dig_no
= int_no
+ (MANT_DIG
- bits
+ 2) / 3 + 2;
1226 neg_exp
= dig_no
- int_no
- exponent
;
1228 /* Construct the denominator. */
1233 if ((neg_exp
& expbit
) != 0)
1240 densize
= ttab
->arraysize
- _FPIO_CONST_OFFSET
;
1241 memcpy (psrc
, &__tens
[ttab
->arrayoff
+ _FPIO_CONST_OFFSET
],
1242 densize
* sizeof (mp_limb_t
));
1246 cy
= __mpn_mul (pdest
, &__tens
[ttab
->arrayoff
1247 + _FPIO_CONST_OFFSET
],
1248 ttab
->arraysize
- _FPIO_CONST_OFFSET
,
1250 densize
+= ttab
->arraysize
- _FPIO_CONST_OFFSET
;
1253 (void) SWAP (psrc
, pdest
);
1259 while (neg_exp
!= 0);
1262 memcpy (den
, num
, densize
* sizeof (mp_limb_t
));
1264 /* Read the fractional digits from the string. */
1265 (void) str_to_mpn (startp
, dig_no
- int_no
, num
, &numsize
, &exponent
1266 #ifndef USE_WIDE_CHAR
1267 , decimal
, decimal_len
, thousands
1271 /* We now have to shift both numbers so that the highest bit in the
1272 denominator is set. In the same process we copy the numerator to
1273 a high place in the array so that the division constructs the wanted
1274 digits. This is done by a "quasi fix point" number representation.
1276 num: ddddddddddd . 0000000000000000000000
1278 den: ddddddddddd n >= m
1282 count_leading_zeros (cnt
, den
[densize
- 1]);
1286 /* Don't call `mpn_shift' with a count of zero since the specification
1287 does not allow this. */
1288 (void) __mpn_lshift (den
, den
, densize
, cnt
);
1289 cy
= __mpn_lshift (num
, num
, numsize
, cnt
);
1291 num
[numsize
++] = cy
;
1294 /* Now we are ready for the division. But it is not necessary to
1295 do a full multi-precision division because we only need a small
1296 number of bits for the result. So we do not use __mpn_divmod
1297 here but instead do the division here by hand and stop whenever
1298 the needed number of bits is reached. The code itself comes
1299 from the GNU MP Library by Torbj\"orn Granlund. */
1307 mp_limb_t d
, n
, quot
;
1312 assert (numsize
== 1 && n
< d
);
1316 udiv_qrnnd (quot
, n
, n
, 0, d
);
1323 cnt = BITS_PER_MP_LIMB; \
1325 count_leading_zeros (cnt, quot); \
1327 if (BITS_PER_MP_LIMB - cnt > MANT_DIG) \
1329 used = MANT_DIG + cnt; \
1330 retval[0] = quot >> (BITS_PER_MP_LIMB - used); \
1331 bits = MANT_DIG + 1; \
1335 /* Note that we only clear the second element. */ \
1336 /* The conditional is determined at compile time. */ \
1337 if (RETURN_LIMB_SIZE > 1) \
1343 else if (bits + BITS_PER_MP_LIMB <= MANT_DIG) \
1344 __mpn_lshift_1 (retval, RETURN_LIMB_SIZE, BITS_PER_MP_LIMB, \
1348 used = MANT_DIG - bits; \
1350 __mpn_lshift_1 (retval, RETURN_LIMB_SIZE, used, quot); \
1352 bits += BITS_PER_MP_LIMB
1356 while (bits
<= MANT_DIG
);
1358 return round_and_return (retval
, exponent
- 1, negative
,
1359 quot
, BITS_PER_MP_LIMB
- 1 - used
,
1360 more_bits
|| n
!= 0);
1364 mp_limb_t d0
, d1
, n0
, n1
;
1371 if (numsize
< densize
)
1375 /* The numerator of the number occupies fewer bits than
1376 the denominator but the one limb is bigger than the
1377 high limb of the numerator. */
1384 exponent
-= BITS_PER_MP_LIMB
;
1387 if (bits
+ BITS_PER_MP_LIMB
<= MANT_DIG
)
1388 __mpn_lshift_1 (retval
, RETURN_LIMB_SIZE
,
1389 BITS_PER_MP_LIMB
, 0);
1392 used
= MANT_DIG
- bits
;
1394 __mpn_lshift_1 (retval
, RETURN_LIMB_SIZE
, used
, 0);
1396 bits
+= BITS_PER_MP_LIMB
;
1408 while (bits
<= MANT_DIG
)
1414 /* QUOT should be either 111..111 or 111..110. We need
1415 special treatment of this rare case as normal division
1416 would give overflow. */
1417 quot
= ~(mp_limb_t
) 0;
1420 if (r
< d1
) /* Carry in the addition? */
1422 add_ssaaaa (n1
, n0
, r
- d0
, 0, 0, d0
);
1425 n1
= d0
- (d0
!= 0);
1430 udiv_qrnnd (quot
, r
, n1
, n0
, d1
);
1431 umul_ppmm (n1
, n0
, d0
, quot
);
1435 if (n1
> r
|| (n1
== r
&& n0
> 0))
1437 /* The estimated QUOT was too large. */
1440 sub_ddmmss (n1
, n0
, n1
, n0
, 0, d0
);
1442 if (r
>= d1
) /* If not carry, test QUOT again. */
1445 sub_ddmmss (n1
, n0
, r
, 0, n1
, n0
);
1451 return round_and_return (retval
, exponent
- 1, negative
,
1452 quot
, BITS_PER_MP_LIMB
- 1 - used
,
1453 more_bits
|| n1
!= 0 || n0
!= 0);
1458 mp_limb_t cy
, dX
, d1
, n0
, n1
;
1462 dX
= den
[densize
- 1];
1463 d1
= den
[densize
- 2];
1465 /* The division does not work if the upper limb of the two-limb
1466 numerator is greater than the denominator. */
1467 if (__mpn_cmp (num
, &den
[densize
- numsize
], numsize
) > 0)
1470 if (numsize
< densize
)
1472 mp_size_t empty
= densize
- numsize
;
1476 exponent
-= empty
* BITS_PER_MP_LIMB
;
1479 if (bits
+ empty
* BITS_PER_MP_LIMB
<= MANT_DIG
)
1481 /* We make a difference here because the compiler
1482 cannot optimize the `else' case that good and
1483 this reflects all currently used FLOAT types
1484 and GMP implementations. */
1485 #if RETURN_LIMB_SIZE <= 2
1486 assert (empty
== 1);
1487 __mpn_lshift_1 (retval
, RETURN_LIMB_SIZE
,
1488 BITS_PER_MP_LIMB
, 0);
1490 for (i
= RETURN_LIMB_SIZE
- 1; i
>= empty
; --i
)
1491 retval
[i
] = retval
[i
- empty
];
1498 used
= MANT_DIG
- bits
;
1499 if (used
>= BITS_PER_MP_LIMB
)
1502 (void) __mpn_lshift (&retval
[used
1503 / BITS_PER_MP_LIMB
],
1506 - used
/ BITS_PER_MP_LIMB
),
1507 used
% BITS_PER_MP_LIMB
);
1508 for (i
= used
/ BITS_PER_MP_LIMB
- 1; i
>= 0; --i
)
1512 __mpn_lshift_1 (retval
, RETURN_LIMB_SIZE
, used
, 0);
1514 bits
+= empty
* BITS_PER_MP_LIMB
;
1516 for (i
= numsize
; i
> 0; --i
)
1517 num
[i
+ empty
] = num
[i
- 1];
1518 MPN_ZERO (num
, empty
+ 1);
1523 assert (numsize
== densize
);
1524 for (i
= numsize
; i
> 0; --i
)
1525 num
[i
] = num
[i
- 1];
1532 while (bits
<= MANT_DIG
)
1535 /* This might over-estimate QUOT, but it's probably not
1536 worth the extra code here to find out. */
1537 quot
= ~(mp_limb_t
) 0;
1542 udiv_qrnnd (quot
, r
, n0
, num
[densize
- 1], dX
);
1543 umul_ppmm (n1
, n0
, d1
, quot
);
1545 while (n1
> r
|| (n1
== r
&& n0
> num
[densize
- 2]))
1549 if (r
< dX
) /* I.e. "carry in previous addition?" */
1556 /* Possible optimization: We already have (q * n0) and (1 * n1)
1557 after the calculation of QUOT. Taking advantage of this, we
1558 could make this loop make two iterations less. */
1560 cy
= __mpn_submul_1 (num
, den
, densize
+ 1, quot
);
1562 if (num
[densize
] != cy
)
1564 cy
= __mpn_add_n (num
, num
, den
, densize
);
1568 n0
= num
[densize
] = num
[densize
- 1];
1569 for (i
= densize
- 1; i
> 0; --i
)
1570 num
[i
] = num
[i
- 1];
1576 for (i
= densize
; num
[i
] == 0 && i
>= 0; --i
)
1578 return round_and_return (retval
, exponent
- 1, negative
,
1579 quot
, BITS_PER_MP_LIMB
- 1 - used
,
1580 more_bits
|| i
>= 0);
1587 #if defined _LIBC && !defined USE_WIDE_CHAR
1588 libc_hidden_def (____STRTOF_INTERNAL
)
1591 /* External user entry point. */
1594 #ifdef weak_function
1597 __STRTOF (nptr
, endptr
, loc
)
1598 const STRING_TYPE
*nptr
;
1599 STRING_TYPE
**endptr
;
1602 return ____STRTOF_INTERNAL (nptr
, endptr
, 0, loc
);
1605 libc_hidden_def (__STRTOF
)
1606 libc_hidden_ver (__STRTOF
, STRTOF
)
1608 weak_alias (__STRTOF
, STRTOF
)
1610 #ifdef LONG_DOUBLE_COMPAT
1611 # if LONG_DOUBLE_COMPAT(libc, GLIBC_2_1)
1612 # ifdef USE_WIDE_CHAR
1613 compat_symbol (libc
, __wcstod_l
, __wcstold_l
, GLIBC_2_1
);
1615 compat_symbol (libc
, __strtod_l
, __strtold_l
, GLIBC_2_1
);
1618 # if LONG_DOUBLE_COMPAT(libc, GLIBC_2_3)
1619 # ifdef USE_WIDE_CHAR
1620 compat_symbol (libc
, wcstod_l
, wcstold_l
, GLIBC_2_3
);
1622 compat_symbol (libc
, strtod_l
, strtold_l
, GLIBC_2_3
);