Generic get_d_2exp failures
Marc Glisse
marc.glisse at inria.fr
Sat Oct 28 15:55:43 UTC 2017
Hello,
the generic path in mpn_get_d sometimes rounds the result away from 0,
while mpz_get_d_2exp relies on a rounding towards 0 to deduce the exponent
from MPN_SIZEINBASE_2EXP. Various ways forward I can think of:
1) work around in mpz_get_d_2exp: if mpn_get_d returns 1 or -1, divide
that number by 2 and increase the exponent by 1. That doesn't give the
same rounding as with IEEE, but at least the output is in the
documented range [.5,1). (mpf has a similar issue IIRC)
2) The rounding occurs in the addition in
weight = 1/MP_BASE_AS_DOUBLE;
d = up[size - 1];
for (i = size - 2; i >= 0; i--)
{
d += up[i] * weight;
weight /= MP_BASE_AS_DOUBLE;
if (weight == 0)
break;
}
we could make that code uglier to make sure there is no rounding up,
maybe comparing d+up[i]*weight to d+oldweight, and trying again by zeroing
out an increasing number of low bits of up[i].
3) Drop the generic path. It hasn't passed the testsuite for a long time,
can't be that important.
Side note: I have seen code that uses MPFR for the sole purpose of
converting GMP types to double (with a controlled rounding).
--
Marc Glisse
More information about the gmp-devel
mailing list