Are these bugs for the GNU Fortran, g++, and clang compilers?

  • Context: Fortran 
  • Thread starter Thread starter Phyty
  • Start date Start date
  • Tags Tags
    Fortran Gfortran
Click For Summary

Discussion Overview

The discussion revolves around the behavior of floating-point arithmetic in the GNU Fortran, g++, and clang compilers, particularly focusing on the calculation of the smallest representable number and the implications of observed results in a loop. Participants explore whether the behavior observed constitutes a bug in the compilers or is a consequence of floating-point representation.

Discussion Character

  • Exploratory
  • Technical explanation
  • Debate/contested

Main Points Raised

  • One participant reports that after 1358 iterations of dividing a small number by 1.73, the value remains constant at approximately 4.94E-324, raising questions about the behavior of floating-point arithmetic.
  • Another participant asserts that the behavior described is not a bug and references floating-point arithmetic principles, suggesting that the results are expected due to the limitations of numerical representation.
  • Additional references to machine epsilon and floating-point information are provided, indicating the precision and limitations of different floating-point types.
  • A detailed Fortran program is shared that demonstrates how to determine the size and smallest values of different floating-point types, further illustrating the characteristics of floating-point representation.

Areas of Agreement / Disagreement

Participants do not reach a consensus. Some argue that the behavior is a result of floating-point arithmetic properties, while others question whether it indicates a bug in the compilers.

Contextual Notes

The discussion highlights the complexities of floating-point arithmetic, including potential misunderstandings about machine precision and representation limits. There is no resolution regarding the classification of the observed behavior as a bug or a feature of the compilers.

Phyty
Messages
1
Reaction score
4
TL;DR
I tried to calculate the smallest number in my machine with gfortran, g++ and clang. I did it by looping as:
a = a / 1.73
In gfortran after 1358 loops, I found that a = 4.9406564584124654E-324, and the value of a never changes from loop 1358 till the end 1500 loop.
It is questionable for me that a non zero real number here a = 4.9406564584124654E-324 divided by a non zero number here 1.73 will keep the same value.
I tried to calculate the smallest number in my machine with gfortran, g++ and clang. I did it by looping as:

a = a / 1.73

In gfortran after 1358 loops, I found that a = 4.9406564584124654E-324, and the value of a never changes from loop 1358 till the end 1500 loop.

It is questionable for me that a non zero real number here a = 4.9406564584124654E-324 divided by a non zero number here 1.73 will keep the same value. The phenomenon seems to the qualification of the ground state energy in a quantum system.

Is this bugs for the gfortran compiler?

I have tested similar codes with gfortran (GNU Fortran 8.2.0 2018, 5.40 201609), g++ (GNU 5.4.0 201609, clang-1000.10.44.4) on Linux (Ubuntu, KylinOS,) and Macintosh. I also tried a = a / 1.7.

All have similar results. A ground state qualification and unchanged value of a.

I guess that it is related to the binary system used by the computer system.

Here with my codes.
 
Technology news on Phys.org
From: https://www.math.utah.edu/~beebe/software/ieee/fpinfo.f ##-##
Fortran: 
      program fpinfo
************************************************************************
*     (Floating-point information)
*     Fortran program to report the size of single-, double-, and
*     quadruple-precision datatypes, and the floating-point precision
*     in bits, the corresponding machine epsilon and smallest
*     floating-point numbers for each.
*     (29-Nov-2001)
************************************************************************

      call ftest
      call dtest
      call qtest

      end      subroutine ftest
      real alog2, fstore
      real x
      real*16 fppow2
      character*36 note

      write (6,'(''single precision:'')')

      x = 1.0e+00
   10 if (fstore(1.0e+00 + x/2.0e+00) .ne. 1.0e+00) then
          x = x / 2.0e+00
          go to 10
      end if

      write (6,'(8x,''sizeof(REAL) =             '',i2)') 4

      write (6,'(8x,''FLT_MANT_DIG =             '',i2)')
     x    (1 + nint(-alog2(x)))

      if (x .eq. real(fppow2(-23))) then
           note = '[IEEE 754 32-bit macheps]'
      else
           note = '[not IEEE 754 conformant]'
      end if
      write (6,'(8x,''machine epsilon =          '', 1pe14.5e4, 1x, a)')
     x     x, note

      x = 1.0e+00
   20 if (fstore(x/2.0e+00) .ne. 0.0e+00) then
          x = x / 2.0e+00
          go to 20
      end if

      if (x .eq. real(fppow2(-149))) then
           note = '[IEEE 754 smallest 32-bit subnormal]'
      else if (x .eq. real(fppow2(-126))) then
           note = '[IEEE 754 smallest 32-bit normal]'
      else
           note = '[not IEEE 754 conformant]'
      end if
      write (6,'(8x,''smallest positive number = '', 1pe14.5e4, 1x, a)')
     x     x, note

      write (6,'()')

      end      subroutine dtest
      double precision dlog2, dstore
      double precision x
      real*16 fppow2
      character*36 note

      write (6,'(''double precision:'')')

      x = 1.0d+00
   10 if (dstore(1.0d+00 + x/2.0d+00) .ne. 1.0d+00) then
          x = x / 2.0d+00
          go to 10
      end if

      write (6,'(8x,''sizeof(DOUBLE PRECISION) = '',i2)') 8

      write (6,'(8x,''DBL_MANT_DIG =             '',i2)')
     x    (1 + nint(-dlog2(x)))

      if (x .eq. dble(fppow2(-52))) then
           note = '[IEEE 754 64-bit macheps]'
      else
           note = '[not IEEE 754 conformant]'
      end if
      write (6,'(8x,''machine epsilon =          '', 1pe14.5e4, 1x, a)')
     x     x, note

      x = 1.0d+00
   20 if (dstore(x/2.0d+00) .ne. 0.0d+00) then
          x = x / 2.0d+00
          go to 20
      end if

      if (x .eq. dble(fppow2(-1074))) then
           note = '[IEEE 754 smallest 64-bit subnormal]'
      else if (x .eq. dble(fppow2(-1022))) then
           note = '[IEEE 754 smallest 64-bit normal]'
      else
           note = '[not IEEE 754 conformant]'
      end if
      write (6,'(8x,''smallest positive number = '', 1pe14.5e4, 1x, a)')
     x     x, note

      write (6,'()')

      end      subroutine qtest
      real*16 qlog2, qstore
      real*16 x
      real*16 fppow2
      character*42 note

      write (6,'(''quadruple precision:'')')

      x = 1.0q+00
   10 if (qstore(1.0q+00 + x/2.0q+00) .ne. 1.0q+00) then
          x = x / 2.0q+00
          go to 10
      end if

      write (6,'(8x,''sizeof(REAL*16) =          '',i2)') 16

      write (6,'(8x,''LDBL_MANT_DIG =           '',i3)')
     x    (1 + nint(-qlog2(x)))

      if (x .eq. fppow2(-52)) then
           note = '[IEEE 754 64-bit macheps]'
      else if (x .eq. fppow2(-63)) then
           note = '[IEEE 754 80-bit macheps]'
      else if (x .eq. fppow2(-112)) then
           note = '[IEEE 754 128-bit macheps]'
      else
           note = '[not IEEE 754 conformant]'
      end if
      write (6,'(8x,''machine epsilon =          '', 1pe14.5e4, 1x, a)')
     x     x, note

      x = 1.0q+00
   20 if (qstore(x/2.0q+00) .ne. 0.0q+00) then
          x = x / 2.0q+00
          go to 20
      end if
      if (x .eq. fppow2(-1074)) then
           note = '[IEEE 754 smallest 64-bit subnormal]'
      else if (x .eq. fppow2(-1022)) then
           note = '[IEEE 754 smallest 64-bit normal]'
      else if (x .eq. fppow2(-16446)) then
           note = '[IEEE 754 smallest 80-bit subnormal]'
      else if (x .eq. fppow2(-16382)) then
           note = '[IEEE 754 smallest 80- and 128-bit normal]'
      else if (x .eq. fppow2(-16494)) then
           note = '[IEEE 754 smallest 128-bit subnormal]'
      else if (x .eq. fppow2(-16382)) then
           note = '[IEEE 754 smallest 128-bit normal]'
      else
           note = '[not IEEE 754 conformant]'
      end if
      write (6,'(8x,''smallest positive number = '', 1pe14.5e4, 1x, a)')
     x     x, note

      write (6,'()')

      end      real function fstore(x)
      real x
      fstore = x
      end      double precision function dstore(x)
      double precision x
      dstore = x
      end      real*16 function qstore(x)
      real*16 x
      qstore = x
      end      real*16 function fppow2(n)
      integer n
      integer p
      real*16 x
      if (n .lt. 0) then
          p = -n
          x = 1.0q+00 / 2.0q+00
      else
          p = n
          x = 2.0q+00
      end if
      fppow2 = 1.0q+00
   10 if (p .gt. 0) then
          p = p - 1
          fppow2 = fppow2 * x
          go to 10
      end if
      end      real function alog2(x)
      real x
      alog2 = alog(x)/alog(2.0e+00)
      end      double precision function dlog2(x)
      double precision x
      dlog2 = dlog(x)/dlog(2.0d+00)
      end      real*16 function qlog2(x)
      real*16 x
      qlog2 = qlog(x)/qlog(2.0q+00)
      end
 

Similar threads

Fortran What are common errors when compiling a simple Fortran code?
  • · Replies 2 ·
Replies
2
Views
3K
Fortran Compiling a Fortran 77 Code on Ubuntu
  • · Replies 3 ·
Replies
3
Views
5K
Fortran Compilation problem with gfortran
  • · Replies 9 ·
Replies
9
Views
9K
Fortran Help Compiling GNU Fortran Code with 2000 Lines
  • · Replies 22 ·
Replies
22
Views
7K
Fortran Help with Debugging Fortran Code | Plato Compiler
  • · Replies 4 ·
Replies
4
Views
3K
Fortran [Fortran] Integer overflow for a real number?
  • · Replies 3 ·
Replies
3
Views
4K
Fortran FORTRAN Problem Runtime Error: End of File
  • · Replies 2 ·
Replies
2
Views
9K
Fortran Can't get Fortran code to give me desired output file
  • · Replies 5 ·
Replies
5
Views
3K
Fortran Why Am I Getting a Parenthesis Error in FORTRAN Code Compilation?
  • · Replies 6 ·
Replies
6
Views
4K
Fortran Error when reading list, Fortran 77
  • · Replies 2 ·
Replies
2
Views
3K