topic I get different results from MKL and LAPACK. in Intel® oneAPI Math Kernel Library

I get different results from MKL and LAPACK.

Naomichi_H_ — Sat, 16 Jun 2018 11:48:24 GMT

Has anybody experienced getting different results from -mkl and -llapack? I've developed a fortran code that uses zgeev as in LWORK=-1 call zgeev('N','N',Ndim,U,Ndim,Eig,VL,1,VR,1,WORKdummy,LWORK,RWORK,Info) LWORK=INT(DBLE(WORKdummy(1))) allocate (WORK(LWORK)) call zgeev('N','N',Ndim,U,Ndim,Eig,VL,1,VR,1,WORK,LWORK,RWORK,INFO) For some parameters (NOT always), I get different results when I linked mkl and lapack. Here are the compiling commands: ifort -O3 -mkl -heap_arrays QWrandomEig.f90 ifort -O3 -heap_arrays -llapack QWrandomEig.f90 I'm attaching the source code as well as the input file and the output files, QWrandomEig_output_mkl.dat QWrandomEig_output_lapack.dat respectively, which are different from each other. From the anticipated behavior, I believe the latter (with lapack) is the correct result. More disturbingly, when I compile the code with -O0 -mkl as in ifort -O0 -heap_arrays -mkl QWrandomEig.f90 I get the following error message when I run it: Intel MKL ERROR: Parameter 2 was incorrect on entry to ZGEHD2. although the code runs till the end and the result QWrandomEig_output_mkl_O0.dat seems basically the same as the one with -O3 -mkl. I do not get this problem for other parameter sets, which is very annoying. Could someone guess what can be wrong with the code (or compiler for that matter)? OS is Mac OS 10.13.4 and the fortran version is ifort version 18.0.2 I appreciate your opinions.

Hi Naomichi H.

Ying_H_Intel — Mon, 25 Jun 2018 06:17:08 GMT

Hi Naomichi H.

I can reproduce the problem even under Linux. But i can't get same result as yours with liblapack. maybe different version.
we will further investigate and update you if any news.
if possible, could you please submit the issue to intel on-line service center: http://supporttickets.intel.com/, where for your private information.

Best Regards,
Ying

Thank you very much for

Naomichi_H_ — Mon, 25 Jun 2018 10:14:31 GMT

Thank you very much for checking! I really appreciate it. I don't know how to check the versions of my mkl and lapack, though I'm sure the mkl is the latest version.

I will submit the issue to intel.

Naomichi

Hi Naomichi,

Sergey_K_Intel1 — Fri, 06 Jul 2018 22:52:01 GMT

Hi Naomichi,

We reproduced the issue reported by you. Many thanks for catching the issue.

As concerns as the deviation in eigenvalues, unlike symmetric eigenvalue problems, eigenvalues and eigenvectors of non-symmetric problems are quite sensitive to round-off errors when the associated sets of eigenvectors are ill-conditioned. You can find many examples in the literature (see for example the book “Matrix algorithm”, Vol 2 “Eigensystems”, pp 37-38 by G.W.Stewart), showing that even tiny perturbations of a matrix can cause serious deviations of eigenvalues from the expected results. Because of this, the reference LAPACK and the MKL LAPACK cannot guarantee the repeatable results in the case of ill-conditioned eigensystems if the user changes compiler or compiler flags or platform or slightly perturbs systems. So all deviations in eigenvalues reported by you and other peoples in the thread are expected. Intel MKL can only guarantee the bit-to-bit correspondence of eigenvalues under the certain conditions. Please visit https://software.intel.com/en-us/articles/introduction-to-the-conditional-numerical-reproducibility-cnr for more details about conditional numerical results.

So one could ask which of the computed eigenvalues are correct (MKL or reference) and which are not. In my opinion the best way to answer the question is to draw the so-called pseudospectra (please visit http://www.cs.ox.ac.uk/pseudospectra/quickstart.html or https://math.nist.gov/MatrixMarket/spectral.html for the definition of pseudospectra, software tools and examples). The pseudospectra is analytical and graphical tool for investigating the spectrum of non-hermitian matrices and operators. Using these resources you can also find a lot of beautiful pictures showing that deviation from expected results can be very high. In other words, if we draw the region on the complex plane corresponding to 10^-16-pseudospectra, any complex number from the region can be considered as an eigenvalues.

So I believe that all results are provided by you are correct.

Thanks

All the best

Sergey

iDear Sergey

Naomichi_H_ — Sat, 07 Jul 2018 05:22:27 GMT

Dear Sergey Thank you for the reply, but I still wonder what the following error message means. >Intel MKL ERROR: Parameter 2 was incorrect on entry to ZGEHD2. which I got by compiling with -O0 option: >ifort -O0 -heap_arrays -mkl QWrandomEig.f90 Why does this message appear? I understand your point of the numerical instabilities. I can only tell you that the result of MKL LAPACK is certainly wrong from the point of view of my research. I was under the impression that it was somehow related to the error message above. best, Naomichi

Dear Naomichi,

Sergey_K_Intel1 — Tue, 10 Jul 2018 01:37:24 GMT

Dear Naomichi,

The error message mentioned by you doesn't have any impact on eigenvalues or eigenvectors. In order to find eigenvalues, the matrix is being reduced to the upper Hessenberg form and the reduction in MKL as in LAPACK is done by two methods: most of work is done by a block version and the rest of work by ZGEHD2. The last routine might issue an error message if all computational work is done by the block version. Of course we will fix the issue with this redundant error message.

By the way, we see the similar eigenvalues with the latest versions of the reference LAPACK. So probably the deviation is in due to changes in the reference LAPACK and we assume that you used an old version of the reference LAPACK for comparison. We would be much obliged to you if you let us know the version of reference LAPACK, compiler version and compiler flags used for building of LAPACK and BLAS binaries.

I also forgot to mention in the previous message, that the regular way to check correctness of non-symmetric eigenproblem is to compute eigenvalues and say the right eigenvectors by calling ZGEEV('N', 'V', ...). Then for each eigenpair (lambda(j) -eigenvalue, v(j) is the right eigenvector) you need to compute || A*v(j) - lambda(j)*v(j) || / || A||. If this norm is small enough for all eigenpairs, ZGEEV works correctly. Would you please check correctness on your system by computing these values?

Thanks in advance

Sergey