asm/float.c

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