GeezerVison is showing github.com — Bug in zheevd · Issue #43 · Reference-LAPACK/lapack
DescriptionThere is an issue with the routine zheevd when compiled with Intel fortran using the optimization -O2 or -O3.
The problem occurs at line 285 of dlaed4:
DELTA( J ) = ( D( J )-D( I ) ) - TAU
which evaluates to zero because of the compiler not respecting the parenthesis. This results in a division-by-zero later on. Including -assume protect_parens fixes the problem.
Here is a simple program which produces the error:
program test
implicit none
integer, parameter :: n=27
integer i,j
integer liwork,lrwork,lwork,info
integer, allocatable :: iwork(:)
real(8), allocatable :: w(:),rwork(:)
complex(8), allocatable :: a(:,:),work(:)
liwork=5*n+3
lrwork=2*n**2+5*n+1
lwork=n**2+2*n
allocate(w(n),a(n,n))
allocate(iwork(liwork),rwork(lrwork),work(lwork))
a(:,:)=0.d0
do i=1,n
a(i,i)=1.d0
do j=i+1,n
a(i,j)=cos(dble(i))
end do
end do
call zheevd('V','U',n,a,n,w,work,lwork,rwork,lrwork,iwork,liwork,info)
print *,'info',info
end program
I'm not sure whether this is a compiler, LAPACK or user issue. I've added -assume protect_parens to our code (elk.sourceforge.net) and made zheev as the default eigenvalue solver for the moment. Activitypoulson commented on Sep 2, 2016 jeffhammond commented on Sep 2, 2016 I cannot reproduce this.
$ ifort -mkl -O2 fma-parens.f90 && ./a.out
info 0
$ ifort -mkl -O3 fma-parens.f90 && ./a.out
info 0
$ ifort -mkl -Ofast -xCORE-AVX2 -fma -fp-model fast=2 fma-parens.f90 && ./a.out
info 0
$ ifort --version
ifort (IFORT) 16.0.2 20160204
Copyright (C) 1985-2016 Intel Corporation. All rights reserved.
What compiler version are you using? jeffhammond commented on Sep 2, 2016 Oh, I see, the reproducer requires me to build LAPACK from source. Sorry.
At least it is good to know that MKL doesn't have any issue here. langou commented on Sep 5, 2016 Thanks Jay and Jeff. As Kay suggested, I have included the -assume protect_parens in the make.inc.ifort file in the INSTALL directory. Hopefully, that should do the trick. Cheers, Julien. jeffhammond commented on Sep 6, 2016 For what it's worth, this page discusses this in quite a bit of detail, including multiple pages on reassociation.
I attached the latest version of the PDF (FP_Consistency_070816.pdf), which was updated for version 17 of the Intel compilers. langou commented on Sep 12, 2016 Here are some comments from Sven Hammarling.
"Perhaps not surprisingly, I would phrase it the other way round that compiler optimizations should not get in the way of code writers.
It is challenging enough to write quality numerical software without having to worry about what the compiler might do to compromise numerical stability. I would hope that users would prefer correct answers over speed.
That is most certainly not to say that code writers should not strive for efficiency, but that compilers should respect the integrity of their code.
Under "Integrity of parentheses" the Fortran standard explicitly says that any expression in parentheses shall be treated as a data entity. The following example is given:
NOTE 7.16
For example, in evaluating the expression A + (B – C) where A, B, and C
are of numeric types, the difference of B and C shall be evaluated
before the addition operation is performed; the processor shall not
evaluate the mathematically equivalent expression (A + B) – C.
So, I still feel strongly that I should not need to use a switch to ensure that parentheses are not discarded. Parentheses being discarded should only be allowed by a positive action, that is a switch that allows such an action." jeffhammond commented on Sep 12, 2016
I would hope that users would prefer correct answers over speed.
It depends on how you define correct. We are talking about floating point, so accuracy is not a boolean but a threshold.
Scientific computing users want the code to run fast and not give incorrect answers. The vast majority of scientific applications are not rendered incorrect by changes in the last few bits of the significand, particularly since no scientific application I have ever seen is rigorous to all the bits anyways.
It so happens that in this particular case, reassociation causes noticeable incorrectness, but this is exceptional. Thus, it makes sense to cure the exception rather than disable productize optimizations (performance-increasing but not correctness-damaging code transformations).
And you can't ask the compiler to know when reassociation is dangerous, because it is strictly data-dependent, and compilation of Fortran is done offline. In a dynamic language like Python or Julia, one might be able to do something better. jeffhammond commented on Sep 12, 2016 In any case, even if the Intel Fortran compiler turns off productive optimizations everywhere because of a small number of reports where it causes a problem [1], there will be many versions of the compiler that contain this feature, and LAPACK will compile incorrectly with all of them. You have a way to solve this problem, if you choose to use it.
- If compilers disabled every optimization that broke at least one code, I suspect no compilers would optimize any code ever.
jeffhammond commented on Jan 12, 2017 I raised these issues with the Intel compiler team and there are open support issues. Will update as appropriate. weslleyspereira commented on May 23, 2023 I couldn't reproduce this issue using ifort version 2021.9.0. Maybe this was solved? I am inclined to close this issue since, as I understand, it is related to compiler optimizations and not to heevd. Reacted by langou and HOMOvsLUMO
MetadataMetadataDevelopmentNo branches or pull requestsIssue actions
|