/* ** Implementation of encode/decode operations. */ #include #include #if defined(__arm) || defined(__hppa__) #define FLT_MINEXP -128 #define FLT_DIGITS 24 #endif #if defined(sequent) || defined(sun) || defined(mips) || defined(hp300) || defined(_IBMR2) || defined(linux) || defined(__386BSD__) || defined(bsdi) || defined(__alpha) || defined(__CYGWIN32__) || defined(__MINGW32__) || defined(__OpenBSD__) #define IEEE #endif #include "gmp.h" #include "gmp-impl.h" #include "flt.h" #ifdef HIGH_BYTE_FIRST #define LL 1 #define HH 0 #else #define LL 0 #define HH 1 #endif #define BASE 4294967296.0 #if defined(__alpha) #define DNBIGIT 1 #else #define DNBIGIT 2 #endif #define SNBIGIT 1 typedef union { float f; int i; } conv; /* ** Encode a mantissa (an Integer) and an exponent (an Int) as a floating point number. ** f = mant * 2^expo ** This operation is not as efficient as it could be, but it's portable. */ float sencode(MP_INT *man, int iexp) { double r; double base = BASE; int i; int size = (int)CONINFO_LARGESIZEU(man->sizeTag); /* Convert bignum to a double */ for(r = 0.0, i = abs(size)-1; i >= 0; i--) r = r * base + man->d[i]; /* Load new exponent */ /*printf("dencode ldexp %g %d\n", r, iexp);*/ if (r != 0.0) /* bug in mips ldexp */ r = ldexp(r, iexp); if (CONINFO_LARGEEXTRA(man->sizeTag)) r = -r; return (float)r; } #undef MINEXP #undef HIGHBIT #undef MSBBIT /* ** Decode a floating point number into a mantissa an an exponent. ** This operation is tricky to do in a portable and efficient way! ** Current implementation is really bad! man = amp_alloc(SNBIGIT); */ #define MINEXP (FLT_MINEXP - FLT_DIGITS - 1) #define HIGHBIT 0x00800000 #define MSBBIT 0x80000000 int sdecode(float f,MP_INT *man) { /* Do some bit fiddling on IEEE */ int high, iexp, sign; int size; conv c; c.f = f; high = c.i; if ((high & ~MSBBIT) == 0) { man->sizeTag = CONSTRW(0,0); man->d[0] = 0; /* !!! not necessary but better debuging */ size = 0; iexp = 0; } else { size = SNBIGIT; iexp = ((high >> 23) & 0xff) + MINEXP; sign = high; /* fprintf(stderr, "decode %g %08x %d\n", f, high, iexp); */ high &= HIGHBIT-1; if (iexp != MINEXP) /* don't add hidden bit to denorms */ high |= HIGHBIT; else { iexp++; /* A denorm, normalize the mantissa */ while (! (high & HIGHBIT)) { high <<= 1; iexp--; } } man->d[0] = high; if (sign < 0) man->sizeTag = CONSTRW((Int)size,1); else man->sizeTag = CONSTRW((Int)size,0); } return iexp; }