ERROR (check_random test 1): bad mpq_set_d results
Paul Zimmermann
Paul.Zimmermann at loria.fr
Tue Nov 25 09:14:33 CET 2008
Dennis,
> ERROR (check_random test 1): bad mpq_set_d results
> 3.098073531878046e-120
> 3.098073531878046e-120
> FAIL: t-get_d
>
> Is there anything I can look into for you ?
that test does the following:
1 - it generates a random double-precision number d
2 - it converts d to a mpq_t q using mpq_set_d
3 - it converts the numerator and denominator of q to double-precision numbers
nd and dd respectively, using mpz_get_d
4 - it divides nd by dd, giving the double-precision number d2
5 - it checks d equals d2
By construction, mpq_set_d produces a denominator of the form 2^e, thus we
should have dd = 2^e, and the division nd/dd in 4 should be exact.
Thus the problem might be due to the fact that mpz_get_d does not yield the
same value (up to scaling by 2^e) than the original d.
I suggest you add the lines marked /*+*/ in tests/mpq/t-get_d.c, recompile,
and run the test again (function check_random).
for (test = 0; test < reps; test++)
{
mpn_random2 (rp, LIMBS_PER_DOUBLE + 1);
d = 0.0;
for (i = LIMBS_PER_DOUBLE - 1; i >= 0; i--)
d = d * MP_BASE_AS_DOUBLE + rp[i];
d = my_ldexp (d, (int) (rp[LIMBS_PER_DOUBLE] % 1000) - 500);
mpq_set_d (q, d);
nd = mpz_get_d (mpq_numref (q));
dd = mpz_get_d (mpq_denref (q));
d2 = nd / dd;
if (d != d2)
{
printf ("ERROR (check_random test %d): bad mpq_set_d results\n", test);
printf ("%.16g\n", d);
printf ("%.16g\n", d2);
printf ("d=%1.16e\n", d); /*+*/
gmp_printf ("num(q)=%Zd\n", mpq_numref (q)); /*+*/
printf ("nd=%1.16e\n", nd); /*+*/
abort ();
}
}
Note that printing 16 digits of d and d2 is not enough, since 17 digits are
needed to distinguish two binary64 numbers in the worst case. Here both
9007199254736896/2^450 = 3.0980735318780457775e-120 and
9007199254736897/2^450 = 3.0980735318780461215e-120 round to
3.098073531878046e-120 to 16 digits. I guess d is one of those, and d2 the
other one.
Paul Zimmermann
More information about the gmp-bugs
mailing list