www.pudn.com > nasm-src.rar > float.c


/* float.c     floating-point constant support for the Netwide Assembler 
 * 
 * The Netwide Assembler is copyright (C) 1996 Simon Tatham and 
 * Julian Hall. All rights reserved. The software is 
 * redistributable under the licence given in the file "Licence" 
 * distributed in the NASM archive. 
 * 
 * initial version 13/ix/96 by Simon Tatham 
 */ 
 
#include  
#include  
#include  
 
#include "nasm.h" 
 
#define TRUE 1 
#define FALSE 0 
 
#define MANT_WORDS 6            /* 64 bits + 32 for accuracy == 96 */ 
#define MANT_DIGITS 28          /* 29 digits don't fit in 96 bits */ 
 
/* 
 * guaranteed top bit of from is set 
 * => we only have to worry about _one_ bit shift to the left 
 */ 
 
static int ieee_multiply(unsigned short *to, unsigned short *from) 
{ 
    unsigned long temp[MANT_WORDS * 2]; 
    int i, j; 
 
    for (i = 0; i < MANT_WORDS * 2; i++) 
        temp[i] = 0; 
 
    for (i = 0; i < MANT_WORDS; i++) 
        for (j = 0; j < MANT_WORDS; j++) { 
            unsigned long n; 
            n = (unsigned long)to[i] * (unsigned long)from[j]; 
            temp[i + j] += n >> 16; 
            temp[i + j + 1] += n & 0xFFFF; 
        } 
 
    for (i = MANT_WORDS * 2; --i;) { 
        temp[i - 1] += temp[i] >> 16; 
        temp[i] &= 0xFFFF; 
    } 
    if (temp[0] & 0x8000) { 
        for (i = 0; i < MANT_WORDS; i++) 
            to[i] = temp[i] & 0xFFFF; 
        return 0; 
    } else { 
        for (i = 0; i < MANT_WORDS; i++) 
            to[i] = (temp[i] << 1) + !!(temp[i + 1] & 0x8000); 
        return -1; 
    } 
} 
 
static void ieee_flconvert(char *string, unsigned short *mant, 
                           long *exponent, efunc error) 
{ 
    char digits[MANT_DIGITS]; 
    char *p, *q, *r; 
    unsigned short mult[MANT_WORDS], bit; 
    unsigned short *m; 
    long tenpwr, twopwr; 
    int extratwos, started, seendot; 
 
    p = digits; 
    tenpwr = 0; 
    started = seendot = FALSE; 
    while (*string && *string != 'E' && *string != 'e') { 
        if (*string == '.') { 
            if (!seendot) 
                seendot = TRUE; 
            else { 
                error(ERR_NONFATAL, 
                      "too many periods in floating-point constant"); 
                return; 
            } 
        } else if (*string >= '0' && *string <= '9') { 
            if (*string == '0' && !started) { 
                if (seendot) 
                    tenpwr--; 
            } else { 
                started = TRUE; 
                if (p < digits + sizeof(digits)) 
                    *p++ = *string - '0'; 
                if (!seendot) 
                    tenpwr++; 
            } 
        } else { 
            error(ERR_NONFATAL, 
                  "floating-point constant: `%c' is invalid character", 
                  *string); 
            return; 
        } 
        string++; 
    } 
    if (*string) { 
        string++;               /* eat the E */ 
        tenpwr += atoi(string); 
    } 
 
    /* 
     * At this point, the memory interval [digits,p) contains a 
     * series of decimal digits zzzzzzz such that our number X 
     * satisfies 
     * 
     * X = 0.zzzzzzz * 10^tenpwr 
     */ 
 
    bit = 0x8000; 
    for (m = mant; m < mant + MANT_WORDS; m++) 
        *m = 0; 
    m = mant; 
    q = digits; 
    started = FALSE; 
    twopwr = 0; 
    while (m < mant + MANT_WORDS) { 
        unsigned short carry = 0; 
        while (p > q && !p[-1]) 
            p--; 
        if (p <= q) 
            break; 
        for (r = p; r-- > q;) { 
            int i; 
 
            i = 2 * *r + carry; 
            if (i >= 10) 
                carry = 1, i -= 10; 
            else 
                carry = 0; 
            *r = i; 
        } 
        if (carry) 
            *m |= bit, started = TRUE; 
        if (started) { 
            if (bit == 1) 
                bit = 0x8000, m++; 
            else 
                bit >>= 1; 
        } else 
            twopwr--; 
    } 
    twopwr += tenpwr; 
 
    /* 
     * At this point the `mant' array contains the first six 
     * fractional places of a base-2^16 real number, which when 
     * multiplied by 2^twopwr and 5^tenpwr gives X. So now we 
     * really do multiply by 5^tenpwr. 
     */ 
 
    if (tenpwr < 0) { 
        for (m = mult; m < mult + MANT_WORDS; m++) 
            *m = 0xCCCC; 
        extratwos = -2; 
        tenpwr = -tenpwr; 
    } else if (tenpwr > 0) { 
        mult[0] = 0xA000; 
        for (m = mult + 1; m < mult + MANT_WORDS; m++) 
            *m = 0; 
        extratwos = 3; 
    } else 
        extratwos = 0; 
    while (tenpwr) { 
        if (tenpwr & 1) 
            twopwr += extratwos + ieee_multiply(mant, mult); 
        extratwos = extratwos * 2 + ieee_multiply(mult, mult); 
        tenpwr >>= 1; 
    } 
 
    /* 
     * Conversion is done. The elements of `mant' contain the first 
     * fractional places of a base-2^16 real number in [0.5,1) 
     * which we can multiply by 2^twopwr to get X. Or, of course, 
     * it contains zero. 
     */ 
    *exponent = twopwr; 
} 
 
/* 
 * Shift a mantissa to the right by i (i < 16) bits. 
 */ 
static void ieee_shr(unsigned short *mant, int i) 
{ 
    unsigned short n = 0, m; 
    int j; 
 
    for (j = 0; j < MANT_WORDS; j++) { 
        m = (mant[j] << (16 - i)) & 0xFFFF; 
        mant[j] = (mant[j] >> i) | n; 
        n = m; 
    } 
} 
 
/* 
 * Round a mantissa off after i words. 
 */ 
static int ieee_round(unsigned short *mant, int i) 
{ 
    if (mant[i] & 0x8000) { 
        do { 
            ++mant[--i]; 
            mant[i] &= 0xFFFF; 
        } while (i > 0 && !mant[i]); 
        return !i && !mant[i]; 
    } 
    return 0; 
} 
 
#define put(a,b) ( (*(a)=(b)), ((a)[1]=(b)>>8) ) 
 
static int to_double(char *str, long sign, unsigned char *result, 
                     efunc error) 
{ 
    unsigned short mant[MANT_WORDS]; 
    long exponent; 
 
    sign = (sign < 0 ? 0x8000L : 0L); 
 
    ieee_flconvert(str, mant, &exponent, error); 
    if (mant[0] & 0x8000) { 
        /* 
         * Non-zero. 
         */ 
        exponent--; 
        if (exponent >= -1022 && exponent <= 1024) { 
            /* 
             * Normalised. 
             */ 
            exponent += 1023; 
            ieee_shr(mant, 11); 
            ieee_round(mant, 4); 
            if (mant[0] & 0x20) /* did we scale up by one? */ 
                ieee_shr(mant, 1), exponent++; 
            mant[0] &= 0xF;     /* remove leading one */ 
            put(result + 6, (exponent << 4) | mant[0] | sign); 
            put(result + 4, mant[1]); 
            put(result + 2, mant[2]); 
            put(result + 0, mant[3]); 
        } else if (exponent < -1022 && exponent >= -1074) { 
            /* 
             * Denormal. 
             */ 
            int shift = -(exponent + 1011); 
            int sh = shift % 16, wds = shift / 16; 
            ieee_shr(mant, sh); 
            if (ieee_round(mant, 4 - wds) 
                || (sh > 0 && (mant[0] & (0x8000 >> (sh - 1))))) { 
                ieee_shr(mant, 1); 
                if (sh == 0) 
                    mant[0] |= 0x8000; 
                exponent++; 
            } 
            put(result + 6, (wds == 0 ? mant[0] : 0) | sign); 
            put(result + 4, (wds <= 1 ? mant[1 - wds] : 0)); 
            put(result + 2, (wds <= 2 ? mant[2 - wds] : 0)); 
            put(result + 0, (wds <= 3 ? mant[3 - wds] : 0)); 
        } else { 
            if (exponent > 0) { 
                error(ERR_NONFATAL, "overflow in floating-point constant"); 
                return 0; 
            } else 
                memset(result, 0, 8); 
        } 
    } else { 
        /* 
         * Zero. 
         */ 
        memset(result, 0, 8); 
    } 
    return 1;                   /* success */ 
} 
 
static int to_float(char *str, long sign, unsigned char *result, 
                    efunc error) 
{ 
    unsigned short mant[MANT_WORDS]; 
    long exponent; 
 
    sign = (sign < 0 ? 0x8000L : 0L); 
 
    ieee_flconvert(str, mant, &exponent, error); 
    if (mant[0] & 0x8000) { 
        /* 
         * Non-zero. 
         */ 
        exponent--; 
        if (exponent >= -126 && exponent <= 128) { 
            /* 
             * Normalised. 
             */ 
            exponent += 127; 
            ieee_shr(mant, 8); 
            ieee_round(mant, 2); 
            if (mant[0] & 0x100)        /* did we scale up by one? */ 
                ieee_shr(mant, 1), exponent++; 
            mant[0] &= 0x7F;    /* remove leading one */ 
            put(result + 2, (exponent << 7) | mant[0] | sign); 
            put(result + 0, mant[1]); 
        } else if (exponent < -126 && exponent >= -149) { 
            /* 
             * Denormal. 
             */ 
            int shift = -(exponent + 118); 
            int sh = shift % 16, wds = shift / 16; 
            ieee_shr(mant, sh); 
            if (ieee_round(mant, 2 - wds) 
                || (sh > 0 && (mant[0] & (0x8000 >> (sh - 1))))) { 
                ieee_shr(mant, 1); 
                if (sh == 0) 
                    mant[0] |= 0x8000; 
                exponent++; 
            } 
            put(result + 2, (wds == 0 ? mant[0] : 0) | sign); 
            put(result + 0, (wds <= 1 ? mant[1 - wds] : 0)); 
        } else { 
            if (exponent > 0) { 
                error(ERR_NONFATAL, "overflow in floating-point constant"); 
                return 0; 
            } else 
                memset(result, 0, 4); 
        } 
    } else { 
        memset(result, 0, 4); 
    } 
    return 1; 
} 
 
static int to_ldoub(char *str, long sign, unsigned char *result, 
                    efunc error) 
{ 
    unsigned short mant[MANT_WORDS]; 
    long exponent; 
 
    sign = (sign < 0 ? 0x8000L : 0L); 
 
    ieee_flconvert(str, mant, &exponent, error); 
    if (mant[0] & 0x8000) { 
        /* 
         * Non-zero. 
         */ 
        exponent--; 
        if (exponent >= -16383 && exponent <= 16384) { 
            /* 
             * Normalised. 
             */ 
            exponent += 16383; 
            if (ieee_round(mant, 4))    /* did we scale up by one? */ 
                ieee_shr(mant, 1), mant[0] |= 0x8000, exponent++; 
            put(result + 8, exponent | sign); 
            put(result + 6, mant[0]); 
            put(result + 4, mant[1]); 
            put(result + 2, mant[2]); 
            put(result + 0, mant[3]); 
        } else if (exponent < -16383 && exponent >= -16446) { 
            /* 
             * Denormal. 
             */ 
            int shift = -(exponent + 16383); 
            int sh = shift % 16, wds = shift / 16; 
            ieee_shr(mant, sh); 
            if (ieee_round(mant, 4 - wds) 
                || (sh > 0 && (mant[0] & (0x8000 >> (sh - 1))))) { 
                ieee_shr(mant, 1); 
                if (sh == 0) 
                    mant[0] |= 0x8000; 
                exponent++; 
            } 
            put(result + 8, sign); 
            put(result + 6, (wds == 0 ? mant[0] : 0)); 
            put(result + 4, (wds <= 1 ? mant[1 - wds] : 0)); 
            put(result + 2, (wds <= 2 ? mant[2 - wds] : 0)); 
            put(result + 0, (wds <= 3 ? mant[3 - wds] : 0)); 
        } else { 
            if (exponent > 0) { 
                error(ERR_NONFATAL, "overflow in floating-point constant"); 
                return 0; 
            } else 
                memset(result, 0, 10); 
        } 
    } else { 
        /* 
         * Zero. 
         */ 
        memset(result, 0, 10); 
    } 
    return 1; 
} 
 
int float_const(char *number, long sign, unsigned char *result, int bytes, 
                efunc error) 
{ 
    if (bytes == 4) 
        return to_float(number, sign, result, error); 
    else if (bytes == 8) 
        return to_double(number, sign, result, error); 
    else if (bytes == 10) 
        return to_ldoub(number, sign, result, error); 
    else { 
        error(ERR_PANIC, "strange value %d passed to float_const", bytes); 
        return 0; 
    } 
}