[Bug fortran/100778] New: Get SIGFPE on simple test with -fpe-trap=invalid and SLP vectorization ON, with gfortran 11.1.0 on x86_64

["gabrielle.hugo at cern dot ch" <gcc-bugzilla@gcc.gnu.org>]

https://gcc.gnu.org/bugzilla/show_bug.cgi?id=100778

            Bug ID: 100778
           Summary: Get SIGFPE on simple test with -fpe-trap=invalid and
                    SLP vectorization ON, with gfortran 11.1.0 on x86_64
           Product: gcc
           Version: 11.1.0
            Status: UNCONFIRMED
          Severity: normal
          Priority: P3
         Component: fortran
          Assignee: unassigned at gcc dot gnu.org
          Reporter: gabrielle.hugo at cern dot ch
  Target Milestone: ---

Hello,

Updating gfortran from 10.3.0 to 11.1.0 leads to SIGFPE in a Fortran 77 project
(all things else equal).

The issue can actually be isolated and reproduced on x86_64 with the equivalent
simple test.f [1]: it runs smoothly when compiled with gfortran 10.3.0, and
leads to SIGFPE with gfortran 11.1.0.

To reproduce:
gfortran -O3 -Wall -Wextra -g -ffpe-trap=invalid test.f -o test
./test

The FPE is thrown from `x = x / norm` (while the branch should not have been
entered).

The issue disappears if additionally compiling with -fno-tree-slp-vectorize, or
if placing a printout in between `x = x / norm` and ` y = y / norm`.
The issue does not appear on ARM (tested on Apple M1).

Thanks in advance


[1] test.f:

      program test
      implicit none  

      double precision input(2)
      input(1) = 0.D+00
      input(2) = 0.D+00

      call normalize(input)

      end


      subroutine normalize(input)
      implicit none

      double precision input(2)
      double precision x, y, norm

      common / coordinates / x, y
      save / coordinates /

      x = input(1)
      y = input(2)
      norm = sqrt(x**2 + y**2)

      if (norm .gt. 1.D+00) then
         x = x / norm
         y = y / norm
         write (*, *) 'Normalized x and y.'
      end if

      write (*, *) 'Ended as expected.' 

      end

[2] gcc versions:

gcc version 10.3.0 (GCC)
Configured with:
/build/dkonst/gcc-clang/build/contrib/gcc-10.3.0/src/gcc/10.3.0/configure
--prefix=/build/dkonst/gcc-clang/lcgcmake-install/gcc/10.3.0/x86_64-centos7
-with-system-zlib --disable-multilib --enable-languages=all --with-gnu-ld
--with-gnu-as
Thread model: posix
Supported LTO compression algorithms: zlib

gcc version 11.1.0 (GCC)
Target: x86_64-pc-linux-gnu
Configured with:
/build/dkonst/gcc-clang-2/build/contrib/gcc-11.1.0/src/gcc/11.1.0/configure
--prefix=/build/dkonst/gcc-clang-2/lcgcmake-install/gcc/11.1.0/x86_64-centos7
-with-system-zlib --disable-multilib --enable-languages=all --with-gnu-ld
--with-gnu-as
Thread model: posix
Supported LTO compression algorithms: zlib

Tested on Centos7 and Fedora33.

[3] objdump -S test
gfortran 10.3.0 case:

      x = input(1)
  4011d9:       f2 0f 10 0f             movsd  (%rdi),%xmm1
      y = input(2)
  4011dd:       f2 0f 10 57 08          movsd  0x8(%rdi),%xmm2
      x = input(1)
  4011e2:       66 0f 10 27             movupd (%rdi),%xmm4
      norm = sqrt(x**2 + y**2)
  4011e6:       66 0f 28 c1             movapd %xmm1,%xmm0
  4011ea:       66 0f 28 da             movapd %xmm2,%xmm3
  4011ee:       f2 0f 59 c1             mulsd  %xmm1,%xmm0
      x = input(1)
  4011f2:       0f 29 25 77 2e 00 00    movaps %xmm4,0x2e77(%rip)        #
404070 <coordinates_>
      norm = sqrt(x**2 + y**2)
  4011f9:       f2 0f 59 da             mulsd  %xmm2,%xmm3
  4011fd:       f2 0f 58 c3             addsd  %xmm3,%xmm0
  401201:       f2 0f 51 c0             sqrtsd %xmm0,%xmm0

      if (norm .gt. 1.D+00) then
  401205:       66 0f 2f 05 53 0e 00    comisd 0xe53(%rip),%xmm0        #
402060 <options.2.0+0x20>
  40120c:       00 
  40120d:       77 51                   ja     401260 <normalize_+0x90>
  40120f:       48 8b 1d 52 0e 00 00    mov    0xe52(%rip),%rbx        # 402068
<options.2.0+0x28>
  401216:       48 89 e5                mov    %rsp,%rbp
         x = x / norm
         y = y / norm
         write (*, *) 'Normalized x and y.'
      end if

gfortran 11.1.0 case:
   x = input(1)
  4011d9:       66 0f 10 0f             movupd (%rdi),%xmm1
  4011dd:       66 0f 28 d9             movapd %xmm1,%xmm3
  4011e1:       66 0f 28 c1             movapd %xmm1,%xmm0
  4011e5:       0f 29 0d 84 2e 00 00    movaps %xmm1,0x2e84(%rip)        #
404070 <coordinates_>
      y = input(2)
      norm = sqrt(x**2 + y**2)
  4011ec:       f2 0f 59 c1             mulsd  %xmm1,%xmm0
  4011f0:       66 0f 15 db             unpckhpd %xmm3,%xmm3
  4011f4:       66 0f 28 d3             movapd %xmm3,%xmm2
  4011f8:       f2 0f 59 d3             mulsd  %xmm3,%xmm2
  4011fc:       f2 0f 58 c2             addsd  %xmm2,%xmm0
  401200:       f2 0f 51 c0             sqrtsd %xmm0,%xmm0

      if (norm .gt. 1.D+00) then
  401204:       66 0f 2f 05 54 0e 00    comisd 0xe54(%rip),%xmm0        #
402060 <options.2.0+0x20>
  40120b:       00 
         x = x / norm
  40120c:       66 0f 28 d0             movapd %xmm0,%xmm2
  401210:       66 0f 14 d2             unpcklpd %xmm2,%xmm2
  401214:       66 0f 5e ca             divpd  %xmm2,%xmm1
      if (norm .gt. 1.D+00) then
  401218:       77 4e                   ja     401268 <normalize_+0x98>
  40121a:       48 8b 1d 47 0e 00 00    mov    0xe47(%rip),%rbx        # 402068
<options.2.0+0x28>
  401221:       48 89 e5                mov    %rsp,%rbp
         y = y / norm
         write (*, *) 'Normalized x and y.'
      end if

Here, divpd is executed before ja.