CODEWITHSUNDEEP
c++matlabeigen
F1sher
F1sher

Reputation: 7310

Eigen Library - complex values in EigenSolver

I want to rewrite to C++ eig function from Matlab:

[V,D] = eig(A,B) produces a diagonal matrix D of generalized
eigenvalues and a full matrix V whose columns are the corresponding
eigenvectors so that A*V = B*V*D.

I get positive results for Matrix4d class:

pair<Matrix4d, Vector4d> eig(const Matrix4d& A, const Matrix4d& B)
{
const Eigen::GeneralizedSelfAdjointEigenSolver<Matrix4d> solver(A, B);

const Matrix4d V = solver.eigenvectors();
const Vector4d D = solver.eigenvalues();

return pair<Matrix4d, Vector4d>(V, D);
}

The problem I face now is that eig in Matlab is capable of using complex numbers - what I need in my algorithm.

How can I change it to receive same/very similar effect with usage of Matrix4cd class instead of Matrix4d?

Upvotes: 0

Views: 1688

Answers (2)

Avi Ginsburg
Avi Ginsburg

Reputation: 10596

If you don't/can't use auto in davidhigh's solution try this 

std::pair<Matrix4cd, Vector4d> eig(const Matrix4cd& A, const Matrix4cd& B)
{
    Eigen::GeneralizedSelfAdjointEigenSolver<Matrix4cd> solver(A, B);

    Matrix4cd V = solver.eigenvectors();
    Vector4d D = solver.eigenvalues();

    return std::make_pair(V, D);
}

Upvotes: 1

davidhigh
davidhigh

Reputation: 15468

Replace the MatrixXd by the complex version of it, Matrix4cd (which is actually a shortcut for Matrix< std::complex < double >, Dynamic, 4 >, see here):

std::pair<Matrix4cd, Vector4d> eig(const Matrix4cd& A, const Matrix4cd& B)
{
    Eigen::GeneralizedSelfAdjointEigenSolver<Matrix4cd> solver(A, B);

    Matrix4cd V = solver.eigenvectors();
    Vector4d D = solver.eigenvalues();

    return std::make_pair(V, D);
}

This is it already. Note however that you cannot pass arbitrary matrices, but the matrix A needs to be self-adjoint (also called hermitian) and B needs to be positive definite.

EDIT: thanks to the comment by @AviGinsburg, I corrected the mistake that the eigenvalues are real (and thus should be mapped to a real Eigen vector). Here is a C++14 implementation which gets around this problem (and so supports my stupidity ;-) ):

auto eig(const Matrix4cd& A, const Matrix4cd& B)
{
    Eigen::GeneralizedSelfAdjointEigenSolver<Matrix4cd> solver(A, B);

    auto V = solver.eigenvectors();
    auto D = solver.eigenvalues();

    return std::make_pair(V, D);
}

Upvotes: 2

Related Questions