[Bug libfortran/48906] New: Wrong rounding results with -m32

[Bug libfortran/48906] New: Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

           Summary: Wrong rounding results with -m32
           Product: gcc
           Version: 4.7.0
            Status: UNCONFIRMED
          Severity: normal
          Priority: P3
         Component: libfortran
        AssignedTo: jvdelisle@gcc.gnu.org
        ReportedBy: jvdelisle@gcc.gnu.org


This PR is to carry forward the remaining issue identified in PR48602 comment
#37

Here is a test case showing the problem.

 print "(ru,g15.2)", .099d0 ! > 0.10E+00 expected 0.10
 print "(rc,g15.1)", .095d0 ! > 0.1E+00 expected 0.1
 print "(rc,g15.2)", .0995d0 ! > 0.10E+00 expected 0.10
 print "(ru,g15.3)", .0999d0 ! > 0.100E+00 expected 0.100

end

$ gfc bug.f03 
$ ./a.out 
       0.10    
        0.1    
       0.10    
      0.100    
$ gfc -O1 bug.f03 
$ ./a.out 
       0.10    
        0.1    
       0.10    
      0.100    
$ gfc -m32 -O1 bug.f03 
$ ./a.out 
       0.10E+00
        0.1E+00
       0.10E+00
      0.100E+00
$ gfc -m32 bug.f03 
$ ./a.out 
       0.10E+00
        0.1E+00
       0.10E+00
      0.100E+00

The problem is rooted in excess precision and that we use sensitive floating
point operations to determine the formatting parameters.  I am seeking a
solution that will avoid using floating point to do this.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #1 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-05-06 21:10:07 UTC ---
Repeating Thomas T. suggestion:

"We can treat FMT_G like FMT_E. After the rounding step, i.e. when the final
value of the variable e in output_float() is known (at the label 'skip'), if e
is within limits for F editing, revert the formatting to F and add the blanks
at the end." This is my plan A and I have a patch started.

Another idea would be to perform the floating point operations and tests with
mpfr.  This will avoid the -m32 issue.  As a side note, mpfr is compiled with
-ffloat-store on 32bit platforms by default.  I will take this idea as plan B
or C. ;)

[Bug libfortran/48906] Wrong rounding results with -m32

[kargl at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

kargl at gcc dot gnu.org changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
                 CC|                            |kargl at gcc dot gnu.org

--- Comment #2 from kargl at gcc dot gnu.org 2011-05-06 21:51:26 UTC ---
(In reply to comment #1)
> Another idea would be to perform the floating point operations and tests with
> mpfr.  This will avoid the -m32 issue.  As a side note, mpfr is compiled with
> -ffloat-store on 32bit platforms by default.  I will take this idea as plan B
> or C. ;)

If the "floating point operations and tests" are part of the runtime
library, then you cannot currently use mpfr because the final executable
is not linked with gmp and mpfr.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #3 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-05-10 12:53:07 UTC ---
Status: I have managed to eliminate the OUTPUT_FLOAT_G macro completely and all
floating point operations.  What remains is to adjust the trailing blanks.  The
patch is simple looking once you get there, though a bit subtle to sort out
initially, especially when I get only 20 or 30 minute intervals to work on this
every two or three days.  ;)

Stay tuned.

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #4 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-05-10 17:16:47 UTC ---
Way to go! I'll be happy to test.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #5 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-05-19 05:13:25 UTC ---
Thomas,

What are the correct results for this test case:


implicit none
integer, parameter :: RT = 8
write(*, "(rc,f11.2,4x,'<RC')") .995_RT
write(*, "(rc,g15.2,'<RC')") .995_RT
write(*, "(rd,f11.2,4x,'<RD')") .995_RT
write(*, "(rd,g15.2,'<RD')") .995_RT
write(*, "(ru,f11.2,4x,'<RU')") .995_RT
write(*, "(ru,g15.2,'<RU')") .995_RT
write(*,'(a)') "12345678901234567890"
end

I have finished re-factoring the code, but have some remaining breakages to
fix.  I just want to confirm what the above should give vs what I am getting at
the moment.  I will post a patch to test after I get all known breakage fixed.

Thanks.

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #6 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-05-19 08:05:32 UTC ---
(In reply to comment #5)

As a general rule, d specifies the number of significant digits in the result,
i.e. the number of digits counting from the first non-zero digit.

So in the example all numbers should be formatted with 2 significant digits,
i.e. 2 decimal digits if the value rounded to 2 significant digits is smaller
than 1, and 1 decimal digit if the value rounded to 2 significant digits is
larger than or equal to 1.

> implicit none
> integer, parameter :: RT = 8
> write(*, "(rc,f11.2,4x,'<RC')") .995_RT
> write(*, "(rc,g15.2,'<RC')") .995_RT
> write(*, "(rd,f11.2,4x,'<RD')") .995_RT
> write(*, "(rd,g15.2,'<RD')") .995_RT
> write(*, "(ru,f11.2,4x,'<RU')") .995_RT
> write(*, "(ru,g15.2,'<RU')") .995_RT
> write(*,'(a)') "12345678901234567890"
> end

       0.99    <RC
       0.99    <RC
       0.99    <RD
       0.99    <RD
       1.00    <RU
        1.0    <RU
12345678901234567890

and
  write(*, "(ru,f11.1,4x,'<RU')") .995_RT
gets
        1.0    <RU

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #7 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-05-19 12:58:34 UTC ---
It was the last case I was having trouble with. I also was getting different
results from different compilers, so I needed to confirm independently. I have
this part working now. Still have some bugs here.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #8 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-05-23 19:00:54 UTC ---
Status:  I am down to about 5 testsuite failures on the patch for this.  There
is a lot if interplay going on, so i will be factoring the code some as part of
the cleanup. Stay tuned.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #9 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-05-26 12:53:30 UTC ---
Created attachment 24364
  --> http://gcc.gnu.org/bugzilla/attachment.cgi?id=24364
Preliminary patch

I am getting two testsuite regressions with this patch. At least one of those I
think is an adjustment needed to the test case.  Possibly the same with the
other, but I am still testing.

After fixing the above, all that remains is a little code cleanup.  Please
begin testing this patch and report back.

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #10 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-05-27 07:15:16 UTC ---
The test cases from above still fail (but with a different result):
    print "(ru,g15.2)", .099d0 ! 1.0E-01 expected 0.10
    print "(rc,g15.1)", .095d0 ! 1.E-01 expected 0.1
    print "(rc,g15.2)", .0995d0 ! 1.0E-01 expected 0.10
    print "(ru,g15.3)", .0999d0 ! 1.00E-01 expected 0.100

Also, G0.d tests fail in the E formatting case, the chosen format is E25.(d)E3
and not the minimal width required for G0.d:
    print "(rc,g0.4,'<')", .1234d-3 ! fixed width, not G0.d format
    print "(rc,g0.4,'<')", .1234d-2 ! fixed width, not G0.d format
    print "(rc,g0.4,'<')", .1234d-1 ! fixed width, not G0.d format
    print "(rc,g0.4,'<')", .1234d0
    print "(rc,g0.4,'<')", .1234d1
    print "(rc,g0.4,'<')", .1234d2
    print "(rc,g0.4,'<')", .1234d3
    print "(rc,g0.4,'<')", .1234d4
    print "(rc,g0.4,'<')", .1234d5 ! fixed width, not G0.d format
    print "(rc,g0.4,'<')", .1234d6 ! fixed width, not G0.d format
    print "(rc,g0.4,'<')", .1234d7 ! fixed width, not G0.d format

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #11 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-05-27 08:55:14 UTC ---
The scale factor is applied to F editing, but it shouldn't. See Fortran 2008:

NOTE 10.20
The scale factor has no effect on output unless the magnitude of the datum to
be edited is outside the range that permits effective use of F editing.

    print "(rc,1p,g15.4e2,'<')", 0.1234 ! 1.2340E-01 expected 0.1234
    print "(rc,g15.4e2,'<')", 0.1234 ! 0.1234

Additionally, there is a strange digit mixup in this case:
    print "(rc,-1p,g15.4e2,'<')", 0.1234 ! 1.0234 expected 0.1234

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #12 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-05-27 13:12:32 UTC ---
The following examples fail:

    print "(rc,g10.2,'<')", 99.5 ! 10. expected 0.10E+03
    print "(rc,g10.2,'<')", 995. ! 1.0E+03 expected 0.10E+04
    print "(rc,g10.3,'<')", 999.5 ! 100. expected 0.100E+04
    print "(rc,g10.3,'<')", 9995. ! 1.0E+04. expected 0.100E+05

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #13 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-05-27 14:32:12 UTC ---
OK, thanks for test cases, these are very helpful.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

Jerry DeLisle <jvdelisle at gcc dot gnu.org> changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
  Attachment #24364|0                           |1
        is obsolete|                            |

--- Comment #14 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-05-30 18:48:27 UTC ---
Created attachment 24399
  --> http://gcc.gnu.org/bugzilla/attachment.cgi?id=24399
Updated preliminary patch

This updated patch resolves the issues identified so far.  Please continue
testing. I have regression tested and the testsuite case gfortran.dg/pr20755.f
fails because with the patch I am ignoring the scale factor. Perhaps I have
been too aggressive with that. Regardless, please continue testing while I
study the scale factor question a bit more.

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #15 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-05-31 06:05:19 UTC ---
I hadn't really thought of that, but now I see what a PITA a scale factor > 0
is going to be (scale factor < 0 is simply padding with zeros and shifting
digits):

We have to ignore the scale factor, until we decide that we use E editing, at
which point we need an additional digit in the result.

For efficiency: Maybe we can "decide that we use E editing" without actually
overwriting the original digits at that point yet, so we can go back when we
actually switch to E editing and re-do the rounding with the extra digit.

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #16 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-05-31 06:28:25 UTC ---
    print "(rc,g11.1)", 2. ! 0.2E+01 expected 2.
    print "(rc,g11.2)", 20. ! 0.20E+02 expected 20.
    print "(rc,g11.3)", 200. ! 0.200E+03 expected 200.
    print "(rc,g11.4)", 2000. ! 0.2000E+04 expected 2000.

    print "(rc,g11.2)", .04 ! 0.40 expected 0.40E-01
    print "(rc,g11.3)", .04 ! 0.400 expected 0.400E-01

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #17 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-05-31 12:30:00 UTC ---
As you can tell, finding the magic formula with this new approach is not a
trivial effort. ;) My fear is that there is an unlimited number of corner
cases.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

Jerry DeLisle <jvdelisle at gcc dot gnu.org> changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
  Attachment #24399|0                           |1
        is obsolete|                            |

--- Comment #18 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-06-01 05:09:12 UTC ---
Created attachment 24406
  --> http://gcc.gnu.org/bugzilla/attachment.cgi?id=24406
New updated patch

This updated patch takes care of Comment #16.

Also, the test case, pr20755.f was originally intended to pass only with
-std=legacy.  This is to mimic g77 which does not ignore the scale factor.  At
least recently, gfortran passes this test case regardless of the -std=.  I am
tempted to delete the test case, since it looks like it goes against current
standards.  Any opinions about this?

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #19 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-06-01 06:26:27 UTC ---
> Also, the test case, pr20755.f was originally intended to pass only with
> -std=legacy.  This is to mimic g77 which does not ignore the scale factor.  At
> least recently, gfortran passes this test case regardless of the -std=.  I am
> tempted to delete the test case, since it looks like it goes against current
> standards.  Any opinions about this?

Now you've lost me.

>From what I understand PR20755 was about G editing resetting the scale factor
in certain cases (to 0 in the example) and not restoring it afterwards (to 1 in
this case). This is a bug both in legacy and in Fortran standard modes. The bug
was fixed and a testcase was created. The result is standard-compliant. Why
would you want to remove the testcase? The point of having a testcase is to
prevent a reappearing of the same bug in the future, so maybe we should keep
it.

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #20 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-06-01 08:12:31 UTC ---
(In reply to comment #18)
> Created attachment 24406 [details]
> New updated patch
> 
> This updated patch takes care of Comment #16.

Unfortunately, now the other testcases fail again (although with still another
results):

print "(ru,g15.2)", .099d0 ! > 0.10E-01 expected 0.10
print "(rc,g15.1)", .095d0 ! > 0.1E-01 expected 0.1
print "(rc,g15.2)", .0995d0 ! > 0.10E-01 expected 0.10
print "(ru,g15.3)", .0999d0 ! > 0.100E-01 expected 0.100

print "(rc,g10.2)", 99.5 ! 0.10E+02 expected 0.10E+03
print "(rc,g10.2)", 995. ! 0.10E+03 expected 0.10E+04
print "(rc,g10.3)", 999.5 ! 0.100E+03 expected 0.100E+04
print "(rc,g10.3)", 9995. ! 0.100E+04 expected 0.100E+05

It seems there is still a bug that doesn't increment the exponent in the
"overflow" case, i.e. when going from 0.9999...E(x) to 0.1000...E(x+1)

The exponent test (decision between F and E editing) must be done after that! I
can not find this in the code right now.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #21 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-06-01 12:26:28 UTC ---
Regarding comment #18, the test case passes regardless 0f the -std=legacy
option specified at the beginning of the source file. I am questioning why that
option is ignored in current gfortran releases. Somewhere along the way,
something changed.  I am also asking, do we want the test case to pass under
-std=legacy only, or pass all the time.

Regarding, not incrementing the exponent, I missed that and will update on next
iteration.  thanks for reviews.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

Jerry DeLisle <jvdelisle at gcc dot gnu.org> changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
  Attachment #24406|0                           |1
        is obsolete|                            |

--- Comment #22 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-06-02 03:57:49 UTC ---
Created attachment 24414
  --> http://gcc.gnu.org/bugzilla/attachment.cgi?id=24414
Fixed revised patch

This patch fixes my omissions on the last patch.  In addition, studying the
test case pr20755.f, I see that it has the wrong number of significant digits
in the expected result.  I have revised it now so that it passes. I have also
added a new test case, to capture the various errors we found here while
developing the patch

Feel free to continue testing and if no new errors show up, I will submit to
the gfortran list for approval.  /me crosses his fingers.  ;)

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #23 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-06-03 12:55:21 UTC ---
Patch submitted to list for approval.

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #24 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-06-06 08:38:13 UTC ---
(In reply to comment #23)
> Patch submitted to list for approval.

For a scale factor 0, we are done. Good work, thank you!

A scale factor != 0 does not work yet, you wrote you are still working on it,
is that correct?

print "(-2pg12.3)", 0.02 ! 0.200E-01 expected 0.002E+01
print "(-1pg12.3)", 0.02 ! 0.200E-01 expected 0.020E+00
print "(0pg12.3)", 0.02 ! 0.200E-01
print "(1pg12.3)", 0.02 ! 0.200E-01 expected 2.000E-02
print "(2pg12.3)", 0.02 ! 0.200E-01 expected 20.00E-03

--- gcc/testsuite/gfortran.dg/pr20755.f    (revision 174320)
+++ gcc/testsuite/gfortran.dg/pr20755.f    (working copy)
@@ -5,8 +5,8 @@
       character*30 s

       write (s,2000) 0.0, 0.02
-      if (s .ne. "    0.00       2.000E-02") call abort
+      if (s .ne. "    0.00       0.200E-01") call abort
       write (s,2000) 0.01, 0.02
-      if (s .ne. "   1.000E-02   2.000E-02") call abort
+      if (s .ne. "   0.100E-01   0.200E-01") call abort
  2000 format (1PG12.3,G12.3)
       end

I don't agree with the changes to this testcase, since the scale factor does
not work yet, and the test was correct before.

--- gcc/testsuite/gfortran.dg/fmt_g0_6.f08    (revision 174320)
+++ gcc/testsuite/gfortran.dg/fmt_g0_6.f08    (working copy)
@@ -57,7 +57,7 @@ contains
             do dec = d, 0, -1
                 lower = 10.0_RT ** (d - 1 - dec) - r * 10.0_RT ** (- dec - 1)
                 upper = 10.0_RT ** (d - dec) - r * 10.0_RT ** (- dec)
-                if (lower <= mag .and. mag < upper) then
+                if (lower < mag .and. mag <= upper) then
                     write(fmt_f, "('R', a, ',F', i0, '.', i0, ',', i0, 'X')")
roundmode, w - n, dec, n
                     exit
                 end if

I don't agree with this change, since it was according to the Fortran standard
before (it does not actually make a difference for the values we currently
check in this test, though.)

Re: [Bug libfortran/48906] Wrong rounding results with -m32

[jerry DeLisle]

On 06/06/2011 01:38 AM, thenlich at users dot sourceforge.net wrote:
> For a scale factor 0, we are done. Good work, thank you!
>
> A scale factor != 0 does not work yet, you wrote you are still working on it,
> is that correct?

I am now.  ;)

>
> print "(-2pg12.3)", 0.02 ! 0.200E-01 expected 0.002E+01
> print "(-1pg12.3)", 0.02 ! 0.200E-01 expected 0.020E+00
> print "(0pg12.3)", 0.02 ! 0.200E-01
> print "(1pg12.3)", 0.02 ! 0.200E-01 expected 2.000E-02
> print "(2pg12.3)", 0.02 ! 0.200E-01 expected 20.00E-03
>
My confusion seems to be when scale factor is to be ignored and when not,  I 
will give the standard another read.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at charter dot net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #25 from jvdelisle at charter dot net 2011-06-06 12:09:48 UTC ---
On 06/06/2011 01:38 AM, thenlich at users dot sourceforge.net wrote:
> For a scale factor 0, we are done. Good work, thank you!
>
> A scale factor != 0 does not work yet, you wrote you are still working on it,
> is that correct?

I am now.  ;)

>
> print "(-2pg12.3)", 0.02 ! 0.200E-01 expected 0.002E+01
> print "(-1pg12.3)", 0.02 ! 0.200E-01 expected 0.020E+00
> print "(0pg12.3)", 0.02 ! 0.200E-01
> print "(1pg12.3)", 0.02 ! 0.200E-01 expected 2.000E-02
> print "(2pg12.3)", 0.02 ! 0.200E-01 expected 20.00E-03
>
My confusion seems to be when scale factor is to be ignored and when not,  I 
will give the standard another read.

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #26 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-06-06 12:27:38 UTC ---
(In reply to comment #25)
> My confusion seems to be when scale factor is to be ignored and when not,  I 
> will give the standard another read.

As it happens, you're not the only one to find (this part of) the Fortran
standard confusing and buggy, so if all goes as planned there will be a
rewrite:

---------------------------------------------------------------------
                                                    11-xxx
To: J3
From: John Reid and Thomas Henlich
Subject: Interp: G editing for reals
Date: 2011 June 5

---------------------------------------------------------------------

NUMBER: F08/xxxx
TITLE: G editing for reals
KEYWORDS: format, G editing
DEFECT TYPE: Erratum
STATUS: J3 consideration in progress

QUESTION:

1. Gw.d editing for a real value that is in the range (0.1,10**d) and 
is not near an integer power of 10 uses F editing to produce exactly 
the same value as 0PEw.d editing would produce. For values in this 
range that are near an integer power of 10, is it intended that F 
editing be used to produce exactly the same value as 0PEw.d editing 
would produce? 

The rules in 10.7.5.2.2 usually have this effect, but the following
examples illustrate exceptions for rounding UP and to ZERO.
  i. 
     print "(ru,g11.2)", -9.95 
  or
     print "(rz,g11.2)", -9.95 
  Here, 0PE11.2 editing would produce the value -0.99E+01, which can be 
  represented as -9.9. However, the standard requires F7.0 to be used, 
  which gives the result -9. Note that the standard requires
     print "(rd,g11.2)", 9.95
  and
     print "(rz,g11.2)", 9.95
  to give the result 9.9. 

  ii.
     print "(ru,0p,g11.2)", -99.5
  The standard requires 0PE11.2 editing to be used, which gives 
  -0.99E+02. This is representable as -99.

  iii. 
     print "(ru,0p,g11.2)", 99.
  The standard requires 0PE11.2 editing to be used, which gives 
  0.99E+02. This is representable as 99.


2. COMPATIBLE and NEAREST modes of rounding differ only when the two
nearest representable values are equidistant from the given value. The
similarity appears not to be represented in the second table. What is
meant by "if the higher value is even"?

3. Why is no account taken of the effects when PROCESSOR_DEFINED 
rounding is in effect? 


ANSWER:

1. Yes, this was the intention and it would be clearer for the standard
to state this directly. It would also be easier for implementers to 
implement. 

2. If the standard is rewritten as proposed in the first question, these
further problems would be resolved. 

3. If the standard is rewritten as proposed in the first question, 
PROCESSOR_DEFINED rounding would be covered. 

EDITS to 10-007r1:

[258:14-20] In 10.7.5.2.2, replace paragraph 4 by
"Otherwise, the method of representation in the output field depends
on the internal value being edited. Let N be the value that would be
output with 0PEw.dEe editing and let s be its exponent part unless N 
is identically 0 in which case let s be 1. Let k be the scale factor 
(10.8.5). Let b be a blank. For Gw.d editing, if s lies in the range
0 <= s <= d, F(w-4).(d-s),4('b') editing is used to represent N in the
output field; otherwise, kPEw.d editing is used to represent the
internal value. For Gw.dEe editing, if s lies in the range 0 <= s <= d, 
F(w-n).(d-s),n('b') editing where n=e+2 is used to represent N in the 
output field; otherwise, kPEw.dEe editing is used to represent the 
internal value." 


SUBMITTED BY: John Reid and Thomas Henlich 

HISTORY: 11-xxx    m195  Submitted

------------------------------------------------------------------------

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at charter dot net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #27 from jvdelisle at charter dot net 2011-06-06 12:36:21 UTC ---
>
> print "(-2pg12.3)", 0.02 ! 0.200E-01 expected 0.002E+01
> print "(-1pg12.3)", 0.02 ! 0.200E-01 expected 0.020E+00
> print "(0pg12.3)", 0.02 ! 0.200E-01
> print "(1pg12.3)", 0.02 ! 0.200E-01 expected 2.000E-02   < Too many significant digits?
> print "(2pg12.3)", 0.02 ! 0.200E-01 expected 20.00E-03  < Too many significant digits?

Should these last two cases by 2.00E-02 and 20.0E-03 ? Otherwise we seem to be 
adding an extra significant digit.

Help me understand this.

Jerry

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #28 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-06-06 12:41:55 UTC ---
I had a look at the code and what I think we should do is the following:

If G editing and a scale factor p != 0 is in effect, split the rounding step
into 2 steps:

First, check if overflow occurs (.9999.. => 1.0000...) but do not actually
overwrite the digits yet.

With this information, we determine the final value of e.

If it turns out we need F editing: Repeat the rounding as above (or remember if
overflow occured the first time), and this time actually overwrite the digits.

If it turns out we need E editing: Set the new value for the number of digits
according to p and repeat the rounding.

It may help to refactor some of the code into a separate function, since we
need to call it twice.

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #29 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-06-06 12:47:58 UTC ---
(In reply to comment #27)
> >
> > print "(-2pg12.3)", 0.02 ! 0.200E-01 expected 0.002E+01
> > print "(-1pg12.3)", 0.02 ! 0.200E-01 expected 0.020E+00
> > print "(0pg12.3)", 0.02 ! 0.200E-01
> > print "(1pg12.3)", 0.02 ! 0.200E-01 expected 2.000E-02   < Too many significant digits?
> > print "(2pg12.3)", 0.02 ! 0.200E-01 expected 20.00E-03  < Too many significant digits?
> 
> Should these last two cases by 2.00E-02 and 20.0E-03 ? Otherwise we seem to be 
> adding an extra significant digit.
> 
> Help me understand this.

Don't worry, it's not your fault... Say goodbye to the world of common sense,
welcome to the world of the Fortran standard! ;-)

"The scale factor k controls the decimal normalization (10.3.2, 10.8.5).  If −d
< k <= 0, the output field contains exactly |k| leading zeros and d − |k|
significant digits after the decimal symbol. If 0 < k < d + 2, the output field
contains exactly k significant digits to the left of the decimal symbol and d −
k + 1 significant digits to the right of the decimal symbol. Other values of k
are not permitted."

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #30 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-06-07 06:56:13 UTC ---
Created attachment 24454
  --> http://gcc.gnu.org/bugzilla/attachment.cgi?id=24454
Patch for scale factor. PUBLIC DOMAIN

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #31 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-06-07 06:56:54 UTC ---
(In reply to comment #30)
> Created attachment 24454 [details]
> Patch for scale factor. PUBLIC DOMAIN

A had a go at this. Feel free to improve.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #32 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-06-07 17:06:03 UTC ---
I managed to get pr20755.f fixed last night with only a minor change and I
found a problem elsewhere for fmt_g0_6.f08 which I am still working on.  I will
have a look at your patch in the next few days.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #33 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-06-10 16:16:19 UTC ---
The last test case I am working is fmt_g0_6.f08. 

The apparent failing case is:

      print "(rc,g15.2)", 0.995000_8

Which is resulting in 0.99 and we expect it to be 1.0.

However, the raw internal representation of 0.995 is not exact and is:

99499999999999999555

This places the value less then the threshold given in the table in the
Standard, giving the correct F format as:

      print "(f11.2,4x)", 0.995000_8  ! 0.99

Instead of 

      print "(f11.1,4x)", 0.995000_8  ! 1.0


Our F formatting in the test case is using the f11.1 and the logic does not
look at the third digit and dutifully rounds the value as it should.  The new G
formatting is using exact integer logic and dutifully does look at the third
digit.

If we do:

print "(rc,g15.2)", 0.9950001_8

The internal value does cross the 0.995 threshold and we do get the 1.0
results.

I think this issues is one of those splitting hairs things.  We could
artificially round the digits string early before we process it, but this would
effectively be rounding twice to get there.  I would prefer to revise the test
case like so:

  call check_all(0.995000000001_RT, 15, 2, 0)

and be done with this.

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #34 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-06-10 16:22:20 UTC ---
Additional note:  The standard states:

"Let N be the magnitude of the internal value"

The internal value is to be used to determine the conversion to F formatting. I
think this adds to my point.

I wonder if the standards committee knew the thresholds selected do not have an
exact binary representation so that the internal values will always be above or
below it.

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #35 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-06-10 16:56:02 UTC ---
(In reply to comment #33)
> The last test case I am working is fmt_g0_6.f08. 
> 
> The apparent failing case is:
> 
>       print "(rc,g15.2)", 0.995000_8
> 
> Which is resulting in 0.99 and we expect it to be 1.0.
> 
> However, the raw internal representation of 0.995 is not exact and is:
> 
> 99499999999999999555
> 
> This places the value less then the threshold given in the table in the
> Standard, giving the correct F format as:
> 
>       print "(f11.2,4x)", 0.995000_8  ! 0.99
> 
> Instead of 
> 
>       print "(f11.1,4x)", 0.995000_8  ! 1.0
> 
> 
> Our F formatting in the test case is using the f11.1 and the logic does not
> look at the third digit and dutifully rounds the value as it should.  The new G
> formatting is using exact integer logic and dutifully does look at the third
> digit.
> 
> If we do:
> 
> print "(rc,g15.2)", 0.9950001_8
> 
> The internal value does cross the 0.995 threshold and we do get the 1.0
> results.
> 
> I think this issues is one of those splitting hairs things.  We could
> artificially round the digits string early before we process it, but this would
> effectively be rounding twice to get there.  I would prefer to revise the test
> case like so:
> 
>   call check_all(0.995000000001_RT, 15, 2, 0)
> 
> and be done with this.

The result is correct, but our expectation was wrong.

0.995000_8 which is internally 0.99499999999999999555... is smaller than 0.995
(exact value) so it must output 0.99, not 1.0.

Please don't waste time with this test case (fmt_g0_6.f08), it is broken. I
will fix this some time.

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #36 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-06-10 17:01:34 UTC ---
(In reply to comment #34)
> Additional note:  The standard states:
> 
> "Let N be the magnitude of the internal value"
> 
> The internal value is to be used to determine the conversion to F formatting. I
> think this adds to my point.
> 
> I wonder if the standards committee knew the thresholds selected do not have an
> exact binary representation so that the internal values will always be above or
> below it.

Not always (e.g. 99.5 is the exact internal value), but it many cases, yes.

This and other problems with G editing will be fixed in the standard soon, if
all goes well.

I helped to prepare a request for this, which has now been submitted:

http://j3-fortran.org/doc/year/2011/11-174.txt

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #37 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-06-10 17:29:54 UTC ---
Updated patch posted for approval:

http://gcc.gnu.org/ml/fortran/2011-06/msg00097.html

Thomas, thanks for working the Standard issues and your testing.

[Bug libfortran/48906] Wrong rounding results with -m32

[thenlich at users dot sourceforge.net]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #38 from Thomas Henlich <thenlich at users dot sourceforge.net> 2011-07-01 09:22:48 UTC ---
The Fortran standards committee has voted to edit the standard:

http://j3-fortran.org/doc/meeting/195/11-174r2.txt

This makes our current approach standard compliant (with the corrigendum to be
published)

This also means an additional todo item for this bug is

"If d is zero, kPEw.0 or kPEw.0Ee editing is used for Gw.0 editing or Gw.0Ee
editing respectively."

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

--- Comment #39 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> 2011-07-02 01:32:42 UTC ---
Very Good, I will work toward this end.

[Bug libfortran/48906] Wrong rounding results with -m32

[dominiq at lps dot ens.fr]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

Dominique d'Humieres <dominiq at lps dot ens.fr> changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
             Status|UNCONFIRMED                 |WAITING
   Last reconfirmed|                            |2013-06-16
     Ever confirmed|0                           |1

--- Comment #40 from Dominique d'Humieres <dominiq at lps dot ens.fr> ---
This PR seems to have been fixed from 4.7.3 up to trunk. Could it be closed as
fixed or is there some issue I am missing?

[Bug libfortran/48906] Wrong rounding results with -m32

[jvdelisle at gcc dot gnu.org]

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=48906

Jerry DeLisle <jvdelisle at gcc dot gnu.org> changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
             Status|WAITING                     |RESOLVED
         Resolution|---                         |FIXED

--- Comment #41 from Jerry DeLisle <jvdelisle at gcc dot gnu.org> ---
Closing this long PR and opening a new one for the cleanup mentioned in comment
#38