In the attached test case (a.cpp), I'm trying to use the scalapack PZGESVD routine, but it crashes for a reason I cannot figure out. I'm calling from C++. If I call PDGESVD on a real (double precision) matrix, then everything works fine. I've attached the screen output from both a real compile (real.png) and a complex compile (complex.png). The only difference between the two compiles is that the real version defines dtype=1 at line 94 and the complex version defines dtype=2 at line 94. Also, the PZLAPRNT routine produces weird results unless you comment out line 118 (the PZGESVD call still crashes though).
I'm using MKL version 11.1.0 and intel MPI version 4.1 on a linux system.
Am I calling the routine incorrectly, or is this a library issue? Thanks in advance.
Gennady Fedorov (Intel) wrote:
ok, we will try to check how it works on our side.
Also, the mkl manual documentation for PZGESVD/PCGESVD is somewhat confusing to me for the following reasons:
1. The c-interface for these routines are do not include the 'rwork' parameter but the fortran interface does. This seems odd.
2. The documentation states that the 's'-parameter and 'rwork'-parameter for PZGESVD and PCGESVD are complex/MKL_Complex8 and double complex/MKL_Complex16 for PCGESVD and PZGESVD, respectively. The stock fortran Scalapack implementation available at nelib states that the 's' and 'rwork' parameters are real and double precision for PCGESVD and PZGESVD, respectively. We have used the MKL PZGESVD/PCGESVD fortran interfaces extensively and get correct answers using real arrays for s/rwork and not complex arrays, which conflicts with the MKL documentation. If you use complex arrays from fortran for s/rwork, you do not get the proper answers in s which has the consecutive singular values stuffed into the real/imaginary parts of the complex s (which is incorrect).
If the c-interface functions are redirecting to the fortran versions under-the-hood, I'm not sure how they could work properly with the given interface unless temporary memory is being generated under-the-hood for rwork and the s-values are converted from real-to-complex.
Attached is a short test with screen output of the possibility that the documentation for PZGESVD/PCGESVD is incorrect as well as the c-interface in mkl_scalapack.h. I removed mkl_scalapack.h and replaced it with an interface for pzgesvd_ that follows the stock netlib pzgesvd implementation. In this case, the error I was originally seeing no longer occurs and the result is correct.