Is this expected?

1

Code:

! test_eqv.f90                                                                                                                                       
program test_eqv                                                                                                                        
implicit none                                                                                                                           
real :: x, y                                                                                                                            
character(len=100) :: str                                                                                                               
x = huge(x)                                                                                                                             
x = x**2 / x**3       ! Evaluates to NaN. Intended.                                                                                                                
write (str, *) x                                                                                                                        
read (str, *) y                                                                                                                         
write (*, *) x, y, str                                                                                                                  
write (*, *) x <= 0.0, y <= 0.0, (x <= 0.0 .eqv. y <= 0.0), (x <= 0.0 .neqv. y <= 0.0)                                                  
                                                                                                                                        
if (x <= 0.0 .eqv. y <= 0.0) then                                                                     
    write (*, *) 'Right.'                                                                                                               
else                                                                                                                                    
    write (*, *) 'Wrong!'                                                                                                               
    stop 1                                                                                                                              
end if                                                                                                                                  
                                                                                                                                        
end program test_eqv

Test:

$ uname -a && gfortran --version && gfortran test_eqv.f90 -Ofast && ./a.out
Linux zX11 6.2.0-32-generic #32~22.04.1-Ubuntu SMP PREEMPT_DYNAMIC Fri Aug 18 10:40:13 UTC 2 x86_64 x86_64 x86_64 GNU/Linux
GNU Fortran (Ubuntu 13.1.0-8ubuntu1~22.04) 13.1.0
Copyright (C) 2023 Free Software Foundation, Inc.
This is free software; see the source for copying conditions.  There is NO
warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.

              NaN              NaN               NaN                                                                                   
 F F F F
 Wrong!
Note: The following floating-point exceptions are signalling: IEEE_INVALID_FLAG
STOP 1

Is this expected? I understand that -Ofast affects the evaluation of floating-point expressions, but does it affect things like eqv as well?

Thanks.

2

The execution of any numeric operation whose result is not defined by the arithmetic used by the processor is prohibited.

You’re sequence x = huge(x); x = x**2 / x**3 is not defined by the arithmetic used by the processor, so technically anything that happens after is undefined behavior. I agree that the observed behavior is surprising though.

1 Like
3

It should not be surprising, the intermediates that are computed are clearly beyond the range of a representable floating point value.

There is also one other potential surprise. Suppose your compiler recognizes the expression and decides to simplify it as

x = 1.0 / huge(x)

That might look like it is an expression that can be evaluated correctly and safely, but it really isn’t. With most floating point formats, including the IEEE ones, that is a denormal number. That is, the exponent is a little too small for the number to be represented exactly. The way the result is treated depends on the compile time and run time options to handle denormal floating point numbers, so you might get anything from segment fault to an exception being signaled to just silently computing the low-precision result with no notice of the precision loss (which in this case would be just a couple of bits).

4

The behavior of the floating point operations is not the surprising part. It’s the subsequent behavior of the logical operations. OP has shown that

somehow manage to all evaluate to false, which seems to be contradictory. I.e. how can false not be equivalent to false, and at the same time not be not equivalent.

But at the same time, perhaps it’s treating the answer to NaN <= 0.0 as ¯\_(ツ)_/¯ and so ¯\_(ツ)_/¯ .eqv. ¯\_(ツ)_/¯ is also ¯\_(ツ)_/¯. And then ¯\_(ツ)_/¯ is just printed as F. Like I said, we’re in UB land here so anything is possible.

3 Likes
5

I experimented a little with this, and a workaround for the .eqv. problem is to define the logical variables

   logical :: qx, qy
   ...
   qx = x <= 0.0
   qy = y <= 0.0

and then use those variables in the various write statements and expressions. This seems to work correctly, even with aggressive compiler options.

$ gfortran --version -ffast-math eqv.f90 && a.out
GNU Fortran (Fink gcc11 11.3.0-1) 11.3.0
Copyright (C) 2021 Free Software Foundation, Inc.
This is free software; see the source for copying conditions.  There is NO
warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.

              NaN              NaN               NaN                                                                                   
 F F F F
 F F T F
 Right.

The first line shows the results of the original expressions (which is incorrect), while the remaining (correct) lines are with the logical variables.

1 Like
6

Thank @RonShepard for the suggestion.

Indeed, as long as I can confirm that the behaviour of gfortran -Ofast regarding .eqv. illustrated in my example is unexpected or questionable, then I will stop worrying and will not strive to circumvent it.

This is quite different from the problem discussed in another thread, where the behaviour of gfortran -Ofast regarding floating point computation is totally expected and I understand why it happens, although I would like workarounds to “fix” it.