Search code examples
c++matrixeigenmatrix-multiplicationcomplex-numbers

Complex Matrix matrix multiplication using Eigen


I am doing complex matrix matrix product as shown in the code here. The Eigen * operator doesn't seem to give the right result. To verify I have also written another routine to find the product. From this it was seen that the imaginary part of the result was fine but the real part is inaccurate. Can someone let me know why?

Eigen::MatrixXcd U2 = Eigen::MatrixXcd::Random(N1, computed_rank);
Eigen::MatrixXcd V2 = Eigen::MatrixXcd::Random(computed_rank, N2);
Eigen::MatrixXcd W2 = Mat::Zero(N1,N2);
for (size_t k = 0; k < computed_rank; k++) {
    W2 += U2.col(k)*V2.row(k);
}
Eigen::MatrixXcd Q = U2*V2;
Eigen::MatrixXcd E2 = Q-W2;
cout << "Err Total: " << E2.norm()/W2.norm() << endl;
cout << "Err Real: " << E2.real().norm()/W2.real().norm() << endl;
cout << "Err Imag: " << E2.imag().norm()/W2.imag().norm() << endl;

which gave the result:

Err Total: 0.84969
Err Real: 1.17859
Err Imag: 1.18274e-16

I have observed that this problem didn't occur when this was the only piece of code in the project. But I have this as part of a bigger project, where it seems to fail.


Solution

  • It turns out that -ffast-math optimization flag of GNU is not compatible with all programs. I have used it while compiling the bigger project and not while doing the one that is in isolation. This is the reason it failed! For more details refer https://gcc.gnu.org/onlinedocs/gcc/Optimize-Options.html.