[Bug-gmp] Bug in mpf_get_str.

[wilco oelen]

Hello,

I'm using gmp-3.0.1 for high precision floating point arithmatic. I think I have
found a bug, regarding the printing of floating point numbers in radix-10 format.

The system properties are as follows:

$ uname -a
FreeBSD ser2.home 3.3-RELEASE FreeBSD 3.3-RELEASE #15: Wed Jan 26
19:27:38 CET 2000     address@hidden:/usr/src/sys/compile/HOME  i386
$ gcc -v
gcc version 2.7.2.3
$ ./config.guess
k63-pc-freebsd3.3

No special options were used for building GMP 3.0.1. I simply did
./configure; make. The library is built in the directory .libs and the
library libgmp.a is linked with my application and the file gmp.h from
the toplevel gmp-3.0.1 directory is included in my C-sources. The
command make check reports no errors.
 

The bug is fairly simple to describe. When the number 1 is printed,
then it is printed as 0.9999999999..., even if the accuracy is much
higher and the number is given precisely as 1. The same holds for
printing 10, 100, 1000 and also for 0.1, 0.01 etc. No problems appear
for 2, 20, 200, 0.2, 0.02 etc.
The problem only appears, when printing with 35-digit accuracy.
Printing with 34 or 36 digit accuracy gives the correct result.
The problem is not affected by the internal accuracy with which
the floating point arithmatic is carried out.

The following C-program demonstrates the problem:
 

#include <stdio.h>
#include <stdlib.h>
 
#include "gmp.h"
 
int main(void)
{
  mp_exp_t e;
  char mant[1000];
 
  mpf_t x, y;
 
  /* An accuracy is taken of 1024 bits, appr. 308 digits. */
  /* The bug is not affected by this accuracy.            */
  mpf_set_default_prec(1024);
 
 
  mpf_init(x);
  mpf_init(y);
 
   /* Set x to exactly 1. mpf_set_ui(x, 1) has the same effect. */
  mpf_set_str(x, "1", 10);
 
  /* Print the number 1 with 35 digit accuracy.        */
  /* We expect '1', but we get 0.999...999 (35 nines). */
  /* When we print at 34-digit accuracy or 36-digit    */
  /* accuracy, then all works (we get '1').            */
  mpf_get_str(mant, &e, 10, 35, x);
  printf("WRONG: 0.%s * 10^%d\n", mant, (int)e);
 
  /* The problem with 35-digit accuracy does not occur */
  /* for 34 or 36 digit accuracy, not with 2, 3 etc.   */
  /* It occurs, however, for 10, 100, etc. and 0.1,    */
  /* 0.01, 0.001 etc. only for 35-digit accuracy.      */
 
  /* Print the number 1 with full precision. */
  mpf_get_str(mant, &e, 10, 0, x);
  printf("FULL1: 0.%s * 10^%d\n", mant, (int)e);  /* This works OK. */
 

  mpf_set_str(y, "1e-300", 10);
  mpf_sub(x, x, y);   /* Now we have 1 - 10^(-300) in x. */
  /* Print the number 1 - 10^(-300) with 35 digit accuracy. */
  /* Although this is less than 1, it is printed as 1 and   */
  /* this is right, because it is rounded.                  */
  /* This number is printed as 1, so we certainly expect    */
  /* the exact one to be printed as 1 and not as 0.999...   */
  mpf_get_str(mant, &e, 10, 35, x);
  printf("RIGHT: 0.%s * 10^%d\n", mant, (int)e);
 
  /* Print the number 1 - 10^(-300) with full precision. */
  /* This works OK. We get 0.999...999 (300 nines).      */
  mpf_get_str(mant, &e, 10, 0, x);
  printf("FULL1-1E-300: 0.%s * 10^%d\n", mant, (int)e);
 
  return 0;
}
 

I tested this program on a Pentium MMX 133 MHz with Redhat Linux and on that machine I got the same error,
with GMP 3.0.1 compiled for that platform. I think the bug is platform/CPU independent.

With regards,

Wilco Oelen