Ciao,
Il Mar, 6 Marzo 2018 11:58 pm, Marc Glisse ha scritto:
> It is also possible to determine it at runtime (comparing 2^n+{0,1,2} for
> instance should work whether we round up, down or to the nearest;
> tests/misc.c compares 2^n+1-2^n to 1, which also looks good), doing it
Well, stealing the tests_dbl_mant_bits code from tests-misc.c, I wrote the
following code for mini-gmp. It seams to work...
Comments?
-----8<----------8<----
static unsigned int
gmp_tests_dbl_mant_bits (void)
{
unsigned int n;
volatile double x, y, d;
n = 1;
x = d = 1.0;
do {
x *= 2.0;
y = x + d;
d = y - x;
} while (d == 1.0 && ++n < sizeof(double) * CHAR_BIT);
return n;
}
double
mpz_get_d (const mpz_t u)
{
static unsigned int c = 0;
int m;
mp_limb_t l;
mp_size_t un;
double x;
double B = 2.0 * (double) GMP_LIMB_HIGHBIT;
un = GMP_ABS (u->_mp_size);
if (un == 0)
return 0.0;
if (c == 0)
c = gmp_tests_dbl_mant_bits ();
l = u->_mp_d[--un];
gmp_clz (m, l);
m = m + c - GMP_LIMB_BITS;
if (m < 0)
l &= GMP_LIMB_MAX << -m;
for (x = l; --un >= 0;)
{
x = B*x;
if (m > 0) {
l = u->_mp_d[un];
m -= GMP_LIMB_BITS;
if (m < 0)
l &= GMP_LIMB_MAX << -m;
x += l;
}
}
if (u->_mp_size < 0)
x = -x;
return x;
}
-----8<----------8<----
Ĝis,
m
--
http://bodrato.it/
_______________________________________________
gmp-devel mailing list
gmp-devel@gmplib.org
https://gmplib.org/mailman/listinfo/gmp-devel
