float.c

   1 /* float.c     floating-point constant support for the Netwide Assembler
   2  *
   3  * The Netwide Assembler is copyright (C) 1996 Simon Tatham and
   4  * Julian Hall. All rights reserved. The software is
   5  * redistributable under the license given in the file "LICENSE"
   6  * distributed in the NASM archive.
   7  *
   8  * initial version 13/ix/96 by Simon Tatham
   9  */
  10
  11 #include "compiler.h"
  12
  13 #include <ctype.h>
  14 #include <stdio.h>
  15 #include <stdlib.h>
  16 #include <string.h>
  17 #include <inttypes.h>
  18
  19 #include "nasm.h"
  20 #include "float.h"
  21
  22 /*
  23  * -----------------
  24  *  local variables
  25  * -----------------
  26  */
  27 static efunc error;
  28 static bool daz = false;        /* denormals as zero */
  29 static enum float_round rc = FLOAT_RC_NEAR;     /* rounding control */
  30
  31 /*
  32  * -----------
  33  *  constants
  34  * -----------
  35  */
  36
  37 /* "A limb is like a digit but bigger */
  38 typedef uint32_t fp_limb;
  39 typedef uint64_t fp_2limb;
  40
  41 #define LIMB_BITS       32
  42 #define LIMB_BYTES      (LIMB_BITS/8)
  43 #define LIMB_TOP_BIT    ((fp_limb)1 << (LIMB_BITS-1))
  44 #define LIMB_MASK       ((fp_limb)(~0))
  45 #define LIMB_ALL_BYTES  ((fp_limb)0x01010101)
  46 #define LIMB_BYTE(x)    ((x)*LIMB_ALL_BYTES)
  47
  48 #if X86_MEMORY
  49 #define put(a,b) (*(uint32_t *)(a) = (b))
  50 #else
  51 #define put(a,b) (((a)[0] = (b)),               \
  52                   ((a)[1] = (b) >> 8),          \
  53                   ((a)[2] = (b) >> 16),         \
  54                   ((a)[3] = (b) >> 24))
  55 #endif
  56
  57 /* 112 bits + 64 bits for accuracy + 16 bits for rounding */
  58 #define MANT_LIMBS 6
  59
  60 /* 52 digits fit in 176 bits because 10^53 > 2^176 > 10^52 */
  61 #define MANT_DIGITS 52
  62
  63 /* the format and the argument list depend on MANT_LIMBS */
  64 #define MANT_FMT "%08x_%08x_%08x_%08x_%08x_%08x"
  65 #define MANT_ARG SOME_ARG(mant, 0)
  66
  67 #define SOME_ARG(a,i) (a)[(i)+0], (a)[(i)+1], (a)[(i)+2], (a)[(i)+3],   \
  68         (a)[(i)+4], (a)[(i)+5]
  69
  70 /*
  71  * ---------------------------------------------------------------------------
  72  *  emit a printf()-like debug message... but only if DEBUG_FLOAT was defined
  73  * ---------------------------------------------------------------------------
  74  */
  75
  76 #ifdef DEBUG_FLOAT
  77 #define dprintf(x) printf x
  78 #else                           /*  */
  79 #define dprintf(x) do { } while (0)
  80 #endif                          /*  */
  81
  82 /*
  83  * ---------------------------------------------------------------------------
  84  *  multiply
  85  * ---------------------------------------------------------------------------
  86  */
  87 static int float_multiply(fp_limb *to, fp_limb *from)
  88 {
  89     fp_2limb temp[MANT_LIMBS * 2];
  90     int i, j;
  91
  92     /*
  93      * guaranteed that top bit of 'from' is set -- so we only have
  94      * to worry about _one_ bit shift to the left
  95      */
  96     dprintf(("%s=" MANT_FMT "\n", "mul1", SOME_ARG(to, 0)));
  97     dprintf(("%s=" MANT_FMT "\n", "mul2", SOME_ARG(from, 0)));
  98
  99     memset(temp, 0, sizeof temp);
 100
 101     for (i = 0; i < MANT_LIMBS; i++) {
 102         for (j = 0; j < MANT_LIMBS; j++) {
 103             fp_2limb n;
 104             n = (fp_2limb) to[i] * (fp_2limb) from[j];
 105             temp[i + j] += n >> LIMB_BITS;
 106             temp[i + j + 1] += (fp_limb)n;
 107         }
 108     }
 109
 110     for (i = MANT_LIMBS * 2; --i;) {
 111         temp[i - 1] += temp[i] >> LIMB_BITS;
 112         temp[i] &= LIMB_MASK;
 113     }
 114
 115     dprintf(("%s=" MANT_FMT "_" MANT_FMT "\n", "temp", SOME_ARG(temp, 0),
 116              SOME_ARG(temp, MANT_LIMBS)));
 117
 118     if (temp[0] & LIMB_TOP_BIT) {
 119         for (i = 0; i < MANT_LIMBS; i++) {
 120             to[i] = temp[i] & LIMB_MASK;
 121         }
 122         dprintf(("%s=" MANT_FMT " (%i)\n", "prod", SOME_ARG(to, 0), 0));
 123         return 0;
 124     } else {
 125         for (i = 0; i < MANT_LIMBS; i++) {
 126             to[i] = (temp[i] << 1) + !!(temp[i + 1] & LIMB_TOP_BIT);
 127         }
 128         dprintf(("%s=" MANT_FMT " (%i)\n", "prod", SOME_ARG(to, 0), -1));
 129         return -1;
 130     }
 131 }
 132
 133 /*
 134  * ---------------------------------------------------------------------------
 135  *  read an exponent; returns INT32_MAX on error
 136  * ---------------------------------------------------------------------------
 137  */
 138 static int32_t read_exponent(const char *string, int32_t max)
 139 {
 140     int32_t i = 0;
 141     bool neg = false;
 142
 143     if (*string == '+') {
 144         string++;
 145     } else if (*string == '-') {
 146         neg = true;
 147         string++;
 148     }
 149     while (*string) {
 150         if (*string >= '0' && *string <= '9') {
 151             i = (i * 10) + (*string - '0');
 152
 153             /*
 154              * To ensure that underflows and overflows are
 155              * handled properly we must avoid wraparounds of
 156              * the signed integer value that is used to hold
 157              * the exponent. Therefore we cap the exponent at
 158              * +/-5000, which is slightly more/less than
 159              * what's required for normal and denormal numbers
 160              * in single, double, and extended precision, but
 161              * sufficient to avoid signed integer wraparound.
 162              */
 163             if (i > max)
 164                 i = max;
 165         } else if (*string == '_') {
 166             /* do nothing */
 167         } else {
 168             error(ERR_NONFATAL|ERR_PASS1,
 169                   "invalid character in floating-point constant %s: '%c'",
 170                   "exponent", *string);
 171             return INT32_MAX;
 172         }
 173         string++;
 174     }
 175
 176     return neg ? -i : i;
 177 }
 178
 179 /*
 180  * ---------------------------------------------------------------------------
 181  *  convert
 182  * ---------------------------------------------------------------------------
 183  */
 184 static bool ieee_flconvert(const char *string, fp_limb *mant,
 185                            int32_t * exponent)
 186 {
 187     char digits[MANT_DIGITS];
 188     char *p, *q, *r;
 189     fp_limb mult[MANT_LIMBS], bit;
 190     fp_limb *m;
 191     int32_t tenpwr, twopwr;
 192     int32_t extratwos;
 193     bool started, seendot, warned;
 194
 195     warned = false;
 196     p = digits;
 197     tenpwr = 0;
 198     started = seendot = false;
 199
 200     while (*string && *string != 'E' && *string != 'e') {
 201         if (*string == '.') {
 202             if (!seendot) {
 203                 seendot = true;
 204             } else {
 205                 error(ERR_NONFATAL|ERR_PASS1,
 206                       "too many periods in floating-point constant");
 207                 return false;
 208             }
 209         } else if (*string >= '0' && *string <= '9') {
 210             if (*string == '0' && !started) {
 211                 if (seendot) {
 212                     tenpwr--;
 213                 }
 214             } else {
 215                 started = true;
 216                 if (p < digits + sizeof(digits)) {
 217                     *p++ = *string - '0';
 218                 } else {
 219                     if (!warned) {
 220                         error(ERR_WARNING|ERR_WARN_FL_TOOLONG|ERR_PASS1,
 221                               "floating-point constant significand contains "
 222                               "more than %i digits", MANT_DIGITS);
 223                         warned = true;
 224                     }
 225                 }
 226                 if (!seendot) {
 227                     tenpwr++;
 228                 }
 229             }
 230         } else if (*string == '_') {
 231             /* do nothing */
 232         } else {
 233             error(ERR_NONFATAL|ERR_PASS1,
 234                   "invalid character in floating-point constant %s: '%c'",
 235                   "significand", *string);
 236             return false;
 237         }
 238         string++;
 239     }
 240
 241     if (*string) {
 242         int32_t e;
 243
 244         string++;               /* eat the E */
 245         e = read_exponent(string, 5000);
 246         if (e == INT32_MAX)
 247             return false;
 248         tenpwr += e;
 249     }
 250
 251     /*
 252      * At this point, the memory interval [digits,p) contains a
 253      * series of decimal digits zzzzzzz, such that our number X
 254      * satisfies X = 0.zzzzzzz * 10^tenpwr.
 255      */
 256     q = digits;
 257     dprintf(("X = 0."));
 258     while (q < p) {
 259         dprintf(("%c", *q + '0'));
 260         q++;
 261     }
 262     dprintf((" * 10^%i\n", tenpwr));
 263
 264     /*
 265      * Now convert [digits,p) to our internal representation.
 266      */
 267     bit = LIMB_TOP_BIT;
 268     for (m = mant; m < mant + MANT_LIMBS; m++) {
 269         *m = 0;
 270     }
 271     m = mant;
 272     q = digits;
 273     started = false;
 274     twopwr = 0;
 275     while (m < mant + MANT_LIMBS) {
 276         fp_limb carry = 0;
 277         while (p > q && !p[-1]) {
 278             p--;
 279         }
 280         if (p <= q) {
 281             break;
 282         }
 283         for (r = p; r-- > q;) {
 284             int32_t i;
 285             i = 2 * *r + carry;
 286             if (i >= 10) {
 287                 carry = 1;
 288                 i -= 10;
 289             } else {
 290                 carry = 0;
 291             }
 292             *r = i;
 293         }
 294         if (carry) {
 295             *m |= bit;
 296             started = true;
 297         }
 298         if (started) {
 299             if (bit == 1) {
 300                 bit = LIMB_TOP_BIT;
 301                 m++;
 302             } else {
 303                 bit >>= 1;
 304             }
 305         } else {
 306             twopwr--;
 307         }
 308     }
 309     twopwr += tenpwr;
 310
 311     /*
 312      * At this point, the 'mant' array contains the first frac-
 313      * tional places of a base-2^16 real number which when mul-
 314      * tiplied by 2^twopwr and 5^tenpwr gives X.
 315      */
 316     dprintf(("X = " MANT_FMT " * 2^%i * 5^%i\n", MANT_ARG, twopwr,
 317              tenpwr));
 318
 319     /*
 320      * Now multiply 'mant' by 5^tenpwr.
 321      */
 322     if (tenpwr < 0) {           /* mult = 5^-1 = 0.2 */
 323         for (m = mult; m < mult + MANT_LIMBS - 1; m++) {
 324             *m = LIMB_BYTE(0xcc);
 325         }
 326         mult[MANT_LIMBS - 1] = LIMB_BYTE(0xcc)+1;
 327         extratwos = -2;
 328         tenpwr = -tenpwr;
 329
 330         /*
 331          * If tenpwr was 1000...000b, then it becomes 1000...000b. See
 332          * the "ANSI C" comment below for more details on that case.
 333          *
 334          * Because we already truncated tenpwr to +5000...-5000 inside
 335          * the exponent parsing code, this shouldn't happen though.
 336          */
 337     } else if (tenpwr > 0) {    /* mult = 5^+1 = 5.0 */
 338         mult[0] = (fp_limb)5 << (LIMB_BITS-3); /* 0xA000... */
 339         for (m = mult + 1; m < mult + MANT_LIMBS; m++) {
 340             *m = 0;
 341         }
 342         extratwos = 3;
 343     } else {
 344         extratwos = 0;
 345     }
 346     while (tenpwr) {
 347         dprintf(("loop=" MANT_FMT " * 2^%i * 5^%i (%i)\n", MANT_ARG,
 348                  twopwr, tenpwr, extratwos));
 349         if (tenpwr & 1) {
 350             dprintf(("mant*mult\n"));
 351             twopwr += extratwos + float_multiply(mant, mult);
 352         }
 353         dprintf(("mult*mult\n"));
 354         extratwos = extratwos * 2 + float_multiply(mult, mult);
 355         tenpwr >>= 1;
 356
 357         /*
 358          * In ANSI C, the result of right-shifting a signed integer is
 359          * considered implementation-specific. To ensure that the loop
 360          * terminates even if tenpwr was 1000...000b to begin with, we
 361          * manually clear the MSB, in case a 1 was shifted in.
 362          *
 363          * Because we already truncated tenpwr to +5000...-5000 inside
 364          * the exponent parsing code, this shouldn't matter; neverthe-
 365          * less it is the right thing to do here.
 366          */
 367         tenpwr &= (uint32_t) - 1 >> 1;
 368     }
 369
 370     /*
 371      * At this point, the 'mant' array contains the first frac-
 372      * tional places of a base-2^16 real number in [0.5,1) that
 373      * when multiplied by 2^twopwr gives X. Or it contains zero
 374      * of course. We are done.
 375      */
 376     *exponent = twopwr;
 377     return true;
 378 }
 379
 380 /*
 381  * ---------------------------------------------------------------------------
 382  *  operations of specific bits
 383  * ---------------------------------------------------------------------------
 384  */
 385
 386 /* Set a bit, using *bigendian* bit numbering (0 = MSB) */
 387 static void set_bit(fp_limb *mant, int bit)
 388 {
 389     mant[bit/LIMB_BITS] |= LIMB_TOP_BIT >> (bit & (LIMB_BITS-1));
 390 }
 391
 392 /* Test a single bit */
 393 static int test_bit(const fp_limb *mant, int bit)
 394 {
 395     return (mant[bit/LIMB_BITS] >> (~bit & (LIMB_BITS-1))) & 1;
 396 }
 397
 398 /* Report if the mantissa value is all zero */
 399 static bool is_zero(const fp_limb *mant)
 400 {
 401     int i;
 402
 403     for (i = 0; i < MANT_LIMBS; i++)
 404         if (mant[i])
 405             return false;
 406
 407     return true;
 408 }
 409
 410 /*
 411  * ---------------------------------------------------------------------------
 412  *  round a mantissa off after i words
 413  * ---------------------------------------------------------------------------
 414  */
 415
 416 #define ROUND_COLLECT_BITS                      \
 417     do {                                        \
 418         m = mant[i] & (2*bit-1);                \
 419         for (j = i+1; j < MANT_LIMBS; j++)      \
 420             m = m | mant[j];                    \
 421     } while (0)
 422
 423 #define ROUND_ABS_DOWN                          \
 424     do {                                        \
 425         mant[i] &= ~(bit-1);                    \
 426         for (j = i+1; j < MANT_LIMBS; j++)      \
 427             mant[j] = 0;                        \
 428         return false;                           \
 429     } while (0)
 430
 431 #define ROUND_ABS_UP                            \
 432     do {                                        \
 433         mant[i] = (mant[i] & ~(bit-1)) + bit;   \
 434         for (j = i+1; j < MANT_LIMBS; j++)      \
 435             mant[j] = 0;                        \
 436         while (i > 0 && !mant[i])               \
 437             ++mant[--i];                        \
 438         return !mant[0];                        \
 439     } while (0)
 440
 441 static bool ieee_round(bool minus, fp_limb *mant, int bits)
 442 {
 443     fp_limb m = 0;
 444     int32_t j;
 445     int i = bits / LIMB_BITS;
 446     int p = bits % LIMB_BITS;
 447     fp_limb bit = LIMB_TOP_BIT >> p;
 448
 449     if (rc == FLOAT_RC_NEAR) {
 450         if (mant[i] & bit) {
 451             mant[i] &= ~bit;
 452             ROUND_COLLECT_BITS;
 453             mant[i] |= bit;
 454             if (m) {
 455                 ROUND_ABS_UP;
 456             } else {
 457                 if (test_bit(mant, bits-1)) {
 458                     ROUND_ABS_UP;
 459                 } else {
 460                     ROUND_ABS_DOWN;
 461                 }
 462             }
 463         } else {
 464             ROUND_ABS_DOWN;
 465         }
 466     } else if (rc == FLOAT_RC_ZERO ||
 467                rc == (minus ? FLOAT_RC_UP : FLOAT_RC_DOWN)) {
 468         ROUND_ABS_DOWN;
 469     } else {
 470         /* rc == (minus ? FLOAT_RC_DOWN : FLOAT_RC_UP) */
 471         /* Round toward +/- infinity */
 472         ROUND_COLLECT_BITS;
 473         if (m) {
 474             ROUND_ABS_UP;
 475         } else {
 476             ROUND_ABS_DOWN;
 477         }
 478     }
 479     return false;
 480 }
 481
 482 /* Returns a value >= 16 if not a valid hex digit */
 483 static unsigned int hexval(char c)
 484 {
 485     unsigned int v = (unsigned char) c;
 486
 487     if (v >= '0' && v <= '9')
 488         return v - '0';
 489     else
 490         return (v|0x20) - 'a' + 10;
 491 }
 492
 493 /* Handle floating-point numbers with radix 2^bits and binary exponent */
 494 static bool ieee_flconvert_bin(const char *string, int bits,
 495                                fp_limb *mant, int32_t *exponent)
 496 {
 497     static const int log2tbl[16] =
 498         { -1, 0, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3 };
 499     fp_limb mult[MANT_LIMBS + 1], *mp;
 500     int ms;
 501     int32_t twopwr;
 502     bool seendot, seendigit;
 503     unsigned char c;
 504     int radix = 1 << bits;
 505     fp_limb v;
 506
 507     twopwr = 0;
 508     seendot = seendigit = false;
 509     ms = 0;
 510     mp = NULL;
 511
 512     memset(mult, 0, sizeof mult);
 513
 514     while ((c = *string++) != '\0') {
 515         if (c == '.') {
 516             if (!seendot)
 517                 seendot = true;
 518             else {
 519                 error(ERR_NONFATAL|ERR_PASS1,
 520                       "too many periods in floating-point constant");
 521                 return false;
 522             }
 523         } else if ((v = hexval(c)) < (unsigned int)radix) {
 524             if (!seendigit && v) {
 525                 int l = log2tbl[v];
 526
 527                 seendigit = true;
 528                 mp = mult;
 529                 ms = (LIMB_BITS-1)-l;
 530
 531                 twopwr = seendot ? twopwr-bits+l : l+1-bits;
 532             }
 533
 534             if (seendigit) {
 535                 if (ms <= 0) {
 536                     *mp |= v >> -ms;
 537                     mp++;
 538                     if (mp > &mult[MANT_LIMBS])
 539                         mp = &mult[MANT_LIMBS]; /* Guard slot */
 540                     ms += LIMB_BITS;
 541                 }
 542                 *mp |= v << ms;
 543                 ms -= bits;
 544
 545                 if (!seendot)
 546                     twopwr += bits;
 547             } else {
 548                 if (seendot)
 549                     twopwr -= bits;
 550             }
 551         } else if (c == 'p' || c == 'P') {
 552             int32_t e;
 553             e = read_exponent(string, 20000);
 554             if (e == INT32_MAX)
 555                 return false;
 556             twopwr += e;
 557             break;
 558         } else if (c == '_') {
 559             /* ignore */
 560         } else {
 561             error(ERR_NONFATAL|ERR_PASS1,
 562                   "floating-point constant: `%c' is invalid character", c);
 563             return false;
 564         }
 565     }
 566
 567     if (!seendigit) {
 568         memset(mant, 0, sizeof mult); /* Zero */
 569         *exponent = 0;
 570     } else {
 571         memcpy(mant, mult, sizeof mult);
 572         *exponent = twopwr;
 573     }
 574
 575     return true;
 576 }
 577
 578 /*
 579  * Shift a mantissa to the right by i bits.
 580  */
 581 static void ieee_shr(fp_limb *mant, int i)
 582 {
 583     fp_limb n, m;
 584     int j = 0;
 585     int sr, sl, offs;
 586
 587     sr = i % LIMB_BITS; sl = LIMB_BITS-sr;
 588     offs = i/LIMB_BITS;
 589
 590     if (sr == 0) {
 591         if (offs)
 592             for (j = MANT_LIMBS-1; j >= offs; j--)
 593                 mant[j] = mant[j-offs];
 594     } else {
 595         n = mant[MANT_LIMBS-1-offs] >> sr;
 596         for (j = MANT_LIMBS-1; j > offs; j--) {
 597             m = mant[j-offs-1];
 598             mant[j] = (m << sl) | n;
 599             n = m >> sr;
 600         }
 601         mant[j--] = n;
 602     }
 603     while (j >= 0)
 604         mant[j--] = 0;
 605 }
 606
 607 /* Produce standard IEEE formats, with implicit or explicit integer
 608    bit; this makes the following assumptions:
 609
 610    - the sign bit is the MSB, followed by the exponent,
 611      followed by the integer bit if present.
 612    - the sign bit plus exponent fit in 16 bits.
 613    - the exponent bias is 2^(n-1)-1 for an n-bit exponent */
 614
 615 struct ieee_format {
 616     int bytes;
 617     int mantissa;               /* Fractional bits in the mantissa */
 618     int explicit;               /* Explicit integer */
 619     int exponent;               /* Bits in the exponent */
 620 };
 621
 622 /*
 623  * The 16- and 128-bit formats are expected to be in IEEE 754r.
 624  * AMD SSE5 uses the 16-bit format.
 625  *
 626  * The 32- and 64-bit formats are the original IEEE 754 formats.
 627  *
 628  * The 80-bit format is x87-specific, but widely used.
 629  *
 630  * The 8-bit format appears to be the consensus 8-bit floating-point
 631  * format.  It is apparently used in graphics applications.
 632  */
 633 static const struct ieee_format ieee_8   = {  1,   3, 0,  4 };
 634 static const struct ieee_format ieee_16  = {  2,  10, 0,  5 };
 635 static const struct ieee_format ieee_32  = {  4,  23, 0,  8 };
 636 static const struct ieee_format ieee_64  = {  8,  52, 0, 11 };
 637 static const struct ieee_format ieee_80  = { 10,  63, 1, 15 };
 638 static const struct ieee_format ieee_128 = { 16, 112, 0, 15 };
 639
 640 /* Types of values we can generate */
 641 enum floats {
 642     FL_ZERO,
 643     FL_DENORMAL,
 644     FL_NORMAL,
 645     FL_INFINITY,
 646     FL_QNAN,
 647     FL_SNAN
 648 };
 649
 650 static int to_float(const char *str, int s, uint8_t * result,
 651                     const struct ieee_format *fmt)
 652 {
 653     fp_limb mant[MANT_LIMBS], *mp, m;
 654     int32_t exponent = 0;
 655     int32_t expmax = 1 << (fmt->exponent - 1);
 656     fp_limb one_mask = LIMB_TOP_BIT >>
 657         ((fmt->exponent+fmt->explicit) % LIMB_BITS);
 658     int one_pos = (fmt->exponent+fmt->explicit)/LIMB_BITS;
 659     int i;
 660     int shift;
 661     enum floats type;
 662     bool ok;
 663     bool minus = s < 0;
 664     int bits = fmt->bytes * 8;
 665
 666     if (str[0] == '_') {
 667         /* Special tokens */
 668
 669         switch (str[2]) {
 670         case 'n':              /* __nan__ */
 671         case 'N':
 672         case 'q':              /* __qnan__ */
 673         case 'Q':
 674             type = FL_QNAN;
 675             break;
 676         case 's':              /* __snan__ */
 677         case 'S':
 678             type = FL_SNAN;
 679             break;
 680         case 'i':              /* __infinity__ */
 681         case 'I':
 682             type = FL_INFINITY;
 683             break;
 684         default:
 685             error(ERR_NONFATAL|ERR_PASS1,
 686                   "internal error: unknown FP constant token `%s'\n", str);
 687             type = FL_QNAN;
 688             break;
 689         }
 690     } else {
 691         if (str[0] == '0') {
 692             switch (str[1]) {
 693             case 'x': case 'X':
 694             case 'h': case 'H':
 695                 ok = ieee_flconvert_bin(str+2, 4, mant, &exponent);
 696                 break;
 697             case 'o': case 'O':
 698             case 'q': case 'Q':
 699                 ok = ieee_flconvert_bin(str+2, 3, mant, &exponent);
 700                 break;
 701             case 'b': case 'B':
 702             case 'y': case 'Y':
 703                 ok = ieee_flconvert_bin(str+2, 1, mant, &exponent);
 704                 break;
 705             case 'd': case 'D':
 706             case 't': case 'T':
 707                 ok = ieee_flconvert(str+2, mant, &exponent);
 708                 break;
 709             default:
 710                 /* Leading zero was just a zero? */
 711                 ok = ieee_flconvert(str, mant, &exponent);
 712                 break;
 713             }
 714         } else if (str[0] == '$') {
 715             ok = ieee_flconvert_bin(str+1, 4, mant, &exponent);
 716         } else {
 717             ok = ieee_flconvert(str, mant, &exponent);
 718         }
 719
 720         if (!ok) {
 721             type = FL_QNAN;
 722         } else if (mant[0] & LIMB_TOP_BIT) {
 723             /*
 724              * Non-zero.
 725              */
 726             exponent--;
 727             if (exponent >= 2 - expmax && exponent <= expmax) {
 728                 type = FL_NORMAL;
 729             } else if (exponent > 0) {
 730                 if (pass0 == 1)
 731                     error(ERR_WARNING|ERR_WARN_FL_OVERFLOW|ERR_PASS1,
 732                           "overflow in floating-point constant");
 733                 type = FL_INFINITY;
 734             } else {
 735                 /* underflow or denormal; the denormal code handles
 736                    actual underflow. */
 737                 type = FL_DENORMAL;
 738             }
 739         } else {
 740             /* Zero */
 741             type = FL_ZERO;
 742         }
 743     }
 744
 745     switch (type) {
 746     case FL_ZERO:
 747     zero:
 748         memset(mant, 0, sizeof mant);
 749         break;
 750
 751     case FL_DENORMAL:
 752     {
 753         shift = -(exponent + expmax - 2 - fmt->exponent)
 754             + fmt->explicit;
 755         ieee_shr(mant, shift);
 756         ieee_round(minus, mant, bits);
 757         if (mant[one_pos] & one_mask) {
 758             /* One's position is set, we rounded up into normal range */
 759             exponent = 1;
 760             if (!fmt->explicit)
 761                 mant[one_pos] &= ~one_mask;     /* remove explicit one */
 762             mant[0] |= exponent << (LIMB_BITS-1 - fmt->exponent);
 763         } else {
 764             if (daz || is_zero(mant)) {
 765                 /* Flush denormals to zero */
 766                 error(ERR_WARNING|ERR_WARN_FL_UNDERFLOW|ERR_PASS1,
 767                       "underflow in floating-point constant");
 768                 goto zero;
 769             } else {
 770                 error(ERR_WARNING|ERR_WARN_FL_DENORM|ERR_PASS1,
 771                       "denormal floating-point constant");
 772             }
 773         }
 774         break;
 775     }
 776
 777     case FL_NORMAL:
 778         exponent += expmax - 1;
 779         ieee_shr(mant, fmt->exponent+fmt->explicit);
 780         ieee_round(minus, mant, bits);
 781         /* did we scale up by one? */
 782         if (test_bit(mant, fmt->exponent+fmt->explicit-1)) {
 783             ieee_shr(mant, 1);
 784             exponent++;
 785             if (exponent >= (expmax << 1)-1) {
 786                     error(ERR_WARNING|ERR_WARN_FL_OVERFLOW|ERR_PASS1,
 787                           "overflow in floating-point constant");
 788                 type = FL_INFINITY;
 789                 goto overflow;
 790             }
 791         }
 792
 793         if (!fmt->explicit)
 794             mant[one_pos] &= ~one_mask; /* remove explicit one */
 795         mant[0] |= exponent << (LIMB_BITS-1 - fmt->exponent);
 796         break;
 797
 798     case FL_INFINITY:
 799     case FL_QNAN:
 800     case FL_SNAN:
 801     overflow:
 802         memset(mant, 0, sizeof mant);
 803         mant[0] = (((fp_limb)1 << fmt->exponent)-1)
 804             << (LIMB_BITS-1 - fmt->exponent);
 805         if (fmt->explicit)
 806             mant[one_pos] |= one_mask;
 807         if (type == FL_QNAN)
 808             set_bit(mant, fmt->exponent+fmt->explicit+1);
 809         else if (type == FL_SNAN)
 810             set_bit(mant, fmt->exponent+fmt->explicit+fmt->mantissa);
 811         break;
 812     }
 813
 814     mant[0] |= minus ? LIMB_TOP_BIT : 0;
 815
 816     m = mant[fmt->bytes/LIMB_BYTES];
 817     for (i = LIMB_BYTES-(fmt->bytes % LIMB_BYTES); i < LIMB_BYTES; i++)
 818         *result++ = m >> (i*8);
 819
 820     for (mp = &mant[fmt->bytes/LIMB_BYTES], i = 0;
 821          i < fmt->bytes; i += LIMB_BYTES) {
 822         m = *--mp;
 823         put(result, m);
 824         result += LIMB_BYTES;
 825     }
 826
 827     return 1;                   /* success */
 828 }
 829
 830 int float_const(const char *number, int sign, uint8_t * result,
 831                 int bytes, efunc err)
 832 {
 833     error = err;
 834
 835     switch (bytes) {
 836     case 1:
 837         return to_float(number, sign, result, &ieee_8);
 838     case 2:
 839         return to_float(number, sign, result, &ieee_16);
 840     case 4:
 841         return to_float(number, sign, result, &ieee_32);
 842     case 8:
 843         return to_float(number, sign, result, &ieee_64);
 844     case 10:
 845         return to_float(number, sign, result, &ieee_80);
 846     case 16:
 847         return to_float(number, sign, result, &ieee_128);
 848     default:
 849         error(ERR_PANIC, "strange value %d passed to float_const", bytes);
 850         return 0;
 851     }
 852 }
 853
 854 /* Set floating-point options */
 855 int float_option(const char *option)
 856 {
 857     if (!nasm_stricmp(option, "daz")) {
 858         daz = true;
 859         return 0;
 860     } else if (!nasm_stricmp(option, "nodaz")) {
 861         daz = false;
 862         return 0;
 863     } else if (!nasm_stricmp(option, "near")) {
 864         rc = FLOAT_RC_NEAR;
 865         return 0;
 866     } else if (!nasm_stricmp(option, "down")) {
 867         rc = FLOAT_RC_DOWN;
 868         return 0;
 869     } else if (!nasm_stricmp(option, "up")) {
 870         rc = FLOAT_RC_UP;
 871         return 0;
 872     } else if (!nasm_stricmp(option, "zero")) {
 873         rc = FLOAT_RC_ZERO;
 874         return 0;
 875     } else if (!nasm_stricmp(option, "default")) {
 876         rc = FLOAT_RC_NEAR;
 877         daz = false;
 878         return 0;
 879     } else {
 880         return -1;              /* Unknown option */
 881     }
 882 }