X-Git-Url: http://drtracing.org/?a=blobdiff_plain;f=sim%2Fcommon%2Fsim-fpu.c;h=74f5fd488c3409c5546faffe28b2d94bac4aa460;hb=611fa2f46c4ea0ca2f4153b6cf7ec31b8486eb28;hp=184f802492e0d9439e725b7ae69e20eacfa0d6ed;hpb=bdfe5c0439d0fdd6db272c2527841bf5640ce5be;p=deliverable%2Fbinutils-gdb.git diff --git a/sim/common/sim-fpu.c b/sim/common/sim-fpu.c index 184f802492..74f5fd488c 100644 --- a/sim/common/sim-fpu.c +++ b/sim/common/sim-fpu.c @@ -2,29 +2,20 @@ of the floating point routines in libgcc1.c for targets without hardware floating point. */ -/* Copyright (C) 1994,1997 Free Software Foundation, Inc. - -This file is free software; you can redistribute it and/or modify it -under the terms of the GNU General Public License as published by the -Free Software Foundation; either version 2, or (at your option) any -later version. - -In addition to the permissions in the GNU General Public License, the -Free Software Foundation gives you unlimited permission to link the -compiled version of this file with other programs, and to distribute -those programs without any restriction coming from the use of this -file. (The General Public License restrictions do apply in other -respects; for example, they cover modification of the file, and -distribution when not linked into another program.) - -This file is distributed in the hope that it will be useful, but -WITHOUT ANY WARRANTY; without even the implied warranty of -MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -General Public License for more details. +/* Copyright 1994-2019 Free Software Foundation, Inc. + +This program is free software; you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation; either version 3 of the License, or +(at your option) any later version. + +This program is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. You should have received a copy of the GNU General Public License -along with this program; see the file COPYING. If not, write to -the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA. */ +along with this program. If not, see . */ /* As a special exception, if you link this library with other files, some of which are compiled with GCC, to produce an executable, @@ -50,38 +41,48 @@ the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA. */ #include "sim-io.h" #include "sim-assert.h" +#ifdef HAVE_STDLIB_H +#include +#endif -/* Debugging support. */ +/* Debugging support. + If digits is -1, then print all digits. */ static void print_bits (unsigned64 x, int msbit, + int digits, sim_fpu_print_func print, void *arg) { unsigned64 bit = LSBIT64 (msbit); int i = 4; - while (bit) + while (bit && digits) { if (i == 0) print (arg, ","); + if ((x & bit)) print (arg, "1"); else print (arg, "0"); bit >>= 1; + + if (digits > 0) + digits--; i = (i + 1) % 4; } } -/* Quick and dirty conversion between a host double and host 64bit int */ +/* Quick and dirty conversion between a host double and host 64bit int. */ -typedef union { +typedef union +{ double d; unsigned64 i; -} sim_fpu_map; +} sim_fpu_map; /* A packed IEEE floating point number. @@ -97,7 +98,7 @@ typedef union { Zero (0 == BIASEDEXP && FRAC == 0): (sign ? "-" : "+") 0.0 - + Infinity (BIASEDEXP == EXPMAX && FRAC == 0): (sign ? "-" : "+") "infinity" @@ -110,11 +111,16 @@ typedef union { */ #define NR_EXPBITS (is_double ? 11 : 8) -#define NR_FRACBITS (is_double ? 52 : 23) +#define NR_FRACBITS (is_double ? 52 : 23) #define SIGNBIT (is_double ? MSBIT64 (0) : MSBIT64 (32)) -#define EXPMAX ((unsigned) (is_double ? 2047 : 255)) -#define EXPBIAS (is_double ? 1023 : 127) +#define EXPMAX32 (255) +#define EXMPAX64 (2047) +#define EXPMAX ((unsigned) (is_double ? EXMPAX64 : EXPMAX32)) + +#define EXPBIAS32 (127) +#define EXPBIAS64 (1023) +#define EXPBIAS (is_double ? EXPBIAS64 : EXPBIAS32) #define QUIET_NAN LSBIT64 (NR_FRACBITS - 1) @@ -128,10 +134,16 @@ typedef union { 64 bit - 32 bit - */ -#define NR_PAD (is_double ? 0 : 30) -#define PADMASK (is_double ? 0 : LSMASK64 (29, 0)) -#define NR_GUARDS ((is_double ? 8 : 7 ) + NR_PAD) +#define NR_PAD32 (30) +#define NR_PAD64 (0) +#define NR_PAD (is_double ? NR_PAD64 : NR_PAD32) +#define PADMASK (is_double ? 0 : LSMASK64 (NR_PAD32 - 1, 0)) + +#define NR_GUARDS32 (7 + NR_PAD32) +#define NR_GUARDS64 (8 + NR_PAD64) +#define NR_GUARDS (is_double ? NR_GUARDS64 : NR_GUARDS32) #define GUARDMASK LSMASK64 (NR_GUARDS - 1, 0) + #define GUARDMSB LSBIT64 (NR_GUARDS - 1) #define GUARDLSB LSBIT64 (NR_PAD) #define GUARDROUND LSMASK64 (NR_GUARDS - 2, 0) @@ -145,6 +157,9 @@ typedef union { #define FRAC32MASK LSMASK64 (63, NR_FRAC_GUARD - 32 + 1) #define NORMAL_EXPMIN (-(EXPBIAS)+1) + +#define NORMAL_EXPMAX32 (EXPBIAS32) +#define NORMAL_EXPMAX64 (EXPBIAS64) #define NORMAL_EXPMAX (EXPBIAS) @@ -163,7 +178,7 @@ typedef union { #define MAX_UINT (is_64bit ? MAX_UINT64 : MAX_UINT32) #define NR_INTBITS (is_64bit ? 64 : 32) -/* Squeese an unpacked sim_fpu struct into a 32/64 bit integer */ +/* Squeeze an unpacked sim_fpu struct into a 32/64 bit integer. */ STATIC_INLINE_SIM_FPU (unsigned64) pack_fpu (const sim_fpu *src, int is_double) @@ -175,22 +190,30 @@ pack_fpu (const sim_fpu *src, switch (src->class) { - /* create a NaN */ + /* Create a NaN. */ case sim_fpu_class_qnan: sign = src->sign; exp = EXPMAX; - /* force fraction to correct class */ + /* Force fraction to correct class. */ fraction = src->fraction; fraction >>= NR_GUARDS; +#ifdef SIM_QUIET_NAN_NEGATED + fraction |= QUIET_NAN - 1; +#else fraction |= QUIET_NAN; +#endif break; case sim_fpu_class_snan: sign = src->sign; exp = EXPMAX; - /* force fraction to correct class */ + /* Force fraction to correct class. */ fraction = src->fraction; fraction >>= NR_GUARDS; +#ifdef SIM_QUIET_NAN_NEGATED + fraction |= QUIET_NAN; +#else fraction &= ~QUIET_NAN; +#endif break; case sim_fpu_class_infinity: sign = src->sign; @@ -215,7 +238,7 @@ pack_fpu (const sim_fpu *src, int nr_shift = NORMAL_EXPMIN - src->normal_exp; if (nr_shift > NR_FRACBITS) { - /* underflow, just make the number zero */ + /* Underflow, just make the number zero. */ sign = src->sign; exp = 0; fraction = 0; @@ -224,7 +247,7 @@ pack_fpu (const sim_fpu *src, { sign = src->sign; exp = 0; - /* Shift by the value */ + /* Shift by the value. */ fraction = src->fraction; fraction >>= NR_GUARDS; fraction >>= nr_shift; @@ -235,7 +258,7 @@ pack_fpu (const sim_fpu *src, /* Infinity */ sign = src->sign; exp = EXPMAX; - fraction = 0; + fraction = 0; } else { @@ -243,7 +266,7 @@ pack_fpu (const sim_fpu *src, sign = src->sign; fraction = src->fraction; /* FIXME: Need to round according to WITH_SIM_FPU_ROUNDING - or some such */ + or some such. */ /* Round to nearest: If the guard bits are the all zero, but the first, then we're half way between two numbers, choose the one which makes the lsb of the answer 0. */ @@ -254,17 +277,17 @@ pack_fpu (const sim_fpu *src, } else { - /* Add a one to the guards to force round to nearest */ + /* Add a one to the guards to force round to nearest. */ fraction += GUARDROUND; } - if ((fraction & IMPLICIT_2)) /* rounding resulted in carry */ + if ((fraction & IMPLICIT_2)) /* Rounding resulted in carry. */ { exp += 1; fraction >>= 1; } fraction >>= NR_GUARDS; /* When exp == EXPMAX (overflow from carry) fraction must - have been made zero */ + have been made zero. */ ASSERT ((exp == EXPMAX) <= ((fraction & ~IMPLICIT_1) == 0)); } break; @@ -276,7 +299,7 @@ pack_fpu (const sim_fpu *src, | (exp << NR_FRACBITS) | LSMASKED64 (fraction, NR_FRACBITS - 1, 0)); - /* trace operation */ + /* Trace operation. */ #if 0 if (is_double) { @@ -290,12 +313,12 @@ pack_fpu (const sim_fpu *src, (long) LSEXTRACTED32 (packed, 23 - 1, 0)); } #endif - + return packed; } -/* Unpack a 32/64 bit integer into a sim_fpu structure */ +/* Unpack a 32/64 bit integer into a sim_fpu structure. */ STATIC_INLINE_SIM_FPU (void) unpack_fpu (sim_fpu *dst, unsigned64 packed, int is_double) { @@ -308,9 +331,10 @@ unpack_fpu (sim_fpu *dst, unsigned64 packed, int is_double) /* Hmm. Looks like 0 */ if (fraction == 0) { - /* tastes like zero */ + /* Tastes like zero. */ dst->class = sim_fpu_class_zero; dst->sign = sign; + dst->normal_exp = 0; } else { @@ -334,7 +358,7 @@ unpack_fpu (sim_fpu *dst, unsigned64 packed, int is_double) /* Huge exponent*/ if (fraction == 0) { - /* Attached to a zero fraction - means infinity */ + /* Attached to a zero fraction - means infinity. */ dst->class = sim_fpu_class_infinity; dst->sign = sign; /* dst->normal_exp = EXPBIAS; */ @@ -342,10 +366,17 @@ unpack_fpu (sim_fpu *dst, unsigned64 packed, int is_double) } else { - /* Non zero fraction, means NaN */ + int qnan; + + /* Non zero fraction, means NaN. */ dst->sign = sign; dst->fraction = (fraction << NR_GUARDS); - if (fraction >= QUIET_NAN) +#ifdef SIM_QUIET_NAN_NEGATED + qnan = (fraction & QUIET_NAN) == 0; +#else + qnan = fraction >= QUIET_NAN; +#endif + if (qnan) dst->class = sim_fpu_class_qnan; else dst->class = sim_fpu_class_snan; @@ -353,14 +384,14 @@ unpack_fpu (sim_fpu *dst, unsigned64 packed, int is_double) } else { - /* Nothing strange about this number */ + /* Nothing strange about this number. */ dst->class = sim_fpu_class_number; dst->sign = sign; dst->fraction = ((fraction << NR_GUARDS) | IMPLICIT_1); dst->normal_exp = exp - EXPBIAS; } - /* trace operation */ + /* Trace operation. */ #if 0 if (is_double) { @@ -392,7 +423,7 @@ unpack_fpu (sim_fpu *dst, unsigned64 packed, int is_double) } -/* Convert a floating point into an integer */ +/* Convert a floating point into an integer. */ STATIC_INLINE_SIM_FPU (int) fpu2i (signed64 *i, const sim_fpu *s, @@ -417,13 +448,13 @@ fpu2i (signed64 *i, *i = MIN_INT; /* FIXME */ return sim_fpu_status_invalid_cvi; } - /* map infinity onto MAX_INT... */ + /* Map infinity onto MAX_INT... */ if (sim_fpu_is_infinity (s)) { *i = s->sign ? MIN_INT : MAX_INT; return sim_fpu_status_invalid_cvi; } - /* it is a number, but a small one */ + /* It is a number, but a small one. */ if (s->normal_exp < 0) { *i = 0; @@ -438,7 +469,7 @@ fpu2i (signed64 *i, return 0; /* exact */ if (is_64bit) /* can't round */ return sim_fpu_status_invalid_cvi; /* must be overflow */ - /* For a 32bit with MAX_INT, rounding is possible */ + /* For a 32bit with MAX_INT, rounding is possible. */ switch (round) { case sim_fpu_round_default: @@ -474,7 +505,7 @@ fpu2i (signed64 *i, *i = s->sign ? MIN_INT : MAX_INT; return sim_fpu_status_invalid_cvi; } - /* normal number shift it into place */ + /* Normal number, shift it into place. */ tmp = s->fraction; shift = (s->normal_exp - (NR_FRAC_GUARD)); if (shift > 0) @@ -492,7 +523,7 @@ fpu2i (signed64 *i, return status; } -/* convert an integer into a floating point */ +/* Convert an integer into a floating point. */ STATIC_INLINE_SIM_FPU (int) i2fpu (sim_fpu *f, signed64 i, int is_64bit) { @@ -501,6 +532,7 @@ i2fpu (sim_fpu *f, signed64 i, int is_64bit) { f->class = sim_fpu_class_zero; f->sign = 0; + f->normal_exp = 0; } else { @@ -508,10 +540,10 @@ i2fpu (sim_fpu *f, signed64 i, int is_64bit) f->sign = (i < 0); f->normal_exp = NR_FRAC_GUARD; - if (f->sign) + if (f->sign) { /* Special case for minint, since there is no corresponding - +ve integer representation for it */ + +ve integer representation for it. */ if (i == MIN_INT) { f->fraction = IMPLICIT_1; @@ -525,9 +557,9 @@ i2fpu (sim_fpu *f, signed64 i, int is_64bit) if (f->fraction >= IMPLICIT_2) { - do + do { - f->fraction >>= 1; + f->fraction = (f->fraction >> 1) | (f->fraction & 1); f->normal_exp += 1; } while (f->fraction >= IMPLICIT_2); @@ -564,7 +596,7 @@ i2fpu (sim_fpu *f, signed64 i, int is_64bit) } -/* Convert a floating point into an integer */ +/* Convert a floating point into an integer. */ STATIC_INLINE_SIM_FPU (int) fpu2u (unsigned64 *u, const sim_fpu *s, int is_64bit) { @@ -581,19 +613,19 @@ fpu2u (unsigned64 *u, const sim_fpu *s, int is_64bit) *u = 0; return 0; } - /* it is a negative number */ + /* It is a negative number. */ if (s->sign) { *u = 0; return 0; } - /* get reasonable MAX_USI_INT... */ + /* Get reasonable MAX_USI_INT... */ if (sim_fpu_is_infinity (s)) { *u = MAX_UINT; return 0; } - /* it is a number, but a small one */ + /* It is a number, but a small one. */ if (s->normal_exp < 0) { *u = 0; @@ -621,7 +653,7 @@ fpu2u (unsigned64 *u, const sim_fpu *s, int is_64bit) return 0; } -/* Convert an unsigned integer into a floating point */ +/* Convert an unsigned integer into a floating point. */ STATIC_INLINE_SIM_FPU (int) u2fpu (sim_fpu *f, unsigned64 u, int is_64bit) { @@ -629,6 +661,7 @@ u2fpu (sim_fpu *f, unsigned64 u, int is_64bit) { f->class = sim_fpu_class_zero; f->sign = 0; + f->normal_exp = 0; } else { @@ -698,6 +731,40 @@ sim_fpu_to64 (unsigned64 *u, } +INLINE_SIM_FPU (void) +sim_fpu_fractionto (sim_fpu *f, + int sign, + int normal_exp, + unsigned64 fraction, + int precision) +{ + int shift = (NR_FRAC_GUARD - precision); + f->class = sim_fpu_class_number; + f->sign = sign; + f->normal_exp = normal_exp; + /* Shift the fraction to where sim-fpu expects it. */ + if (shift >= 0) + f->fraction = (fraction << shift); + else + f->fraction = (fraction >> -shift); + f->fraction |= IMPLICIT_1; +} + + +INLINE_SIM_FPU (unsigned64) +sim_fpu_tofraction (const sim_fpu *d, + int precision) +{ + /* We have NR_FRAC_GUARD bits, we want only PRECISION bits. */ + int shift = (NR_FRAC_GUARD - precision); + unsigned64 fraction = (d->fraction & ~IMPLICIT_1); + if (shift >= 0) + return fraction >> shift; + else + return fraction << -shift; +} + + /* Rounding */ STATIC_INLINE_SIM_FPU (int) @@ -760,7 +827,7 @@ do_normal_underflow (sim_fpu *f, /* Round a number using NR_GUARDS. - Will return the rounded number or F->FRACTION == 0 when underflow */ + Will return the rounded number or F->FRACTION == 0 when underflow. */ STATIC_INLINE_SIM_FPU (int) do_normal_round (sim_fpu *f, @@ -802,7 +869,7 @@ do_normal_round (sim_fpu *f, break; } f->fraction &= ~guardmask; - /* round if needed, handle resulting overflow */ + /* Round if needed, handle resulting overflow. */ if ((status & sim_fpu_status_rounded)) { f->fraction += fraclsb; @@ -833,7 +900,7 @@ do_round (sim_fpu *f, return 0; break; case sim_fpu_class_snan: - /* Quieten a SignalingNaN */ + /* Quieten a SignalingNaN. */ f->class = sim_fpu_class_qnan; return sim_fpu_status_invalid_snan; break; @@ -855,7 +922,7 @@ do_round (sim_fpu *f, && !(denorm & sim_fpu_denorm_zero)) { status = do_normal_round (f, shift + NR_GUARDS, round); - if (f->fraction == 0) /* rounding underflowed */ + if (f->fraction == 0) /* Rounding underflowed. */ { status |= do_normal_underflow (f, is_double, round); } @@ -867,7 +934,7 @@ do_round (sim_fpu *f, before rounding, some after! */ if (status & sim_fpu_status_inexact) status |= sim_fpu_status_underflow; - /* Flag that resultant value has been denormalized */ + /* Flag that resultant value has been denormalized. */ f->class = sim_fpu_class_denorm; } else if ((denorm & sim_fpu_denorm_underflow_inexact)) @@ -893,7 +960,7 @@ do_round (sim_fpu *f, /* f->class = sim_fpu_class_zero; */ status |= do_normal_underflow (f, is_double, round); else if (f->normal_exp > NORMAL_EXPMAX) - /* oops! rounding caused overflow */ + /* Oops! rounding caused overflow. */ status |= do_normal_overflow (f, is_double, round); } ASSERT ((f->class == sim_fpu_class_number @@ -992,13 +1059,13 @@ sim_fpu_add (sim_fpu *f, /* use exp of larger */ if (shift >= NR_FRAC_GUARD) { - /* left has much bigger magnitute */ + /* left has much bigger magnitude */ *f = *l; return sim_fpu_status_inexact; } if (shift <= - NR_FRAC_GUARD) { - /* right has much bigger magnitute */ + /* right has much bigger magnitude */ *f = *r; return sim_fpu_status_inexact; } @@ -1010,7 +1077,7 @@ sim_fpu_add (sim_fpu *f, if (rfraction & LSMASK64 (shift - 1, 0)) { status |= sim_fpu_status_inexact; - rfraction |= LSBIT64 (shift); /* stick LSBit */ + rfraction |= LSBIT64 (shift); /* Stick LSBit. */ } rfraction >>= shift; } @@ -1020,7 +1087,7 @@ sim_fpu_add (sim_fpu *f, if (lfraction & LSMASK64 (- shift - 1, 0)) { status |= sim_fpu_status_inexact; - lfraction |= LSBIT64 (- shift); /* stick LSBit */ + lfraction |= LSBIT64 (- shift); /* Stick LSBit. */ } lfraction >>= -shift; } @@ -1029,7 +1096,7 @@ sim_fpu_add (sim_fpu *f, f->normal_exp = r->normal_exp; } - /* perform the addition */ + /* Perform the addition. */ if (l->sign) lfraction = - lfraction; if (r->sign) @@ -1045,7 +1112,7 @@ sim_fpu_add (sim_fpu *f, /* sign? */ f->class = sim_fpu_class_number; - if ((signed64) f->fraction >= 0) + if (((signed64) f->fraction) >= 0) f->sign = 0; else { @@ -1053,7 +1120,7 @@ sim_fpu_add (sim_fpu *f, f->fraction = - f->fraction; } - /* normalize it */ + /* Normalize it. */ if ((f->fraction & IMPLICIT_2)) { f->fraction = (f->fraction >> 1) | (f->fraction & 1); @@ -1145,13 +1212,13 @@ sim_fpu_sub (sim_fpu *f, /* use exp of larger */ if (shift >= NR_FRAC_GUARD) { - /* left has much bigger magnitute */ + /* left has much bigger magnitude */ *f = *l; return sim_fpu_status_inexact; } if (shift <= - NR_FRAC_GUARD) { - /* right has much bigger magnitute */ + /* right has much bigger magnitude */ *f = *r; f->sign = !r->sign; return sim_fpu_status_inexact; @@ -1164,7 +1231,7 @@ sim_fpu_sub (sim_fpu *f, if (rfraction & LSMASK64 (shift - 1, 0)) { status |= sim_fpu_status_inexact; - rfraction |= LSBIT64 (shift); /* stick LSBit */ + rfraction |= LSBIT64 (shift); /* Stick LSBit. */ } rfraction >>= shift; } @@ -1174,7 +1241,7 @@ sim_fpu_sub (sim_fpu *f, if (lfraction & LSMASK64 (- shift - 1, 0)) { status |= sim_fpu_status_inexact; - lfraction |= LSBIT64 (- shift); /* stick LSBit */ + lfraction |= LSBIT64 (- shift); /* Stick LSBit. */ } lfraction >>= -shift; } @@ -1183,7 +1250,7 @@ sim_fpu_sub (sim_fpu *f, f->normal_exp = r->normal_exp; } - /* perform the subtraction */ + /* Perform the subtraction. */ if (l->sign) lfraction = - lfraction; if (!r->sign) @@ -1199,7 +1266,7 @@ sim_fpu_sub (sim_fpu *f, /* sign? */ f->class = sim_fpu_class_number; - if ((signed64) f->fraction >= 0) + if (((signed64) f->fraction) >= 0) f->sign = 0; else { @@ -1207,7 +1274,7 @@ sim_fpu_sub (sim_fpu *f, f->fraction = - f->fraction; } - /* normalize it */ + /* Normalize it. */ if ((f->fraction & IMPLICIT_2)) { f->fraction = (f->fraction >> 1) | (f->fraction & 1); @@ -1284,7 +1351,7 @@ sim_fpu_mul (sim_fpu *f, return 0; } /* Calculate the mantissa by multiplying both 64bit numbers to get a - 128 bit number */ + 128 bit number. */ { unsigned64 low; unsigned64 high; @@ -1308,14 +1375,14 @@ sim_fpu_mul (sim_fpu *f, res2 += ((ps_hh__ >> 32) & 0xffffffff) + pp_hh; high = res2; low = res0; - + f->normal_exp = l->normal_exp + r->normal_exp; f->sign = l->sign ^ r->sign; f->class = sim_fpu_class_number; /* Input is bounded by [1,2) ; [2^60,2^61) Output is bounded by [1,4) ; [2^120,2^122) */ - + /* Adjust the exponent according to where the decimal point ended up in the high 64 bit word. In the source the decimal point was at NR_FRAC_GUARD. */ @@ -1327,15 +1394,7 @@ sim_fpu_mul (sim_fpu *f, ASSERT (high >= LSBIT64 ((NR_FRAC_GUARD * 2) - 64)); ASSERT (LSBIT64 (((NR_FRAC_GUARD + 1) * 2) - 64) < IMPLICIT_1); -#if 0 - printf ("\n"); - print_bits (high, 63, (sim_fpu_print_func*)fprintf, stdout); - printf (";"); - print_bits (low, 63, (sim_fpu_print_func*)fprintf, stdout); - printf ("\n"); -#endif - - /* normalize */ + /* Normalize. */ do { f->normal_exp--; @@ -1346,13 +1405,6 @@ sim_fpu_mul (sim_fpu *f, } while (high < IMPLICIT_1); -#if 0 - print_bits (high, 63, (sim_fpu_print_func*)fprintf, stdout); - printf (";"); - print_bits (low, 63, (sim_fpu_print_func*)fprintf, stdout); - printf ("\n"); -#endif - ASSERT (high >= IMPLICIT_1 && high < IMPLICIT_2); if (low != 0) { @@ -1439,7 +1491,7 @@ sim_fpu_div (sim_fpu *f, } /* Calculate the mantissa by multiplying both 64bit numbers to get a - 128 bit number */ + 128 bit number. */ { /* quotient = ( ( numerator / denominator) x 2^(numerator exponent - denominator exponent) @@ -1463,8 +1515,8 @@ sim_fpu_div (sim_fpu *f, f->normal_exp--; } ASSERT (numerator >= denominator); - - /* Gain extra precision, already used one spare bit */ + + /* Gain extra precision, already used one spare bit. */ numerator <<= NR_SPARE; denominator <<= NR_SPARE; @@ -1482,17 +1534,7 @@ sim_fpu_div (sim_fpu *f, numerator <<= 1; } -#if 0 - printf ("\n"); - print_bits (quotient, 63, (sim_fpu_print_func*)fprintf, stdout); - printf ("\n"); - print_bits (numerator, 63, (sim_fpu_print_func*)fprintf, stdout); - printf ("\n"); - print_bits (denominator, 63, (sim_fpu_print_func*)fprintf, stdout); - printf ("\n"); -#endif - - /* discard (but save) the extra bits */ + /* Discard (but save) the extra bits. */ if ((quotient & LSMASK64 (NR_SPARE -1, 0))) quotient = (quotient >> NR_SPARE) | 1; else @@ -1502,7 +1544,7 @@ sim_fpu_div (sim_fpu *f, ASSERT (f->fraction >= IMPLICIT_1 && f->fraction < IMPLICIT_2); if (numerator != 0) { - f->fraction |= 1; /* stick remaining bits */ + f->fraction |= 1; /* Stick remaining bits. */ return sim_fpu_status_inexact; } else @@ -1512,81 +1554,298 @@ sim_fpu_div (sim_fpu *f, INLINE_SIM_FPU (int) -sim_fpu_neg (sim_fpu *f, +sim_fpu_rem (sim_fpu *f, + const sim_fpu *l, const sim_fpu *r) { + if (sim_fpu_is_snan (l)) + { + *f = *l; + f->class = sim_fpu_class_qnan; + return sim_fpu_status_invalid_snan; + } if (sim_fpu_is_snan (r)) { *f = *r; f->class = sim_fpu_class_qnan; return sim_fpu_status_invalid_snan; } + if (sim_fpu_is_qnan (l)) + { + *f = *l; + f->class = sim_fpu_class_qnan; + return 0; + } if (sim_fpu_is_qnan (r)) { *f = *r; + f->class = sim_fpu_class_qnan; return 0; } - *f = *r; - f->sign = !r->sign; - return 0; + if (sim_fpu_is_infinity (l)) + { + *f = sim_fpu_qnan; + return sim_fpu_status_invalid_irx; + } + if (sim_fpu_is_zero (r)) + { + *f = sim_fpu_qnan; + return sim_fpu_status_invalid_div0; + } + if (sim_fpu_is_zero (l)) + { + *f = *l; + return 0; + } + if (sim_fpu_is_infinity (r)) + { + *f = *l; + return 0; + } + { + sim_fpu n, tmp; + + /* Remainder is calculated as l-n*r, where n is l/r rounded to the + nearest integer. The variable n is rounded half even. */ + + sim_fpu_div (&n, l, r); + sim_fpu_round_64 (&n, 0, 0); + + if (n.normal_exp < -1) /* If n looks like zero just return l. */ + { + *f = *l; + return 0; + } + else if (n.class == sim_fpu_class_number + && n.normal_exp <= (NR_FRAC_GUARD)) /* If not too large round. */ + do_normal_round (&n, (NR_FRAC_GUARD) - n.normal_exp, sim_fpu_round_near); + + /* Mark 0's as zero so multiply can detect zero. */ + if (n.fraction == 0) + n.class = sim_fpu_class_zero; + + /* Calculate n*r. */ + sim_fpu_mul (&tmp, &n, r); + sim_fpu_round_64 (&tmp, 0, 0); + + /* Finally calculate l-n*r. */ + sim_fpu_sub (f, l, &tmp); + + return 0; + } } INLINE_SIM_FPU (int) -sim_fpu_abs (sim_fpu *f, +sim_fpu_max (sim_fpu *f, + const sim_fpu *l, const sim_fpu *r) { + if (sim_fpu_is_snan (l)) + { + *f = *l; + f->class = sim_fpu_class_qnan; + return sim_fpu_status_invalid_snan; + } if (sim_fpu_is_snan (r)) { *f = *r; f->class = sim_fpu_class_qnan; return sim_fpu_status_invalid_snan; } + if (sim_fpu_is_qnan (l)) + { + *f = *l; + return 0; + } if (sim_fpu_is_qnan (r)) { *f = *r; return 0; } - *f = *r; - f->sign = 0; - return 0; + if (sim_fpu_is_infinity (l)) + { + if (sim_fpu_is_infinity (r) + && l->sign == r->sign) + { + *f = sim_fpu_qnan; + return sim_fpu_status_invalid_isi; + } + if (l->sign) + *f = *r; /* -inf < anything */ + else + *f = *l; /* +inf > anything */ + return 0; + } + if (sim_fpu_is_infinity (r)) + { + if (r->sign) + *f = *l; /* anything > -inf */ + else + *f = *r; /* anything < +inf */ + return 0; + } + if (l->sign > r->sign) + { + *f = *r; /* -ve < +ve */ + return 0; + } + if (l->sign < r->sign) + { + *f = *l; /* +ve > -ve */ + return 0; + } + ASSERT (l->sign == r->sign); + if (l->normal_exp > r->normal_exp + || (l->normal_exp == r->normal_exp + && l->fraction > r->fraction)) + { + /* |l| > |r| */ + if (l->sign) + *f = *r; /* -ve < -ve */ + else + *f = *l; /* +ve > +ve */ + return 0; + } + else + { + /* |l| <= |r| */ + if (l->sign) + *f = *l; /* -ve > -ve */ + else + *f = *r; /* +ve < +ve */ + return 0; + } } INLINE_SIM_FPU (int) -sim_fpu_inv (sim_fpu *f, +sim_fpu_min (sim_fpu *f, + const sim_fpu *l, const sim_fpu *r) { + if (sim_fpu_is_snan (l)) + { + *f = *l; + f->class = sim_fpu_class_qnan; + return sim_fpu_status_invalid_snan; + } if (sim_fpu_is_snan (r)) { *f = *r; f->class = sim_fpu_class_qnan; return sim_fpu_status_invalid_snan; } + if (sim_fpu_is_qnan (l)) + { + *f = *l; + return 0; + } if (sim_fpu_is_qnan (r)) { *f = *r; - f->class = sim_fpu_class_qnan; + return 0; + } + if (sim_fpu_is_infinity (l)) + { + if (sim_fpu_is_infinity (r) + && l->sign == r->sign) + { + *f = sim_fpu_qnan; + return sim_fpu_status_invalid_isi; + } + if (l->sign) + *f = *l; /* -inf < anything */ + else + *f = *r; /* +inf > anthing */ return 0; } if (sim_fpu_is_infinity (r)) { - *f = sim_fpu_zero; - f->sign = r->sign; + if (r->sign) + *f = *r; /* anything > -inf */ + else + *f = *l; /* anything < +inf */ return 0; } - if (sim_fpu_is_zero (r)) + if (l->sign > r->sign) { - f->class = sim_fpu_class_infinity; - f->sign = r->sign; - return sim_fpu_status_invalid_div0; + *f = *l; /* -ve < +ve */ + return 0; + } + if (l->sign < r->sign) + { + *f = *r; /* +ve > -ve */ + return 0; + } + ASSERT (l->sign == r->sign); + if (l->normal_exp > r->normal_exp + || (l->normal_exp == r->normal_exp + && l->fraction > r->fraction)) + { + /* |l| > |r| */ + if (l->sign) + *f = *l; /* -ve < -ve */ + else + *f = *r; /* +ve > +ve */ + return 0; + } + else + { + /* |l| <= |r| */ + if (l->sign) + *f = *r; /* -ve > -ve */ + else + *f = *l; /* +ve < +ve */ + return 0; + } +} + + +INLINE_SIM_FPU (int) +sim_fpu_neg (sim_fpu *f, + const sim_fpu *r) +{ + if (sim_fpu_is_snan (r)) + { + *f = *r; + f->class = sim_fpu_class_qnan; + return sim_fpu_status_invalid_snan; + } + if (sim_fpu_is_qnan (r)) + { + *f = *r; + return 0; } *f = *r; - f->normal_exp = - r->normal_exp; + f->sign = !r->sign; + return 0; +} + + +INLINE_SIM_FPU (int) +sim_fpu_abs (sim_fpu *f, + const sim_fpu *r) +{ + *f = *r; + f->sign = 0; + if (sim_fpu_is_snan (r)) + { + f->class = sim_fpu_class_qnan; + return sim_fpu_status_invalid_snan; + } return 0; } +INLINE_SIM_FPU (int) +sim_fpu_inv (sim_fpu *f, + const sim_fpu *r) +{ + return sim_fpu_div (f, &sim_fpu_one, r); +} + + INLINE_SIM_FPU (int) sim_fpu_sqrt (sim_fpu *f, const sim_fpu *r) @@ -1605,6 +1864,7 @@ sim_fpu_sqrt (sim_fpu *f, { f->class = sim_fpu_class_zero; f->sign = r->sign; + f->normal_exp = 0; return 0; } if (sim_fpu_is_infinity (r)) @@ -1635,20 +1895,20 @@ sim_fpu_sqrt (sim_fpu *f, * * Developed at SunPro, a Sun Microsystems, Inc. business. * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice + * software is freely granted, provided that this notice * is preserved. * ==================================================== */ - + /* __ieee754_sqrt(x) * Return correctly rounded sqrt. * ------------------------------------------ * | Use the hardware sqrt if you have one | * ------------------------------------------ - * Method: - * Bit by bit method using integer arithmetic. (Slow, but portable) + * Method: + * Bit by bit method using integer arithmetic. (Slow, but portable) * 1. Normalization - * Scale x to y in [1,4) with even powers of 2: + * Scale x to y in [1,4) with even powers of 2: * find an integer k such that 1 <= (y=x*2^(2k)) < 4, then * sqrt(x) = 2^k * sqrt(y) - @@ -1668,9 +1928,9 @@ sim_fpu_sqrt (sim_fpu *f, * i+1 2 * s = 2*q , and y = 2 * ( y - q ). (1) * i i i i - * - * To compute q from q , one checks whether - * i+1 i + * + * To compute q from q , one checks whether + * i+1 i * * -(i+1) 2 * (q + 2 ) <= y. (2) @@ -1679,13 +1939,13 @@ sim_fpu_sqrt (sim_fpu *f, * If (2) is false, then q = q ; otherwise q = q + 2 . * i+1 i i+1 i * - * With some algebric manipulation, it is not difficult to see - * that (2) is equivalent to + * With some algebraic manipulation, it is not difficult to see + * that (2) is equivalent to * -(i+1) * s + 2 <= y (3) * i i * - * The advantage of (3) is that s and y can be computed by + * The advantage of (3) is that s and y can be computed by * i i * the following recurrence formula: * if (3) is false @@ -1701,15 +1961,15 @@ sim_fpu_sqrt (sim_fpu *f, * -i -(i+1) * s = s + 2 , y = y - s - 2 (5) * i+1 i i+1 i i - * + * - - -(i+1) - - NOTE: y = 2 (y - s - 2 ) + - NOTE: y = 2 (y - s - 2 ) - i+1 i i - - * One may easily use induction to prove (4) and (5). + * One may easily use induction to prove (4) and (5). * Note. Since the left hand side of (3) contain only i+2 bits, - * it does not necessary to do a full (53-bit) comparison + * it does not necessary to do a full (53-bit) comparison * in (3). * 3. Final rounding * After generating the 53 bits result, we compute one more bit. @@ -1719,19 +1979,19 @@ sim_fpu_sqrt (sim_fpu *f, * The rounding mode can be detected by checking whether * huge + tiny is equal to huge, and whether huge - tiny is * equal to huge for some floating point number "huge" and "tiny". - * + * * Special cases: * sqrt(+-0) = +-0 ... exact * sqrt(inf) = inf * sqrt(-ve) = NaN ... with invalid signal - * sqrt(NaN) = NaN ... with invalid signal for signaling NaN + * sqrt(NaN) = NaN ... with invalid signal for signalling NaN * * Other methods : see the appended file at the end of the program below. *--------------- */ { - /* generate sqrt(x) bit by bit */ + /* Generate sqrt(x) bit by bit. */ unsigned64 y; unsigned64 q; unsigned64 s; @@ -1742,7 +2002,7 @@ sim_fpu_sqrt (sim_fpu *f, y = r->fraction; f->normal_exp = (r->normal_exp >> 1); /* exp = [exp/2] */ - /* odd exp, double x to make it even */ + /* Odd exp, double x to make it even. */ ASSERT (y >= IMPLICIT_1 && y < IMPLICIT_4); if ((r->normal_exp & 1)) { @@ -1754,7 +2014,7 @@ sim_fpu_sqrt (sim_fpu *f, b = IMPLICIT_1; q = 0; s = 0; - + while (b) { unsigned64 t = s + b; @@ -1772,7 +2032,7 @@ sim_fpu_sqrt (sim_fpu *f, f->fraction = q; if (y != 0) { - f->fraction |= 1; /* stick remaining bits */ + f->fraction |= 1; /* Stick remaining bits. */ return sim_fpu_status_inexact; } else @@ -1876,7 +2136,17 @@ INLINE_SIM_FPU (double) sim_fpu_2d (const sim_fpu *s) { sim_fpu_map val; - val.i = pack_fpu (s, 1); + if (sim_fpu_is_snan (s)) + { + /* gag SNaN's */ + sim_fpu n = *s; + n.class = sim_fpu_class_qnan; + val.i = pack_fpu (&n, 1); + } + else + { + val.i = pack_fpu (s, 1); + } return val.d; } @@ -1991,6 +2261,38 @@ sim_fpu_is_denorm (const sim_fpu *d) } } + +INLINE_SIM_FPU (int) +sim_fpu_sign (const sim_fpu *d) +{ + return d->sign; +} + + +INLINE_SIM_FPU (int) +sim_fpu_exp (const sim_fpu *d) +{ + return d->normal_exp; +} + + +INLINE_SIM_FPU (unsigned64) +sim_fpu_fraction (const sim_fpu *d) +{ + return d->fraction; +} + + +INLINE_SIM_FPU (unsigned64) +sim_fpu_guard (const sim_fpu *d, int is_double) +{ + unsigned64 rv; + unsigned64 guardmask = LSMASK64 (NR_GUARDS - 1, 0); + rv = (d->fraction & guardmask) >> NR_PAD; + return rv; +} + + INLINE_SIM_FPU (int) sim_fpu_is (const sim_fpu *d) { @@ -2001,8 +2303,10 @@ sim_fpu_is (const sim_fpu *d) case sim_fpu_class_snan: return SIM_FPU_IS_SNAN; case sim_fpu_class_infinity: - return SIM_FPU_IS_NINF; - return SIM_FPU_IS_PINF; + if (d->sign) + return SIM_FPU_IS_NINF; + else + return SIM_FPU_IS_PINF; case sim_fpu_class_number: if (d->sign) return SIM_FPU_IS_NNUMBER; @@ -2206,8 +2510,26 @@ sim_fpu_gt (int *is, /* A number of useful constants */ -const sim_fpu sim_fpu_zero = { sim_fpu_class_zero, }; -const sim_fpu sim_fpu_qnan = { sim_fpu_class_qnan, }; +#if EXTERN_SIM_FPU_P +const sim_fpu sim_fpu_zero = { + sim_fpu_class_zero, 0, 0, 0 +}; +const sim_fpu sim_fpu_qnan = { + sim_fpu_class_qnan, 0, 0, 0 +}; +const sim_fpu sim_fpu_one = { + sim_fpu_class_number, 0, IMPLICIT_1, 0 +}; +const sim_fpu sim_fpu_two = { + sim_fpu_class_number, 0, IMPLICIT_1, 1 +}; +const sim_fpu sim_fpu_max32 = { + sim_fpu_class_number, 0, LSMASK64 (NR_FRAC_GUARD, NR_GUARDS32), NORMAL_EXPMAX32 +}; +const sim_fpu sim_fpu_max64 = { + sim_fpu_class_number, 0, LSMASK64 (NR_FRAC_GUARD, NR_GUARDS64), NORMAL_EXPMAX64 +}; +#endif /* For debugging */ @@ -2216,18 +2538,27 @@ INLINE_SIM_FPU (void) sim_fpu_print_fpu (const sim_fpu *f, sim_fpu_print_func *print, void *arg) +{ + sim_fpu_printn_fpu (f, print, -1, arg); +} + +INLINE_SIM_FPU (void) +sim_fpu_printn_fpu (const sim_fpu *f, + sim_fpu_print_func *print, + int digits, + void *arg) { print (arg, "%s", f->sign ? "-" : "+"); switch (f->class) { case sim_fpu_class_qnan: print (arg, "0."); - print_bits (f->fraction, NR_FRAC_GUARD - 1, print, arg); + print_bits (f->fraction, NR_FRAC_GUARD - 1, digits, print, arg); print (arg, "*QuietNaN"); break; case sim_fpu_class_snan: print (arg, "0."); - print_bits (f->fraction, NR_FRAC_GUARD - 1, print, arg); + print_bits (f->fraction, NR_FRAC_GUARD - 1, digits, print, arg); print (arg, "*SignalNaN"); break; case sim_fpu_class_zero: @@ -2239,8 +2570,8 @@ sim_fpu_print_fpu (const sim_fpu *f, case sim_fpu_class_number: case sim_fpu_class_denorm: print (arg, "1."); - print_bits (f->fraction, NR_FRAC_GUARD - 1, print, arg); - print (arg, "*2^%+-5d", f->normal_exp); + print_bits (f->fraction, NR_FRAC_GUARD - 1, digits, print, arg); + print (arg, "*2^%+d", f->normal_exp); ASSERT (f->fraction >= IMPLICIT_1); ASSERT (f->fraction < IMPLICIT_2); } @@ -2253,7 +2584,7 @@ sim_fpu_print_status (int status, void *arg) { int i = 1; - char *prefix = ""; + const char *prefix = ""; while (status >= i) { switch ((sim_fpu_status) (status & i)) @@ -2288,27 +2619,24 @@ sim_fpu_print_status (int status, case sim_fpu_status_invalid_sqrt: print (arg, "%sSQRT", prefix); break; + case sim_fpu_status_invalid_irx: + print (arg, "%sIRX", prefix); break; case sim_fpu_status_inexact: print (arg, "%sX", prefix); break; - break; case sim_fpu_status_overflow: print (arg, "%sO", prefix); break; - break; case sim_fpu_status_underflow: print (arg, "%sU", prefix); break; - break; case sim_fpu_status_invalid_div0: print (arg, "%s/", prefix); break; - break; case sim_fpu_status_rounded: print (arg, "%sR", prefix); break; - break; } i <<= 1; prefix = ",";