Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Silence error while solving the eigenvalues of a 2-by-2 matrix. #158

Open
BangShiuh opened this issue Aug 1, 2023 · 5 comments
Open

Silence error while solving the eigenvalues of a 2-by-2 matrix. #158

BangShiuh opened this issue Aug 1, 2023 · 5 comments

Comments

@BangShiuh
Copy link

BangShiuh commented Aug 1, 2023

The below code failed to obtain the eigenvalues for a 2-by-2 sparse matrix. At least an error message should pop up. (version 1.01)

#include <Eigen/Core>
#include <Eigen/SparseCore>
#include <Spectra/GenEigsSolver.h>
#include <Spectra/MatOp/SparseGenMatProd.h>
#include <iostream>
#include <Eigen/Sparse>
#include <Eigen/SparseLU>

using namespace Spectra;

double solver(Eigen::SparseMatrix<double> M, SortRule rule, int nKrylov)
{
    // Construct matrix operation object using the wrapper class SparseGenMatProd
    SparseGenMatProd<double> op(M);

    // Construct eigen solver object, requesting the largest three eigenvalues
    // how to set the size of Krylov subspace? from 5 to 10
    GenEigsSolver<SparseGenMatProd<double>> eigs(op, 1, nKrylov);

    // Initialize and compute
    eigs.init();
    int nconv = eigs.compute(rule);

    // Retrieve results
    Eigen::VectorXcd evalues;
    Eigen::MatrixXcd evecs;
    if(eigs.info() == CompInfo::Successful) {
        evalues = eigs.eigenvalues();
        evecs = eigs.eigenvectors();
    }

    // max eigenvalue
    std::complex<double> eigmax = evalues.coeff(0,0);

    return eigmax.real();
}

int main()
{
    // Define the matrix A
    Eigen::SparseMatrix<double> M(2, 2);
    M.insert(0, 0) = 1.0;
    M.insert(1, 1) = 2.0;
    M.insert(0, 1) = 0;
    M.insert(1, 0) = 0;

    double eigmax = solver(M,SortRule::LargestReal,2);
    std::cout << "Max real Eigenvalues found:\n" << eigmax << std::endl;

    return 0;
}
@yixuan
Copy link
Owner

yixuan commented Aug 15, 2023

An exception nev must satisfy 1 <= nev <= n - 2, n is the size of matrix should be thrown. Did you capture the exception somewhere else?

@BangShiuh
Copy link
Author

I am using the visual studio 2022 compiler, which does not show any message.

@yixuan
Copy link
Owner

yixuan commented Aug 16, 2023

It depends on how the compiler handles the exception, but indeed the program will terminate due to the error. You can explicitly use try-catch statements to detect such exceptions.

@BangShiuh
Copy link
Author

BangShiuh commented Aug 16, 2023

include/Spectra/JDSymEigsBase.h

void check_argument() const
    {
        if (m_number_eigenvalues < 1 || m_number_eigenvalues > m_matrix_operator.cols() - 1)
            throw std::invalid_argument("nev must satisfy 1 <= nev <= n - 1, n is the size of matrix");
    }

It seems that it is n-1 in the code.

@yixuan
Copy link
Owner

yixuan commented Aug 16, 2023

I think you were using GenEigsSolver, so this is the code: https://github.com/yixuan/spectra/blob/master/include/Spectra/GenEigsBase.h#L336. As a result, for any 2x2 matrix, GenEigsSolver will raise an error.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants