1 /* float.c floating-point constant support for the Netwide Assembler
3 * The Netwide Assembler is copyright (C) 1996 Simon Tatham and
4 * Julian Hall. All rights reserved. The software is
5 * redistributable under the licence given in the file "Licence"
6 * distributed in the NASM archive.
8 * initial version 13/ix/96 by Simon Tatham
22 #define MANT_WORDS 10 /* 112 bits + 48 for accuracy == 160 */
23 #define MANT_DIGITS 49 /* 50 digits don't fit in 160 bits */
26 * guaranteed top bit of from is set
27 * => we only have to worry about _one_ bit shift to the left
30 static int ieee_multiply(uint16_t *to
, uint16_t *from
)
32 uint32_t temp
[MANT_WORDS
* 2];
35 for (i
= 0; i
< MANT_WORDS
* 2; i
++)
38 for (i
= 0; i
< MANT_WORDS
; i
++)
39 for (j
= 0; j
< MANT_WORDS
; j
++) {
41 n
= (uint32_t)to
[i
] * (uint32_t)from
[j
];
42 temp
[i
+ j
] += n
>> 16;
43 temp
[i
+ j
+ 1] += n
& 0xFFFF;
46 for (i
= MANT_WORDS
* 2; --i
;) {
47 temp
[i
- 1] += temp
[i
] >> 16;
50 if (temp
[0] & 0x8000) {
51 memcpy(to
, temp
, 2*MANT_WORDS
);
54 for (i
= 0; i
< MANT_WORDS
; i
++)
55 to
[i
] = (temp
[i
] << 1) + !!(temp
[i
+ 1] & 0x8000);
60 static int hexval(char c
)
62 if (c
>= '0' && c
<= '9')
64 else if (c
>= 'a' && c
<= 'f')
70 static void ieee_flconvert_hex(char *string
, uint16_t *mant
,
71 int32_t *exponent
, efunc error
)
73 static const int log2tbl
[16] =
74 { -1, 0, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3 };
75 uint16_t mult
[MANT_WORDS
+1], *mp
;
78 int seendot
, seendigit
;
82 seendot
= seendigit
= 0;
86 memset(mult
, 0, sizeof mult
);
88 while ((c
= *string
++) != '\0') {
94 "too many periods in floating-point constant");
97 } else if (isxdigit(c
)) {
100 if (!seendigit
&& v
) {
107 twopwr
= seendot
? twopwr
-4+l
: l
-3;
114 if (mp
> &mult
[MANT_WORDS
])
115 mp
= &mult
[MANT_WORDS
]; /* Guard slot */
127 } else if (c
== 'p' || c
== 'P') {
128 twopwr
+= atoi(string
);
132 "floating-point constant: `%c' is invalid character",
139 memset(mant
, 0, 2*MANT_WORDS
); /* Zero */
142 memcpy(mant
, mult
, 2*MANT_WORDS
);
147 static void ieee_flconvert(char *string
, uint16_t *mant
,
148 int32_t *exponent
, efunc error
)
150 char digits
[MANT_DIGITS
];
152 uint16_t mult
[MANT_WORDS
], bit
;
154 int32_t tenpwr
, twopwr
;
155 int extratwos
, started
, seendot
;
157 if (string
[0] == '0' && (string
[1] == 'x' || string
[1] == 'X')) {
158 ieee_flconvert_hex(string
+2, mant
, exponent
, error
);
164 started
= seendot
= FALSE
;
165 while (*string
&& *string
!= 'E' && *string
!= 'e') {
166 if (*string
== '.') {
171 "too many periods in floating-point constant");
174 } else if (*string
>= '0' && *string
<= '9') {
175 if (*string
== '0' && !started
) {
180 if (p
< digits
+ sizeof(digits
))
181 *p
++ = *string
- '0';
187 "floating-point constant: `%c' is invalid character",
194 string
++; /* eat the E */
195 tenpwr
+= atoi(string
);
199 * At this point, the memory interval [digits,p) contains a
200 * series of decimal digits zzzzzzz such that our number X
203 * X = 0.zzzzzzz * 10^tenpwr
207 for (m
= mant
; m
< mant
+ MANT_WORDS
; m
++)
213 while (m
< mant
+ MANT_WORDS
) {
215 while (p
> q
&& !p
[-1])
219 for (r
= p
; r
-- > q
;) {
230 *m
|= bit
, started
= TRUE
;
242 * At this point the `mant' array contains the first six
243 * fractional places of a base-2^16 real number, which when
244 * multiplied by 2^twopwr and 5^tenpwr gives X. So now we
245 * really do multiply by 5^tenpwr.
249 for (m
= mult
; m
< mult
+ MANT_WORDS
; m
++)
253 } else if (tenpwr
> 0) {
255 for (m
= mult
+ 1; m
< mult
+ MANT_WORDS
; m
++)
262 twopwr
+= extratwos
+ ieee_multiply(mant
, mult
);
263 extratwos
= extratwos
* 2 + ieee_multiply(mult
, mult
);
268 * Conversion is done. The elements of `mant' contain the first
269 * fractional places of a base-2^16 real number in [0.5,1)
270 * which we can multiply by 2^twopwr to get X. Or, of course,
277 * Shift a mantissa to the right by i (i < 16) bits.
279 static void ieee_shr(uint16_t *mant
, int i
)
284 for (j
= 0; j
< MANT_WORDS
; j
++) {
285 m
= (mant
[j
] << (16 - i
)) & 0xFFFF;
286 mant
[j
] = (mant
[j
] >> i
) | n
;
292 * Round a mantissa off after i words.
294 static int ieee_round(uint16_t *mant
, int i
)
296 if (mant
[i
] & 0x8000) {
300 } while (i
> 0 && !mant
[i
]);
301 return !i
&& !mant
[i
];
306 #define put(a,b) ( (*(a)=(b)), ((a)[1]=(b)>>8) )
308 /* Set a bit, using *bigendian* bit numbering (0 = MSB) */
309 static void set_bit(uint16_t *mant
, int bit
)
311 mant
[bit
>> 4] |= 1 << (~bit
& 15);
314 /* Produce standard IEEE formats, with implicit "1" bit; this makes
315 the following assumptions:
317 - the sign bit is the MSB, followed by the exponent.
318 - the sign bit plus exponent fit in 16 bits.
319 - the exponent bias is 2^(n-1)-1 for an n-bit exponent */
323 int mantissa
; /* Bits in the mantissa */
324 int exponent
; /* Bits in the exponent */
327 static const struct ieee_format ieee_16
= { 1, 10, 5 };
328 static const struct ieee_format ieee_32
= { 2, 23, 8 };
329 static const struct ieee_format ieee_64
= { 4, 52, 11 };
330 static const struct ieee_format ieee_128
= { 8, 112, 15 };
332 /* Produce all the standard IEEE formats: 16, 32, 64, and 128 bits */
333 static int to_float(char *str
, int32_t sign
, uint8_t *result
,
334 const struct ieee_format
*fmt
, efunc error
)
336 uint16_t mant
[MANT_WORDS
], *mp
;
338 int32_t expmax
= 1 << (fmt
->exponent
-1);
339 uint16_t implicit_one
= 0x8000 >> fmt
->exponent
;
342 sign
= (sign
< 0 ? 0x8000L
: 0L);
345 /* NaN or Infinity */
346 int32_t expmask
= (1 << fmt
->exponent
)-1;
348 memset(mant
, 0, sizeof mant
);
349 mant
[0] = expmask
<< (15-fmt
->exponent
); /* Exponent: all bits one */
352 case 'n': /* __nan__ */
354 case 'q': /* __qnan__ */
356 set_bit(mant
, fmt
->exponent
+1); /* Highest bit in mantissa */
358 case 's': /* __snan__ */
360 set_bit(mant
, fmt
->exponent
+fmt
->mantissa
); /* Last bit */
362 case 'i': /* __infinity__ */
367 ieee_flconvert(str
, mant
, &exponent
, error
);
368 if (mant
[0] & 0x8000) {
373 if (exponent
>= 2-expmax
&& exponent
<= expmax
) {
378 ieee_shr(mant
, fmt
->exponent
);
379 ieee_round(mant
, fmt
->words
);
380 /* did we scale up by one? */
381 if (mant
[0] & (implicit_one
<< 1)) {
386 mant
[0] &= (implicit_one
-1); /* remove leading one */
387 mant
[0] |= exponent
<< (15 - fmt
->exponent
);
388 } else if (exponent
< 2-expmax
&&
389 exponent
>= 2-expmax
-fmt
->mantissa
) {
393 int shift
= -(exponent
+ expmax
-2-fmt
->exponent
);
394 int sh
= shift
% 16, wds
= shift
/ 16;
396 if (ieee_round(mant
, fmt
->words
- wds
)
397 || (sh
> 0 && (mant
[0] & (0x8000 >> (sh
- 1))))) {
405 for (i
= fmt
->words
-1; i
>= wds
; i
--)
406 mant
[i
] = mant
[i
-wds
];
412 error(ERR_NONFATAL
, "overflow in floating-point constant");
415 memset(mant
, 0, 2*fmt
->words
);
420 memset(mant
, 0, 2*fmt
->words
);
426 for (mp
= &mant
[fmt
->words
], i
= 0; i
< fmt
->words
; i
++) {
432 return 1; /* success */
435 /* 80-bit format with 64-bit mantissa *including an explicit integer 1*
436 and 15-bit exponent. */
437 static int to_ldoub(char *str
, int32_t sign
, uint8_t *result
,
440 uint16_t mant
[MANT_WORDS
];
443 sign
= (sign
< 0 ? 0x8000L
: 0L);
446 uint16_t is_snan
= 0, is_qnan
= 0x8000;
462 put(result
+ 0, is_snan
);
465 put(result
+ 6, is_qnan
);
466 put(result
+ 8, 0x7fff|sign
);
470 ieee_flconvert(str
, mant
, &exponent
, error
);
471 if (mant
[0] & 0x8000) {
476 if (exponent
>= -16383 && exponent
<= 16384) {
481 if (ieee_round(mant
, 4)) /* did we scale up by one? */
482 ieee_shr(mant
, 1), mant
[0] |= 0x8000, exponent
++;
483 put(result
+ 0, mant
[3]);
484 put(result
+ 2, mant
[2]);
485 put(result
+ 4, mant
[1]);
486 put(result
+ 6, mant
[0]);
487 put(result
+ 8, exponent
| sign
);
488 } else if (exponent
< -16383 && exponent
>= -16446) {
492 int shift
= -(exponent
+ 16383);
493 int sh
= shift
% 16, wds
= shift
/ 16;
495 if (ieee_round(mant
, 4 - wds
)
496 || (sh
> 0 && (mant
[0] & (0x8000 >> (sh
- 1))))) {
502 put(result
+ 0, (wds
<= 3 ? mant
[3 - wds
] : 0));
503 put(result
+ 2, (wds
<= 2 ? mant
[2 - wds
] : 0));
504 put(result
+ 4, (wds
<= 1 ? mant
[1 - wds
] : 0));
505 put(result
+ 6, (wds
== 0 ? mant
[0] : 0));
506 put(result
+ 8, sign
);
509 error(ERR_NONFATAL
, "overflow in floating-point constant");
524 put(result
+ 8, sign
);
529 int float_const(char *number
, int32_t sign
, uint8_t *result
, int bytes
,
534 return to_float(number
, sign
, result
, &ieee_16
, error
);
536 return to_float(number
, sign
, result
, &ieee_32
, error
);
538 return to_float(number
, sign
, result
, &ieee_64
, error
);
540 return to_ldoub(number
, sign
, result
, error
);
542 return to_float(number
, sign
, result
, &ieee_128
, error
);
544 error(ERR_PANIC
, "strange value %d passed to float_const", bytes
);