[Mplapack-devel] possibly a bug in Cheev
Status: Pre-Alpha
Brought to you by:
nakatamaho
From: Akira S. <aki...@ni...> - 2012-07-26 10:21:31
|
Dear Dr. Maho NAKATA, I am a postdoc working on quantum information and the matrices I have to handle are quite often Hermitian matrices with some eigenvalues degenerate or close to each other. I am trying several libraries as well as writing my own library for matrix calculations required in my research. I would like to report some issue as I think this is possibly a bug of MPACK. Or, it might be a problem in my code or the environment. I tried Cheev of MPACK with GMP but could not get correct results for some simple matrices. For example, it does not diagonalize the 3 x 3 matrix 1 0 1 0 0 0 1 0 1 and gives me "#Rlaev2 Checkpoint 13 Not checked". It works fine for diag(1,0,0) and a matrix filled with 1's among others. I am using mpack-0.7.0 on Fedora 15 where GMP 4.3.2 exists. Thanks for your attentions. Regards, Akira SAITOH, NII, Japan =========================================== [saitoh@localhost test_cheev]$ g++ test.cpp -L/usr/local/lib -lmlapack_gmp -lmblas_gmp -lmblas_gmp_ref -lmlapack_gmp_ref -lgmp -lmpfr -lgmpxx -lgomp -g /usr/bin/ld: warning: libgmp.so.3, needed by /usr/lib/gcc/x86_64-redhat-linux/4.6.3/../../../../lib64/libmpfr.so, may conflict with libgmp.so.10 [saitoh@localhost test_cheev]$ LD_LIBRARY_PATH=/usr/local/lib ./a.out 1 0 1 0 0 0 1 0 1 A = 1 0 1 0 0 0 1 0 1 #Rlaev2 Checkpoint 13 Not checked [saitoh@localhost test_cheev]$ LD_LIBRARY_PATH=/usr/local/lib ./a.out 1 0 0 0 0 0 0 0 0 A = 1 0 0 0 0 0 0 0 0 Eigenvalues: 0 0 1 Unitary matrix: 0+i*0 0+i*0 1+i*0 1+i*0 0+i*0 0+i*0 0+i*0 1+i*0 0+i*0 =========================================== test.cpp: =========================================== #include <gmpxx.h> #include <mpack/mpack_config.h> #include <mpack/mpc_class.h> #include <mpack/mlapack_gmp.h> #include <cstdlib> #include <iostream> int main (int argc, char* argv[]) { mpf_set_default_prec(512); mpackint n = 3; mpc_class *A = new mpc_class [n * n]; mpackint lda = n; mpf_class *w = new mpf_class [n]; mpc_class *work = new mpc_class [2 * n - 1]; mpackint lwork = 2 * n - 1; mpf_class *rwork = new mpf_class [3 * n - 2]; mpackint info; for (int i = 1; i < argc && i < n * n + 1; i++) A[i-1] = ::atof(argv[i]); std::cout << "A = " << std::endl; for (int i = 0; i < n; i ++) { for (int j = 0; j < n; j++) std::cout << A[i + n * j].real().get_d() << " "; std::cout << std::endl; } Cheev ("V", "U", n, A, lda, w, work, lwork, rwork, &info); std::cout << "Eigenvalues: " << std::endl; for (int i = 0; i < n; i++) std::cout << w[i].get_d() << " "; std::cout << std::endl; std::cout << "Unitary matrix:" << std::endl; for (int i = 0; i < n; i ++) { for (int j = 0; j < n; j++) std::cout << A[i + n * j].real().get_d() << "+i*" << A[i + n * j].imag().get_d() << " "; std::cout << std::endl; } delete [] A, w, work, rwork; return 0; } =========================================== |