newlib/libc/stdlib/strtod.c

   1 /*
   2 FUNCTION
   3         <<strtod>>, <<strtof>>, <<strtold>>, <<strtod_l>>, <<strtof_l>>, <<strtold_l>>---string to double or float
   4
   5 INDEX
   6         strtod
   7
   8 INDEX
   9         strtof
  10
  11 INDEX
  12         strtold
  13
  14 INDEX
  15         strtod_l
  16
  17 INDEX
  18         strtof_l
  19
  20 INDEX
  21         strtold_l
  22
  23 INDEX
  24         _strtod_r
  25
  26 SYNOPSIS
  27         #include <stdlib.h>
  28         double strtod(const char *restrict <[str]>, char **restrict <[tail]>);
  29         float strtof(const char *restrict <[str]>, char **restrict <[tail]>);
  30         long double strtold(const char *restrict <[str]>,
  31                             char **restrict <[tail]>);
  32
  33         #include <stdlib.h>
  34         double strtod_l(const char *restrict <[str]>, char **restrict <[tail]>,
  35                         locale_t <[locale]>);
  36         float strtof_l(const char *restrict <[str]>, char **restrict <[tail]>,
  37                        locale_t <[locale]>);
  38         long double strtold_l(const char *restrict <[str]>,
  39                               char **restrict <[tail]>,
  40                               locale_t <[locale]>);
  41
  42         double _strtod_r(void *<[reent]>,
  43                          const char *restrict <[str]>, char **restrict <[tail]>);
  44
  45 DESCRIPTION
  46         <<strtod>>, <<strtof>>, <<strtold>> parse the character string
  47         <[str]>, producing a substring which can be converted to a double,
  48         float, or long double value, respectively.  The substring converted
  49         is the longest initial subsequence of <[str]>, beginning with the
  50         first non-whitespace character, that has one of these formats:
  51         .[+|-]<[digits]>[.[<[digits]>]][(e|E)[+|-]<[digits]>]
  52         .[+|-].<[digits]>[(e|E)[+|-]<[digits]>]
  53         .[+|-](i|I)(n|N)(f|F)[(i|I)(n|N)(i|I)(t|T)(y|Y)]
  54         .[+|-](n|N)(a|A)(n|N)[<(>[<[hexdigits]>]<)>]
  55         .[+|-]0(x|X)<[hexdigits]>[.[<[hexdigits]>]][(p|P)[+|-]<[digits]>]
  56         .[+|-]0(x|X).<[hexdigits]>[(p|P)[+|-]<[digits]>]
  57         The substring contains no characters if <[str]> is empty, consists
  58         entirely of whitespace, or if the first non-whitespace
  59         character is something other than <<+>>, <<->>, <<.>>, or a
  60         digit, and cannot be parsed as infinity or NaN. If the platform
  61         does not support NaN, then NaN is treated as an empty substring.
  62         If the substring is empty, no conversion is done, and
  63         the value of <[str]> is stored in <<*<[tail]>>>.  Otherwise,
  64         the substring is converted, and a pointer to the final string
  65         (which will contain at least the terminating null character of
  66         <[str]>) is stored in <<*<[tail]>>>.  If you want no
  67         assignment to <<*<[tail]>>>, pass a null pointer as <[tail]>.
  68
  69         This implementation returns the nearest machine number to the
  70         input decimal string.  Ties are broken by using the IEEE
  71         round-even rule.  However, <<strtof>> is currently subject to
  72         double rounding errors.
  73
  74         <<strtod_l>>, <<strtof_l>>, <<strtold_l>> are like <<strtod>>,
  75         <<strtof>>, <<strtold>> but perform the conversion based on the
  76         locale specified by the locale object locale.  If <[locale]> is
  77         LC_GLOBAL_LOCALE or not a valid locale object, the behaviour is
  78         undefined.
  79
  80         The alternate function <<_strtod_r>> is a reentrant version.
  81         The extra argument <[reent]> is a pointer to a reentrancy structure.
  82
  83 RETURNS
  84         These functions return the converted substring value, if any.  If
  85         no conversion could be performed, 0 is returned.  If the correct
  86         value is out of the range of representable values, plus or minus
  87         <<HUGE_VAL>> (<<HUGE_VALF>>, <<HUGE_VALL>>) is returned, and
  88         <<ERANGE>> is stored in errno. If the correct value would cause
  89         underflow, 0 is returned and <<ERANGE>> is stored in errno.
  90
  91 PORTABILITY
  92 <<strtod>> is ANSI.
  93 <<strtof>>, <<strtold>> are C99.
  94 <<strtod_l>>, <<strtof_l>>, <<strtold_l>> are GNU extensions.
  95
  96 Supporting OS subroutines required: <<close>>, <<fstat>>, <<isatty>>,
  97 <<lseek>>, <<read>>, <<sbrk>>, <<write>>.
  98 */
  99
 100 /****************************************************************
 101
 102 The author of this software is David M. Gay.
 103
 104 Copyright (C) 1998-2001 by Lucent Technologies
 105 All Rights Reserved
 106
 107 Permission to use, copy, modify, and distribute this software and
 108 its documentation for any purpose and without fee is hereby
 109 granted, provided that the above copyright notice appear in all
 110 copies and that both that the copyright notice and this
 111 permission notice and warranty disclaimer appear in supporting
 112 documentation, and that the name of Lucent or any of its entities
 113 not be used in advertising or publicity pertaining to
 114 distribution of the software without specific, written prior
 115 permission.
 116
 117 LUCENT DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE,
 118 INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS.
 119 IN NO EVENT SHALL LUCENT OR ANY OF ITS ENTITIES BE LIABLE FOR ANY
 120 SPECIAL, INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
 121 WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER
 122 IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION,
 123 ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF
 124 THIS SOFTWARE.
 125
 126 ****************************************************************/
 127
 128 /* Please send bug reports to David M. Gay (dmg at acm dot org,
 129  * with " at " changed at "@" and " dot " changed to ".").      */
 130
 131 /* Original file gdtoa-strtod.c Modified 06-21-2006 by Jeff Johnston to work within newlib.  */
 132
 133 #define _GNU_SOURCE
 134 #include <_ansi.h>
 135 #include <errno.h>
 136 #include <stdlib.h>
 137 #include <string.h>
 138 #include "mprec.h"
 139 #include "gdtoa.h"
 140 #include "../locale/setlocale.h"
 141
 142 /* #ifndef NO_FENV_H */
 143 /* #include <fenv.h> */
 144 /* #endif */
 145
 146 #include "locale.h"
 147
 148 #ifdef IEEE_Arith
 149 #ifndef NO_IEEE_Scale
 150 #define Avoid_Underflow
 151 #undef tinytens
 152 /* The factor of 2^106 in tinytens[4] helps us avoid setting the underflow */
 153 /* flag unnecessarily.  It leads to a song and dance at the end of strtod. */
 154 static const double tinytens[] = { 1e-16, 1e-32,
 155 #ifdef _DOUBLE_IS_32BITS
 156                                     0.0, 0.0, 0.0
 157 #else
 158                                     1e-64, 1e-128,
 159                                     9007199254740992. * 9007199254740992.e-256
 160 #endif
 161                                   };
 162
 163 #endif
 164 #endif
 165
 166 #ifdef Honor_FLT_ROUNDS
 167 #define Rounding rounding
 168 #undef Check_FLT_ROUNDS
 169 #define Check_FLT_ROUNDS
 170 #else
 171 #define Rounding Flt_Rounds
 172 #endif
 173
 174 #ifdef IEEE_MC68k
 175 #define _0 0
 176 #define _1 1
 177 #else
 178 #define _0 1
 179 #define _1 0
 180 #endif
 181
 182 #ifdef Avoid_Underflow /*{*/
 183  static double
 184 sulp (U x,
 185         int scale)
 186 {
 187         U u;
 188         double rv;
 189         int i;
 190
 191         rv = ulp(dval(x));
 192         if (!scale || (i = 2*P + 1 - ((dword0(x) & Exp_mask) >> Exp_shift)) <= 0)
 193                 return rv; /* Is there an example where i <= 0 ? */
 194         dword0(u) = Exp_1 + ((__int32_t)i << Exp_shift);
 195 #ifndef _DOUBLE_IS_32BITS
 196         dword1(u) = 0;
 197 #endif
 198         return rv * u.d;
 199         }
 200 #endif /*}*/
 201
 202
 203 #ifndef NO_HEX_FP
 204
 205 static void
 206 ULtod (__ULong *L,
 207         __ULong *bits,
 208         Long exp,
 209         int k)
 210 {
 211         switch(k & STRTOG_Retmask) {
 212           case STRTOG_NoNumber:
 213           case STRTOG_Zero:
 214                 L[0] = L[1] = 0;
 215                 break;
 216
 217           case STRTOG_Denormal:
 218                 L[_1] = bits[0];
 219                 L[_0] = bits[1];
 220                 break;
 221
 222           case STRTOG_Normal:
 223           case STRTOG_NaNbits:
 224                 L[_1] = bits[0];
 225                 L[_0] = (bits[1] & ~0x100000) | ((exp + 0x3ff + 52) << 20);
 226                 break;
 227
 228           case STRTOG_Infinite:
 229                 L[_0] = 0x7ff00000;
 230                 L[_1] = 0;
 231                 break;
 232
 233           case STRTOG_NaN:
 234                 L[_0] = 0x7fffffff;
 235                 L[_1] = (__ULong)-1;
 236           }
 237         if (k & STRTOG_Neg)
 238                 L[_0] |= 0x80000000L;
 239 }
 240 #endif /* !NO_HEX_FP */
 241
 242 double
 243 _strtod_l (struct _reent *ptr, const char *__restrict s00, char **__restrict se,
 244            locale_t loc)
 245 {
 246 #ifdef Avoid_Underflow
 247         int scale;
 248 #endif
 249         int bb2, bb5, bbe, bd2, bd5, bbbits, bs2, c, decpt, dsign,
 250                  e, e1, esign, i, j, k, nd, nd0, nf, nz, nz0, sign;
 251         const char *s, *s0, *s1;
 252         double aadj, adj;
 253         U aadj1, rv, rv0;
 254         Long L;
 255         __ULong y, z;
 256         _Bigint *bb = NULL, *bb1, *bd = NULL, *bd0, *bs = NULL, *delta = NULL;
 257 #ifdef Avoid_Underflow
 258         __ULong Lsb, Lsb1;
 259 #endif
 260 #ifdef SET_INEXACT
 261         int inexact, oldinexact;
 262 #endif
 263 #ifdef Honor_FLT_ROUNDS
 264         int rounding;
 265 #endif
 266 #ifdef __HAVE_LOCALE_INFO__
 267         const char *decimal_point = __get_numeric_locale(loc)->decimal_point;
 268         const int dec_len = strlen (decimal_point);
 269 #else
 270         const char *decimal_point = ".";
 271         const int dec_len = 1;
 272 #endif
 273
 274         delta = bs = bd = NULL;
 275         sign = nz0 = nz = decpt = 0;
 276         dval(rv) = 0.;
 277         for(s = s00;;s++) switch(*s) {
 278                 case '-':
 279                         sign = 1;
 280                         /* no break */
 281                 case '+':
 282                         if (*++s)
 283                                 goto break2;
 284                         /* no break */
 285                 case 0:
 286                         goto ret0;
 287                 case '\t':
 288                 case '\n':
 289                 case '\v':
 290                 case '\f':
 291                 case '\r':
 292                 case ' ':
 293                         continue;
 294                 default:
 295                         goto break2;
 296                 }
 297  break2:
 298         if (*s == '0') {
 299 #ifndef NO_HEX_FP
 300                 {
 301                 static const FPI fpi = { 53, 1-1023-53+1, 2046-1023-53+1, 1, SI };
 302                 Long exp;
 303                 __ULong bits[2];
 304                 switch(s[1]) {
 305                   case 'x':
 306                   case 'X':
 307                         /* If the number is not hex, then the parse of
 308                            0 is still valid.  */
 309                         s00 = s + 1;
 310                         {
 311 #if defined(FE_DOWNWARD) && defined(FE_TONEAREST) && defined(FE_TOWARDZERO) && defined(FE_UPWARD)
 312                         FPI fpi1 = fpi;
 313                         switch(fegetround()) {
 314                           case FE_TOWARDZERO:   fpi1.rounding = 0; break;
 315                           case FE_UPWARD:       fpi1.rounding = 2; break;
 316                           case FE_DOWNWARD:     fpi1.rounding = 3;
 317                           }
 318 #else
 319 #define fpi1 fpi
 320 #endif
 321                         switch((i = gethex(ptr, &s, &fpi1, &exp, &bb, sign, loc)) & STRTOG_Retmask) {
 322                           case STRTOG_NoNumber:
 323                                 s = s00;
 324                                 sign = 0;
 325                                 /* FALLTHROUGH */
 326                           case STRTOG_Zero:
 327                                 break;
 328                           default:
 329                                 if (bb) {
 330                                         copybits(bits, fpi.nbits, bb);
 331                                         Bfree(ptr,bb);
 332                                         }
 333                                 ULtod(rv.i, bits, exp, i);
 334 #ifndef NO_ERRNO
 335                                 /* try to avoid the bug of testing an 8087 register value */
 336                                 if ((dword0(rv)&Exp_mask) == 0)
 337                                     errno = ERANGE;
 338 #endif
 339                           }}
 340                         goto ret;
 341                   }
 342                 }
 343 #endif
 344                 nz0 = 1;
 345                 while(*++s == '0') ;
 346                 if (!*s)
 347                         goto ret;
 348                 }
 349         s0 = s;
 350         y = z = 0;
 351         for(nd = nf = 0; (c = *s) >= '0' && c <= '9'; nd++, s++)
 352                 if (nd < 9)
 353                         y = 10*y + c - '0';
 354                 else
 355                         z = 10*z + c - '0';
 356         nd0 = nd;
 357         if (strncmp (s, decimal_point, dec_len) == 0)
 358                 {
 359                 decpt = 1;
 360                 c = *(s += dec_len);
 361                 if (!nd) {
 362                         for(; c == '0'; c = *++s)
 363                                 nz++;
 364                         if (c > '0' && c <= '9') {
 365                                 s0 = s;
 366                                 nf += nz;
 367                                 nz = 0;
 368                                 goto have_dig;
 369                                 }
 370                         goto dig_done;
 371                         }
 372                 for(; c >= '0' && c <= '9'; c = *++s) {
 373  have_dig:
 374                         nz++;
 375                         if (c -= '0') {
 376                                 nf += nz;
 377                                 for(i = 1; i < nz; i++)
 378                                         if (nd++ < 9)
 379                                                 y *= 10;
 380                                         else if (nd <= DBL_DIG + 1)
 381                                                 z *= 10;
 382                                 if (nd++ < 9)
 383                                         y = 10*y + c;
 384                                 else if (nd <= DBL_DIG + 1)
 385                                         z = 10*z + c;
 386                                 nz = 0;
 387                                 }
 388                         }
 389                 }
 390  dig_done:
 391         e = 0;
 392         if (c == 'e' || c == 'E') {
 393                 if (!nd && !nz && !nz0) {
 394                         goto ret0;
 395                         }
 396                 s00 = s;
 397                 esign = 0;
 398                 switch(c = *++s) {
 399                         case '-':
 400                                 esign = 1;
 401                         case '+':
 402                                 c = *++s;
 403                         }
 404                 if (c >= '0' && c <= '9') {
 405                         while(c == '0')
 406                                 c = *++s;
 407                         if (c > '0' && c <= '9') {
 408                                 L = c - '0';
 409                                 s1 = s;
 410                                 while((c = *++s) >= '0' && c <= '9')
 411                                         L = 10*L + c - '0';
 412                                 if (s - s1 > 8 || L > 19999)
 413                                         /* Avoid confusion from exponents
 414                                          * so large that e might overflow.
 415                                          */
 416                                         e = 19999; /* safe for 16 bit ints */
 417                                 else
 418                                         e = (int)L;
 419                                 if (esign)
 420                                         e = -e;
 421                                 }
 422                         else
 423                                 e = 0;
 424                         }
 425                 else
 426                         s = s00;
 427                 }
 428         if (!nd) {
 429                 if (!nz && !nz0) {
 430 #ifdef INFNAN_CHECK
 431                         /* Check for Nan and Infinity */
 432                         __ULong bits[2];
 433                         static const FPI fpinan =       /* only 52 explicit bits */
 434                                 { 52, 1-1023-53+1, 2046-1023-53+1, 1, SI };
 435                         if (!decpt)
 436                          switch(c) {
 437                           case 'i':
 438                           case 'I':
 439                                 if (match(&s,"nf")) {
 440                                         --s;
 441                                         if (!match(&s,"inity"))
 442                                                 ++s;
 443                                         dword0(rv) = 0x7ff00000;
 444 #ifndef _DOUBLE_IS_32BITS
 445                                         dword1(rv) = 0;
 446 #endif /*!_DOUBLE_IS_32BITS*/
 447                                         goto ret;
 448                                         }
 449                                 break;
 450                           case 'n':
 451                           case 'N':
 452                                 if (match(&s, "an")) {
 453 #ifndef No_Hex_NaN
 454                                         if (*s == '(' /*)*/
 455                                          && hexnan(&s, &fpinan, bits)
 456                                                         == STRTOG_NaNbits) {
 457                                                 dword0(rv) = 0x7ff00000 | bits[1];
 458 #ifndef _DOUBLE_IS_32BITS
 459                                                 dword1(rv) = bits[0];
 460 #endif /*!_DOUBLE_IS_32BITS*/
 461                                                 }
 462                                         else {
 463 #endif
 464                                                 dval(rv) = nan ("");
 465 #ifndef No_Hex_NaN
 466                                                 }
 467 #endif
 468                                         goto ret;
 469                                         }
 470                           }
 471 #endif /* INFNAN_CHECK */
 472  ret0:
 473                         s = s00;
 474                         sign = 0;
 475                         }
 476                 goto ret;
 477                 }
 478         e1 = e -= nf;
 479
 480         /* Now we have nd0 digits, starting at s0, followed by a
 481          * decimal point, followed by nd-nd0 digits.  The number we're
 482          * after is the integer represented by those digits times
 483          * 10**e */
 484
 485         if (!nd0)
 486                 nd0 = nd;
 487         k = nd < DBL_DIG + 1 ? nd : DBL_DIG + 1;
 488         dval(rv) = y;
 489         if (k > 9) {
 490 #ifdef SET_INEXACT
 491                 if (k > DBL_DIG)
 492                         oldinexact = get_inexact();
 493 #endif
 494                 dval(rv) = tens[k - 9] * dval(rv) + z;
 495                 }
 496         bd0 = 0;
 497         if (nd <= DBL_DIG
 498 #ifndef RND_PRODQUOT
 499 #ifndef Honor_FLT_ROUNDS
 500                 && Flt_Rounds == 1
 501 #endif
 502 #endif
 503                         ) {
 504                 if (!e)
 505                         goto ret;
 506                 if (e > 0) {
 507                         if (e <= Ten_pmax) {
 508 #ifdef VAX
 509                                 goto vax_ovfl_check;
 510 #else
 511 #ifdef Honor_FLT_ROUNDS
 512                                 /* round correctly FLT_ROUNDS = 2 or 3 */
 513                                 if (sign) {
 514                                         dval(rv) = -dval(rv);
 515                                         sign = 0;
 516                                         }
 517 #endif
 518                                 /* rv = */ rounded_product(dval(rv), tens[e]);
 519                                 goto ret;
 520 #endif
 521                                 }
 522                         i = DBL_DIG - nd;
 523                         if (e <= Ten_pmax + i) {
 524                                 /* A fancier test would sometimes let us do
 525                                  * this for larger i values.
 526                                  */
 527 #ifdef Honor_FLT_ROUNDS
 528                                 /* round correctly FLT_ROUNDS = 2 or 3 */
 529                                 if (sign) {
 530                                         dval(rv) = -dval(rv);
 531                                         sign = 0;
 532                                         }
 533 #endif
 534                                 e -= i;
 535                                 dval(rv) *= tens[i];
 536 #ifdef VAX
 537                                 /* VAX exponent range is so narrow we must
 538                                  * worry about overflow here...
 539                                  */
 540  vax_ovfl_check:
 541                                 dword0(rv) -= P*Exp_msk1;
 542                                 /* rv = */ rounded_product(dval(rv), tens[e]);
 543                                 if ((dword0(rv) & Exp_mask)
 544                                  > Exp_msk1*(DBL_MAX_EXP+Bias-1-P))
 545                                         goto ovfl;
 546                                 dword0(rv) += P*Exp_msk1;
 547 #else
 548                                 /* rv = */ rounded_product(dval(rv), tens[e]);
 549 #endif
 550                                 goto ret;
 551                                 }
 552                         }
 553 #ifndef Inaccurate_Divide
 554                 else if (e >= -Ten_pmax) {
 555 #ifdef Honor_FLT_ROUNDS
 556                         /* round correctly FLT_ROUNDS = 2 or 3 */
 557                         if (sign) {
 558                                 dval(rv) = -dval(rv);
 559                                 sign = 0;
 560                                 }
 561 #endif
 562                         /* rv = */ rounded_quotient(dval(rv), tens[-e]);
 563                         goto ret;
 564                         }
 565 #endif
 566                 }
 567         e1 += nd - k;
 568
 569 #ifdef IEEE_Arith
 570 #ifdef SET_INEXACT
 571         inexact = 1;
 572         if (k <= DBL_DIG)
 573                 oldinexact = get_inexact();
 574 #endif
 575 #ifdef Avoid_Underflow
 576         scale = 0;
 577 #endif
 578 #ifdef Honor_FLT_ROUNDS
 579         if ((rounding = Flt_Rounds) >= 2) {
 580                 if (sign)
 581                         rounding = rounding == 2 ? 0 : 2;
 582                 else
 583                         if (rounding != 2)
 584                                 rounding = 0;
 585                 }
 586 #endif
 587 #endif /*IEEE_Arith*/
 588
 589         /* Get starting approximation = rv * 10**e1 */
 590
 591         if (e1 > 0) {
 592                 if ( (i = e1 & 15) !=0)
 593                         dval(rv) *= tens[i];
 594                 if (e1 &= ~15) {
 595                         if (e1 > DBL_MAX_10_EXP) {
 596  ovfl:
 597 #ifndef NO_ERRNO
 598                                 _REENT_ERRNO(ptr) = ERANGE;
 599 #endif
 600                                 /* Can't trust HUGE_VAL */
 601 #ifdef IEEE_Arith
 602 #ifdef Honor_FLT_ROUNDS
 603                                 switch(rounding) {
 604                                   case 0: /* toward 0 */
 605                                   case 3: /* toward -infinity */
 606                                         dword0(rv) = Big0;
 607 #ifndef _DOUBLE_IS_32BITS
 608                                         dword1(rv) = Big1;
 609 #endif /*!_DOUBLE_IS_32BITS*/
 610                                         break;
 611                                   default:
 612                                         dword0(rv) = Exp_mask;
 613 #ifndef _DOUBLE_IS_32BITS
 614                                         dword1(rv) = 0;
 615 #endif /*!_DOUBLE_IS_32BITS*/
 616                                   }
 617 #else /*Honor_FLT_ROUNDS*/
 618                                 dword0(rv) = Exp_mask;
 619 #ifndef _DOUBLE_IS_32BITS
 620                                 dword1(rv) = 0;
 621 #endif /*!_DOUBLE_IS_32BITS*/
 622 #endif /*Honor_FLT_ROUNDS*/
 623 #ifdef SET_INEXACT
 624                                 /* set overflow bit */
 625                                 dval(rv0) = 1e300;
 626                                 dval(rv0) *= dval(rv0);
 627 #endif
 628 #else /*IEEE_Arith*/
 629                                 dword0(rv) = Big0;
 630 #ifndef _DOUBLE_IS_32BITS
 631                                 dword1(rv) = Big1;
 632 #endif /*!_DOUBLE_IS_32BITS*/
 633 #endif /*IEEE_Arith*/
 634                                 if (bd0)
 635                                         goto retfree;
 636                                 goto ret;
 637                                 }
 638                         e1 >>= 4;
 639                         for(j = 0; e1 > 1; j++, e1 >>= 1)
 640                                 if (e1 & 1)
 641                                         dval(rv) *= bigtens[j];
 642                 /* The last multiplication could overflow. */
 643                         dword0(rv) -= P*Exp_msk1;
 644                         dval(rv) *= bigtens[j];
 645                         if ((z = dword0(rv) & Exp_mask)
 646                          > Exp_msk1*(DBL_MAX_EXP+Bias-P))
 647                                 goto ovfl;
 648                         if (z > Exp_msk1*(DBL_MAX_EXP+Bias-1-P)) {
 649                                 /* set to largest number */
 650                                 /* (Can't trust DBL_MAX) */
 651                                 dword0(rv) = Big0;
 652 #ifndef _DOUBLE_IS_32BITS
 653                                 dword1(rv) = Big1;
 654 #endif /*!_DOUBLE_IS_32BITS*/
 655                                 }
 656                         else
 657                                 dword0(rv) += P*Exp_msk1;
 658                         }
 659                 }
 660         else if (e1 < 0) {
 661                 e1 = -e1;
 662                 if ( (i = e1 & 15) !=0)
 663                         dval(rv) /= tens[i];
 664                 if (e1 >>= 4) {
 665                         if (e1 >= 1 << n_bigtens)
 666                                 goto undfl;
 667 #ifdef Avoid_Underflow
 668                         if (e1 & Scale_Bit)
 669                                 scale = 2*P;
 670                         for(j = 0; e1 > 0; j++, e1 >>= 1)
 671                                 if (e1 & 1)
 672                                         dval(rv) *= tinytens[j];
 673                         if (scale && (j = 2*P + 1 - ((dword0(rv) & Exp_mask)
 674                                                 >> Exp_shift)) > 0) {
 675                                 /* scaled rv is denormal; zap j low bits */
 676                                 if (j >= 32) {
 677 #ifndef _DOUBLE_IS_32BITS
 678                                         dword1(rv) = 0;
 679 #endif /*!_DOUBLE_IS_32BITS*/
 680                                         if (j >= 53)
 681                                          dword0(rv) = (P+2)*Exp_msk1;
 682                                         else
 683                                          dword0(rv) &= 0xffffffff << (j-32);
 684                                         }
 685 #ifndef _DOUBLE_IS_32BITS
 686                                 else
 687                                         dword1(rv) &= 0xffffffff << j;
 688 #endif /*!_DOUBLE_IS_32BITS*/
 689                                 }
 690 #else
 691                         for(j = 0; e1 > 1; j++, e1 >>= 1)
 692                                 if (e1 & 1)
 693                                         dval(rv) *= tinytens[j];
 694                         /* The last multiplication could underflow. */
 695                         dval(rv0) = dval(rv);
 696                         dval(rv) *= tinytens[j];
 697                         if (!dval(rv)) {
 698                                 dval(rv) = 2.*dval(rv0);
 699                                 dval(rv) *= tinytens[j];
 700 #endif
 701                                 if (!dval(rv)) {
 702  undfl:
 703                                         dval(rv) = 0.;
 704 #ifndef NO_ERRNO
 705                                         _REENT_ERRNO(ptr) = ERANGE;
 706 #endif
 707                                         if (bd0)
 708                                                 goto retfree;
 709                                         goto ret;
 710                                         }
 711 #ifndef Avoid_Underflow
 712 #ifndef _DOUBLE_IS_32BITS
 713                                 dword0(rv) = Tiny0;
 714                                 dword1(rv) = Tiny1;
 715 #else
 716                                 dword0(rv) = Tiny1;
 717 #endif /*_DOUBLE_IS_32BITS*/
 718                                 /* The refinement below will clean
 719                                  * this approximation up.
 720                                  */
 721                                 }
 722 #endif
 723                         }
 724                 }
 725
 726         /* Now the hard part -- adjusting rv to the correct value.*/
 727
 728         /* Put digits into bd: true value = bd * 10^e */
 729
 730         bd0 = s2b(ptr, s0, nd0, nd, y);
 731         if (bd0 == NULL)
 732                 goto ovfl;
 733
 734         for(;;) {
 735                 bd = Balloc(ptr,bd0->_k);
 736                 if (bd == NULL)
 737                         goto ovfl;
 738                 Bcopy(bd, bd0);
 739                 bb = d2b(ptr,dval(rv), &bbe, &bbbits);  /* rv = bb * 2^bbe */
 740                 if (bb == NULL)
 741                         goto ovfl;
 742                 bs = i2b(ptr,1);
 743                 if (bs == NULL)
 744                         goto ovfl;
 745
 746                 if (e >= 0) {
 747                         bb2 = bb5 = 0;
 748                         bd2 = bd5 = e;
 749                         }
 750                 else {
 751                         bb2 = bb5 = -e;
 752                         bd2 = bd5 = 0;
 753                         }
 754                 if (bbe >= 0)
 755                         bb2 += bbe;
 756                 else
 757                         bd2 -= bbe;
 758                 bs2 = bb2;
 759 #ifdef Honor_FLT_ROUNDS
 760                 if (rounding != 1)
 761                         bs2++;
 762 #endif
 763 #ifdef Avoid_Underflow
 764                 Lsb = LSB;
 765                 Lsb1 = 0;
 766                 j = bbe - scale;
 767                 i = j + bbbits - 1;     /* logb(rv) */
 768                 j = P + 1 - bbbits;
 769                 if (i < Emin) { /* denormal */
 770                         i = Emin - i;
 771                         j -= i;
 772                         if (i < 32)
 773                                 Lsb <<= i;
 774                         else
 775                                 Lsb1 = Lsb << (i-32);
 776                         }
 777 #else /*Avoid_Underflow*/
 778 #ifdef Sudden_Underflow
 779 #ifdef IBM
 780                 j = 1 + 4*P - 3 - bbbits + ((bbe + bbbits - 1) & 3);
 781 #else
 782                 j = P + 1 - bbbits;
 783 #endif
 784 #else /*Sudden_Underflow*/
 785                 j = bbe;
 786                 i = j + bbbits - 1;     /* logb(rv) */
 787                 if (i < Emin)   /* denormal */
 788                         j += P - Emin;
 789                 else
 790                         j = P + 1 - bbbits;
 791 #endif /*Sudden_Underflow*/
 792 #endif /*Avoid_Underflow*/
 793                 bb2 += j;
 794                 bd2 += j;
 795 #ifdef Avoid_Underflow
 796                 bd2 += scale;
 797 #endif
 798                 i = bb2 < bd2 ? bb2 : bd2;
 799                 if (i > bs2)
 800                         i = bs2;
 801                 if (i > 0) {
 802                         bb2 -= i;
 803                         bd2 -= i;
 804                         bs2 -= i;
 805                         }
 806                 if (bb5 > 0) {
 807                         bs = pow5mult(ptr, bs, bb5);
 808                         if (bs == NULL)
 809                                 goto ovfl;
 810                         bb1 = mult(ptr, bs, bb);
 811                         if (bb1 == NULL)
 812                                 goto ovfl;
 813                         Bfree(ptr, bb);
 814                         bb = bb1;
 815                         }
 816                 if (bb2 > 0) {
 817                         bb = lshift(ptr, bb, bb2);
 818                         if (bb == NULL)
 819                                 goto ovfl;
 820                         }
 821                 if (bd5 > 0) {
 822                         bd = pow5mult(ptr, bd, bd5);
 823                         if (bd == NULL)
 824                                 goto ovfl;
 825                         }
 826                 if (bd2 > 0) {
 827                         bd = lshift(ptr, bd, bd2);
 828                         if (bd == NULL)
 829                                 goto ovfl;
 830                         }
 831                 if (bs2 > 0) {
 832                         bs = lshift(ptr, bs, bs2);
 833                         if (bs == NULL)
 834                                 goto ovfl;
 835                         }
 836                 delta = diff(ptr, bb, bd);
 837                 if (delta == NULL)
 838                         goto ovfl;
 839                 dsign = delta->_sign;
 840                 delta->_sign = 0;
 841                 i = cmp(delta, bs);
 842 #ifdef Honor_FLT_ROUNDS
 843                 if (rounding != 1) {
 844                         if (i < 0) {
 845                                 /* Error is less than an ulp */
 846                                 if (!delta->_x[0] && delta->_wds <= 1) {
 847                                         /* exact */
 848 #ifdef SET_INEXACT
 849                                         inexact = 0;
 850 #endif
 851                                         break;
 852                                         }
 853                                 if (rounding) {
 854                                         if (dsign) {
 855                                                 adj = 1.;
 856                                                 goto apply_adj;
 857                                                 }
 858                                         }
 859                                 else if (!dsign) {
 860                                         adj = -1.;
 861                                         if (!dword1(rv)
 862                                             && !(dword0(rv) & Frac_mask)) {
 863                                                 y = dword0(rv) & Exp_mask;
 864 #ifdef Avoid_Underflow
 865                                                 if (!scale || y > 2*P*Exp_msk1)
 866 #else
 867                                                 if (y)
 868 #endif
 869                                                   {
 870                                                   delta = lshift(ptr, delta,Log2P);
 871                                                   if (cmp(delta, bs) <= 0)
 872                                                         adj = -0.5;
 873                                                   }
 874                                                 }
 875  apply_adj:
 876 #ifdef Avoid_Underflow
 877                                         if (scale && (y = dword0(rv) & Exp_mask)
 878                                                 <= 2*P*Exp_msk1)
 879                                           dword0(adj) += (2*P+1)*Exp_msk1 - y;
 880 #else
 881 #ifdef Sudden_Underflow
 882                                         if ((dword0(rv) & Exp_mask) <=
 883                                                         P*Exp_msk1) {
 884                                                 dword0(rv) += P*Exp_msk1;
 885                                                 dval(rv) += adj*ulp(dval(rv));
 886                                                 dword0(rv) -= P*Exp_msk1;
 887                                                 }
 888                                         else
 889 #endif /*Sudden_Underflow*/
 890 #endif /*Avoid_Underflow*/
 891                                         dval(rv) += adj*ulp(dval(rv));
 892                                         }
 893                                 break;
 894                                 }
 895                         adj = ratio(delta, bs);
 896                         if (adj < 1.)
 897                                 adj = 1.;
 898                         if (adj <= 0x7ffffffe) {
 899                                 /* adj = rounding ? ceil(adj) : floor(adj); */
 900                                 y = adj;
 901                                 if (y != adj) {
 902                                         if (!((rounding>>1) ^ dsign))
 903                                                 y++;
 904                                         adj = y;
 905                                         }
 906                                 }
 907 #ifdef Avoid_Underflow
 908                         if (scale && (y = dword0(rv) & Exp_mask) <= 2*P*Exp_msk1)
 909                                 dword0(adj) += (2*P+1)*Exp_msk1 - y;
 910 #else
 911 #ifdef Sudden_Underflow
 912                         if ((dword0(rv) & Exp_mask) <= P*Exp_msk1) {
 913                                 dword0(rv) += P*Exp_msk1;
 914                                 adj *= ulp(dval(rv));
 915                                 if (dsign)
 916                                         dval(rv) += adj;
 917                                 else
 918                                         dval(rv) -= adj;
 919                                 dword0(rv) -= P*Exp_msk1;
 920                                 goto cont;
 921                                 }
 922 #endif /*Sudden_Underflow*/
 923 #endif /*Avoid_Underflow*/
 924                         adj *= ulp(dval(rv));
 925                         if (dsign) {
 926                                 if (dword0(rv) == Big0 && dword1(rv) == Big1)
 927                                         goto ovfl;
 928                                 dval(rv) += adj;
 929                         else
 930                                 dval(rv) -= adj;
 931                         goto cont;
 932                         }
 933 #endif /*Honor_FLT_ROUNDS*/
 934
 935                 if (i < 0) {
 936                         /* Error is less than half an ulp -- check for
 937                          * special case of mantissa a power of two.
 938                          */
 939                         if (dsign || dword1(rv) || dword0(rv) & Bndry_mask
 940 #ifdef IEEE_Arith
 941 #ifdef Avoid_Underflow
 942                          || (dword0(rv) & Exp_mask) <= (2*P+1)*Exp_msk1
 943 #else
 944                          || (dword0(rv) & Exp_mask) <= Exp_msk1
 945 #endif
 946 #endif
 947                                 ) {
 948 #ifdef SET_INEXACT
 949                                 if (!delta->x[0] && delta->wds <= 1)
 950                                         inexact = 0;
 951 #endif
 952                                 break;
 953                                 }
 954                         if (!delta->_x[0] && delta->_wds <= 1) {
 955                                 /* exact result */
 956 #ifdef SET_INEXACT
 957                                 inexact = 0;
 958 #endif
 959                                 break;
 960                                 }
 961                         delta = lshift(ptr,delta,Log2P);
 962                         if (cmp(delta, bs) > 0)
 963                                 goto drop_down;
 964                         break;
 965                         }
 966                 if (i == 0) {
 967                         /* exactly half-way between */
 968                         if (dsign) {
 969                                 if ((dword0(rv) & Bndry_mask1) == Bndry_mask1
 970                                  &&  dword1(rv) == (
 971 #ifdef Avoid_Underflow
 972                         (scale && (y = dword0(rv) & Exp_mask) <= 2*P*Exp_msk1)
 973                 ? (0xffffffff & (0xffffffff << (2*P+1-(y>>Exp_shift)))) :
 974 #endif
 975                                                    0xffffffff)) {
 976                                         /*boundary case -- increment exponent*/
 977                                         if (dword0(rv) == Big0 && dword1(rv) == Big1)
 978                                                 goto ovfl;
 979                                         dword0(rv) = (dword0(rv) & Exp_mask)
 980                                                 + Exp_msk1
 981 #ifdef IBM
 982                                                 | Exp_msk1 >> 4
 983 #endif
 984                                                 ;
 985 #ifndef _DOUBLE_IS_32BITS
 986                                         dword1(rv) = 0;
 987 #endif /*!_DOUBLE_IS_32BITS*/
 988 #ifdef Avoid_Underflow
 989                                         dsign = 0;
 990 #endif
 991                                         break;
 992                                         }
 993                                 }
 994                         else if (!(dword0(rv) & Bndry_mask) && !dword1(rv)) {
 995  drop_down:
 996                                 /* boundary case -- decrement exponent */
 997 #ifdef Sudden_Underflow /*{{*/
 998                                 L = dword0(rv) & Exp_mask;
 999 #ifdef IBM
1000                                 if (L <  Exp_msk1)
1001 #else
1002 #ifdef Avoid_Underflow
1003                                 if (L <= (scale ? (2*P+1)*Exp_msk1 : Exp_msk1))
1004 #else
1005                                 if (L <= Exp_msk1)
1006 #endif /*Avoid_Underflow*/
1007 #endif /*IBM*/
1008                                         goto undfl;
1009                                 L -= Exp_msk1;
1010 #else /*Sudden_Underflow}{*/
1011 #ifdef Avoid_Underflow
1012                                 if (scale) {
1013                                         L = dword0(rv) & Exp_mask;
1014                                         if (L <= (2*P+1)*Exp_msk1) {
1015                                                 if (L > (P+2)*Exp_msk1)
1016                                                         /* round even ==> */
1017                                                         /* accept rv */
1018                                                         break;
1019                                                 /* rv = smallest denormal */
1020                                                 goto undfl;
1021                                                 }
1022                                         }
1023 #endif /*Avoid_Underflow*/
1024                                 L = (dword0(rv) & Exp_mask) - Exp_msk1;
1025 #endif /*Sudden_Underflow}*/
1026                                 dword0(rv) = L | Bndry_mask1;
1027 #ifndef _DOUBLE_IS_32BITS
1028                                 dword1(rv) = 0xffffffff;
1029 #endif /*!_DOUBLE_IS_32BITS*/
1030 #ifdef IBM
1031                                 goto cont;
1032 #else
1033                                 break;
1034 #endif
1035                                 }
1036 #ifndef ROUND_BIASED
1037 #ifdef Avoid_Underflow
1038                         if (Lsb1) {
1039                                 if (!(dword0(rv) & Lsb1))
1040                                         break;
1041                                 }
1042                         else if (!(dword1(rv) & Lsb))
1043                                 break;
1044 #else
1045                         if (!(dword1(rv) & LSB))
1046                                 break;
1047 #endif
1048 #endif
1049                         if (dsign)
1050 #ifdef Avoid_Underflow
1051                                 dval(rv) += sulp(rv, scale);
1052 #else
1053                                 dval(rv) += ulp(dval(rv));
1054 #endif
1055 #ifndef ROUND_BIASED
1056                         else {
1057 #ifdef Avoid_Underflow
1058                                 dval(rv) -= sulp(rv, scale);
1059 #else
1060                                 dval(rv) -= ulp(dval(rv));
1061 #endif
1062 #ifndef Sudden_Underflow
1063                                 if (!dval(rv))
1064                                         goto undfl;
1065 #endif
1066                                 }
1067 #ifdef Avoid_Underflow
1068                         dsign = 1 - dsign;
1069 #endif
1070 #endif
1071                         break;
1072                         }
1073                 if ((aadj = ratio(delta, bs)) <= 2.) {
1074                         if (dsign)
1075                                 aadj = dval(aadj1) = 1.;
1076                         else if (dword1(rv) || dword0(rv) & Bndry_mask) {
1077 #ifndef Sudden_Underflow
1078                                 if (dword1(rv) == Tiny1 && !dword0(rv))
1079                                         goto undfl;
1080 #endif
1081                                 aadj = 1.;
1082                                 dval(aadj1) = -1.;
1083                                 }
1084                         else {
1085                                 /* special case -- power of FLT_RADIX to be */
1086                                 /* rounded down... */
1087
1088                                 if (aadj < 2./FLT_RADIX)
1089                                         aadj = 1./FLT_RADIX;
1090                                 else
1091                                         aadj *= 0.5;
1092                                 dval(aadj1) = -aadj;
1093                                 }
1094                         }
1095                 else {
1096                         aadj *= 0.5;
1097                         dval(aadj1) = dsign ? aadj : -aadj;
1098 #ifdef Check_FLT_ROUNDS
1099                         switch(Rounding) {
1100                                 case 2: /* towards +infinity */
1101                                         dval(aadj1) -= 0.5;
1102                                         break;
1103                                 case 0: /* towards 0 */
1104                                 case 3: /* towards -infinity */
1105                                         dval(aadj1) += 0.5;
1106                                 }
1107 #else
1108                         if (Flt_Rounds == 0)
1109                                 dval(aadj1) += 0.5;
1110 #endif /*Check_FLT_ROUNDS*/
1111                         }
1112                 y = dword0(rv) & Exp_mask;
1113
1114                 /* Check for overflow */
1115
1116                 if (y == Exp_msk1*(DBL_MAX_EXP+Bias-1)) {
1117                         dval(rv0) = dval(rv);
1118                         dword0(rv) -= P*Exp_msk1;
1119                         adj = dval(aadj1) * ulp(dval(rv));
1120                         dval(rv) += adj;
1121                         if ((dword0(rv) & Exp_mask) >=
1122                                         Exp_msk1*(DBL_MAX_EXP+Bias-P)) {
1123                                 if (dword0(rv0) == Big0 && dword1(rv0) == Big1)
1124                                         goto ovfl;
1125                                 dword0(rv) = Big0;
1126 #ifndef _DOUBLE_IS_32BITS
1127                                 dword1(rv) = Big1;
1128 #endif /*!_DOUBLE_IS_32BITS*/
1129                                 goto cont;
1130                                 }
1131                         else
1132                                 dword0(rv) += P*Exp_msk1;
1133                         }
1134                 else {
1135 #ifdef Avoid_Underflow
1136                         if (scale && y <= 2*P*Exp_msk1) {
1137                                 if (aadj <= 0x7fffffff) {
1138                                         if ((z = aadj) == 0)
1139                                                 z = 1;
1140                                         aadj = z;
1141                                         dval(aadj1) = dsign ? aadj : -aadj;
1142                                         }
1143                                 dword0(aadj1) += (2*P+1)*Exp_msk1 - y;
1144                                 }
1145                         adj = dval(aadj1) * ulp(dval(rv));
1146                         dval(rv) += adj;
1147 #else
1148 #ifdef Sudden_Underflow
1149                         if ((dword0(rv) & Exp_mask) <= P*Exp_msk1) {
1150                                 dval(rv0) = dval(rv);
1151                                 dword0(rv) += P*Exp_msk1;
1152                                 adj = dval(aadj1) * ulp(dval(rv));
1153                                 dval(rv) += adj;
1154 #ifdef IBM
1155                                 if ((dword0(rv) & Exp_mask) <  P*Exp_msk1)
1156 #else
1157                                 if ((dword0(rv) & Exp_mask) <= P*Exp_msk1)
1158 #endif
1159                                         {
1160                                         if (dword0(rv0) == Tiny0
1161                                          && dword1(rv0) == Tiny1)
1162                                                 goto undfl;
1163 #ifndef _DOUBLE_IS_32BITS
1164                                         dword0(rv) = Tiny0;
1165                                         dword1(rv) = Tiny1;
1166 #else
1167                                         dword0(rv) = Tiny1;
1168 #endif /*_DOUBLE_IS_32BITS*/
1169                                         goto cont;
1170                                         }
1171                                 else
1172                                         dword0(rv) -= P*Exp_msk1;
1173                                 }
1174                         else {
1175                                 adj = dval(aadj1) * ulp(dval(rv));
1176                                 dval(rv) += adj;
1177                                 }
1178 #else /*Sudden_Underflow*/
1179                         /* Compute adj so that the IEEE rounding rules will
1180                          * correctly round rv + adj in some half-way cases.
1181                          * If rv * ulp(rv) is denormalized (i.e.,
1182                          * y <= (P-1)*Exp_msk1), we must adjust aadj to avoid
1183                          * trouble from bits lost to denormalization;
1184                          * example: 1.2e-307 .
1185                          */
1186                         if (y <= (P-1)*Exp_msk1 && aadj > 1.) {
1187                                 dval(aadj1) = (double)(int)(aadj + 0.5);
1188                                 if (!dsign)
1189                                         dval(aadj1) = -dval(aadj1);
1190                                 }
1191                         adj = dval(aadj1) * ulp(dval(rv));
1192                         dval(rv) += adj;
1193 #endif /*Sudden_Underflow*/
1194 #endif /*Avoid_Underflow*/
1195                         }
1196                 z = dword0(rv) & Exp_mask;
1197 #ifndef SET_INEXACT
1198 #ifdef Avoid_Underflow
1199                 if (!scale)
1200 #endif
1201                 if (y == z) {
1202                         /* Can we stop now? */
1203 #ifndef _DOUBLE_IS_32BITS
1204                         /* If FE_INVALID floating point exceptions are
1205                            enabled, a conversion to a 32 bit value is
1206                            dangerous.  A positive double value can result
1207                            in a negative 32 bit int, thus raising SIGFPE.
1208                            To avoid this, always convert into 64 bit here. */
1209                         __int64_t L = (__int64_t)aadj;
1210 #else
1211                         L = (Long)aadj;
1212 #endif
1213                         aadj -= L;
1214                         /* The tolerances below are conservative. */
1215                         if (dsign || dword1(rv) || dword0(rv) & Bndry_mask) {
1216                                 if (aadj < .4999999 || aadj > .5000001)
1217                                         break;
1218                                 }
1219                         else if (aadj < .4999999/FLT_RADIX)
1220                                 break;
1221                         }
1222 #endif
1223  cont:
1224                 Bfree(ptr,bb);
1225                 Bfree(ptr,bd);
1226                 Bfree(ptr,bs);
1227                 Bfree(ptr,delta);
1228                 }
1229 #ifdef SET_INEXACT
1230         if (inexact) {
1231                 if (!oldinexact) {
1232                         dword0(rv0) = Exp_1 + (70 << Exp_shift);
1233 #ifndef _DOUBLE_IS_32BITS
1234                         dword1(rv0) = 0;
1235 #endif /*!_DOUBLE_IS_32BITS*/
1236                         dval(rv0) += 1.;
1237                         }
1238                 }
1239         else if (!oldinexact)
1240                 clear_inexact();
1241 #endif
1242 #ifdef Avoid_Underflow
1243         if (scale) {
1244                 dword0(rv0) = Exp_1 - 2*P*Exp_msk1;
1245 #ifndef _DOUBLE_IS_32BITS
1246                 dword1(rv0) = 0;
1247 #endif /*!_DOUBLE_IS_32BITS*/
1248                 dval(rv) *= dval(rv0);
1249 #ifndef NO_ERRNO
1250                 /* try to avoid the bug of testing an 8087 register value */
1251                 if ((dword0(rv) & Exp_mask) == 0)
1252                         _REENT_ERRNO(ptr) = ERANGE;
1253 #endif
1254                 }
1255 #endif /* Avoid_Underflow */
1256 #ifdef SET_INEXACT
1257         if (inexact && !(dword0(rv) & Exp_mask)) {
1258                 /* set underflow bit */
1259                 dval(rv0) = 1e-300;
1260                 dval(rv0) *= dval(rv0);
1261                 }
1262 #endif
1263  retfree:
1264         Bfree(ptr,bb);
1265         Bfree(ptr,bd);
1266         Bfree(ptr,bs);
1267         Bfree(ptr,bd0);
1268         Bfree(ptr,delta);
1269  ret:
1270         if (se)
1271                 *se = (char *)s;
1272         return sign ? -dval(rv) : dval(rv);
1273 }
1274
1275 double
1276 _strtod_r (struct _reent *ptr,
1277         const char *__restrict s00,
1278         char **__restrict se)
1279 {
1280   return _strtod_l (ptr, s00, se, __get_current_locale ());
1281 }
1282
1283 #ifndef _REENT_ONLY
1284
1285 double
1286 strtod_l (const char *__restrict s00, char **__restrict se, locale_t loc)
1287 {
1288   return _strtod_l (_REENT, s00, se, loc);
1289 }
1290
1291 double
1292 strtod (const char *__restrict s00, char **__restrict se)
1293 {
1294   return _strtod_l (_REENT, s00, se, __get_current_locale ());
1295 }
1296
1297 float
1298 strtof_l (const char *__restrict s00, char **__restrict se, locale_t loc)
1299 {
1300   double val = _strtod_l (_REENT, s00, se, loc);
1301   if (isnan (val))
1302     return signbit (val) ? -nanf ("") : nanf ("");
1303   float retval = (float) val;
1304 #ifndef NO_ERRNO
1305   if (isinf (retval) && !isinf (val))
1306     _REENT_ERRNO(_REENT) = ERANGE;
1307 #endif
1308   return retval;
1309 }
1310
1311 /*
1312  * These two functions are not quite correct as they return true for
1313  * zero, however they are 'good enough' for the test in strtof below
1314  * as we only need to know whether the double test is false when
1315  * the float test is true.
1316  */
1317 static inline int
1318 isdenorm(double d)
1319 {
1320     U u;
1321     dval(u) = d;
1322     return (dword0(u) & Exp_mask) == 0;
1323 }
1324
1325 static inline int
1326 isdenormf(float f)
1327 {
1328     union { float f; __uint32_t i; } u;
1329     u.f = f;
1330     return (u.i & 0x7f800000) == 0;
1331 }
1332
1333 float
1334 strtof (const char *__restrict s00,
1335         char **__restrict se)
1336 {
1337   double val = _strtod_l (_REENT, s00, se, __get_current_locale ());
1338   if (isnan (val))
1339     return signbit (val) ? -nanf ("") : nanf ("");
1340   float retval = (float) val;
1341 #ifndef NO_ERRNO
1342   if ((isinf (retval) && !isinf (val)) || (isdenormf(retval) && !isdenorm(val)))
1343     _REENT_ERRNO(_REENT) = ERANGE;
1344 #endif
1345   return retval;
1346 }
1347
1348 #endif