Bug 33469 - Default formats for real input are not precise enough
Summary: Default formats for real input are not precise enough
Status: RESOLVED FIXED
Alias: None
Product: gcc
Classification: Unclassified
Component: libfortran (show other bugs)
Version: 4.3.0
: P3 normal
Target Milestone: 4.3.0
Assignee: Francois-Xavier Coudert
URL:
Keywords: wrong-code
Depends on:
Blocks:
 
Reported: 2007-09-18 09:45 UTC by Dominique d'Humieres
Modified: 2007-10-03 19:34 UTC (History)
2 users (show)

See Also:
Host:
Target:
Build:
Known to work:
Known to fail:
Last reconfirmed: 2007-09-19 09:55:15


Attachments
test case (435 bytes, application/octet-stream)
2007-09-18 09:47 UTC, Dominique d'Humieres
Details

Note You need to log in before you can comment on or make changes to this bug.
Description Dominique d'Humieres 2007-09-18 09:45:54 UTC
Following a discussion on IRC with FX Coudert, I think the number of digits printed by default in formatted output should be increased by one. Otherwise the logic of the test case large_real_kind_form_io_2.f90 is flawed as shown by the attached code. Its output on AMD64 is

 real(4)
 default         808
 1PG20.6        1881
 1PG20.7         808
 1PG20.8           0

 real(8)
 default         1778
 1PG30.14        1978
 1PG30.15        1778
 1PG30.16           0

 real(10)
 default          916
 1PG60.18        1892
 1PG60.19         916
 1PG60.20           0

where the integers are the number of failures to read back a number. g95 gives:

 real(4)
 default 0
 1PG20.6 1881
 1PG20.7 808
 1PG20.8 0

 real(8)
 default  0
 1PG30.14 1978
 1PG30.15 1778
 1PG30.16 0

 real(10)
 default  0
 1PG60.18 1892
 1PG60.19 916
 1PG60.20 0
Comment 1 Dominique d'Humieres 2007-09-18 09:47:19 UTC
Created attachment 14219 [details]
test case
Comment 2 Tobias Burnus 2007-09-18 17:56:57 UTC
A compiler comparison:

 real(4)  gfortran g95 NAGf95 ifort sunf95 openf95
 default   808       0    808   808     0       0
 1PG20.6  1881    1881   1881  1881  1881    1881
 1PG20.7   808     808    808   808   808     808
 1PG20.8     0       0      0     0     0       0

 real(8)
 default  1778       0      0  1778     0       0
 1PG30.14 1978    1978   1978  1978  1978    1978
 1PG30.15 1778    1778   1778  1778  1778    1778
 1PG30.16    0       0      0     0   890       0
                                      ^^^ ?
Comment 3 Francois-Xavier Coudert 2007-09-19 09:54:42 UTC
(In reply to comment #0)
> Following a discussion on IRC with FX Coudert, I think the number of digits
> printed by default in formatted output should be increased by one.

Dominique, what is the output of your test program on ppc-darwin?

I otherwise agree with you: until (or unless) we get more clever formatted I/O, we should increase the default format. I'm waiting for the real(16) output and will then submit a patch.
Comment 4 Dominique d'Humieres 2007-09-19 12:38:45 UTC
Subject: Re:  Add one digit to the default formatted output

> Dominique, what is the output of your test program on ppc-darwin?
 real(4)
 default         808
 1PG20.6        1881
 1PG20.7         808
 1PG20.8           0

 real(8)
 default         1778
 1PG30.14        1978
 1PG30.15        1778
 1PG30.16           0

 real(16)
 default         1999
 1PG60.32        1999
 1PG60.33        1999
 1PG60.34        1999
 1PG60.35        1999
 1PG60.40        1999

but I have no way to know if the problem comes from NEAREST() 
or from the I/O.
Comment 5 Francois-Xavier Coudert 2007-10-01 23:28:40 UTC
Hum, the real(16) case looks like loads of fun to come!

$ cat k.f90
  integer, parameter :: k = 16
  character(80) :: buf
  real(k) :: xk, yk

  xk = 1.0_k - epsilon(xk)
  write (buf,'(1PG60.40)') xk
  read (buf,*) yk
  write (*,'(1PG60.40)') xk
  write (*,'(1PG60.40)') yk
  if (xk /= yk) print *, "Mismatch"
  print *

  xk = 1.0_k + epsilon(xk)
  write (buf,'(1PG60.40)') xk
  read (buf,*) yk
  write (*,'(1PG60.40)') xk
  write (*,'(1PG60.40)') yk
  if (xk /= yk) print *, "Mismatch"
  print *

  end
$ gfortran k.f90 && ./a.out
              1.9999999999999999999999999999999753481000    
              4.9303806576313237838233035330174139355000E-32
 Mismatch

              1.0000000000000000000000000000000246519000    
              1.0000000000000000000000000000000246519000    

(Also, in the case of real(16), you need to add an E4 to the format anyway.)
Comment 6 Francois-Xavier Coudert 2007-10-02 00:37:27 UTC
(In reply to comment #5)
> Hum, the real(16) case looks like loads of fun to come!

Part of it is simply a libc bug. There are numbers close to 1.0 and -1.0 that the darwin libc can't output properly:

$ cat k2.c 
#include <stdio.h>
int main (void)
{
  long double x;

  x = 0.99999999999999998L;
  printf ("%Lg\n", x);
  printf ("%.60Lg\n", x);

  x = -0.999999999999999999999999999999975L;
  printf ("%Lg\n", x);
  printf ("%.60Lg\n", x);
}
$ gcc k2.c && ./a.out
2
1.99999999999999997999999999999999548766360788798414493987417
-2
-1.99999999999999999999999999999997534809671184338108088348233

I've reported this to Apple (it is now bug #5516762). Otherwise, as far as I can tell, the default format width seems wide enough for real(kind=16).
Comment 7 Dominique d'Humieres 2007-10-02 08:29:14 UTC
Subject: Re:  Default formats for real input are not
 precise enough

> Part of it is simply a libc bug. There are numbers close to 1.0 and -1.0 that
> the darwin libc can't output properly:

Nice catch! The problem affects the vicinity of all the powers of 2, try
1.99999999999999998_16, 3.99999999999999998_16, ...
Comment 8 Francois-Xavier Coudert 2007-10-02 23:28:02 UTC
Subject: Bug 33469

Author: fxcoudert
Date: Tue Oct  2 23:27:51 2007
New Revision: 128967

URL: http://gcc.gnu.org/viewcvs?root=gcc&view=rev&rev=128967
Log:
	PR libfortran/33469

	* io/write.c (write_real): Widen the default formats.

	* gfortran.dg/default_format_1.f90: New test.
	* gfortran.dg/default_format_2.f90: New test.
	* gfortran.dg/namelist_print_1.f: Adjust expected output.
	* gfortran.dg/real_const_3.f90: Adjust expected output.

Added:
    trunk/gcc/testsuite/gfortran.dg/default_format_1.f90
    trunk/gcc/testsuite/gfortran.dg/default_format_2.f90
Modified:
    trunk/gcc/testsuite/ChangeLog
    trunk/gcc/testsuite/gfortran.dg/namelist_print_1.f
    trunk/gcc/testsuite/gfortran.dg/real_const_3.f90
    trunk/libgfortran/ChangeLog
    trunk/libgfortran/io/write.c

Comment 9 Francois-Xavier Coudert 2007-10-02 23:30:50 UTC
As I said, the default format is wide enough for powerpc-darwin. It's widely possible that the new testcase gfortran.dg/default_format_2.f90 fails there, though, due to the Apple printf() bug. Dominique, if you confirm that the new testcase fails, I'll XFAIL it. Thanks!
Comment 10 Dominique d'Humieres 2007-10-03 19:30:05 UTC
gcc/testsuite/gfortran.dg/default_format_1.f90 passes the test if I replace

        if (y /= x) res = res + 1

by

        z=0.0_k
        if (abs(x-y)>nearest(z,1.0_k)) res = res + 1

in the two places of test_r8, i.e. the tests fail for TINY() because the result differs from the original value by ~+/-5.0E-324: the smallest denormalized number.  I have noticed that this happens for an odd number of iterations of NEAREST(), but I did not figure out if the culprit was the WRITE() or the READ(). I have also seen the same problem with xlf and g95.

I suggest that the test includes this change and does not xfail on Darwin.

Comment 11 Dominique d'Humieres 2007-10-03 19:34:36 UTC
BTW I have forgotten to explain why I have to use an auxiliary variable 'z': if I usenearest(0.0_8,1.0_8); I get

default_format_1_db.f90:70.29:

        if (abs(x-y)>nearest(0.0_8,1.0_8)) print *, x, y, x-y
                            1
Error: Result of NEAREST underflows its kind at (1)
default_format_1_db.f90:84.29:

        if (abs(x-y)>nearest(0.0_8,1.0_8)) print *, x, y, x-y
                            1
Error: Result of NEAREST underflows its kind at (1)
default_format_1_db.f90:93.25:

  use test_default_format
                        1
Fatal Error: Can't open module file 'test_default_format.mod' for reading at (1): No such file or directory

Producing a hard error for an underflow is indeed an incredibely clever idea!-(

Reminds me of a compiler which crunched NaNs and Infs with delight, but aborted on underflow, very handy for waves simulations!