Fix GeneralizedEigenSolver::info() and Asserts

This commit is contained in:
Arthur 2022-08-25 22:05:04 +00:00 committed by Rasmus Munk Larsen
parent 714678fc6c
commit a7c1cac18b
2 changed files with 53 additions and 17 deletions

View File

@ -121,8 +121,8 @@ template<typename MatrixType_> class GeneralizedEigenSolver
: m_eivec(),
m_alphas(),
m_betas(),
m_valuesOkay(false),
m_vectorsOkay(false),
m_computeEigenvectors(false),
m_isInitialized(false),
m_realQZ()
{}
@ -136,8 +136,8 @@ template<typename MatrixType_> class GeneralizedEigenSolver
: m_eivec(size, size),
m_alphas(size),
m_betas(size),
m_valuesOkay(false),
m_vectorsOkay(false),
m_computeEigenvectors(false),
m_isInitialized(false),
m_realQZ(size),
m_tmp(size)
{}
@ -158,8 +158,8 @@ template<typename MatrixType_> class GeneralizedEigenSolver
: m_eivec(A.rows(), A.cols()),
m_alphas(A.cols()),
m_betas(A.cols()),
m_valuesOkay(false),
m_vectorsOkay(false),
m_computeEigenvectors(false),
m_isInitialized(false),
m_realQZ(A.cols()),
m_tmp(A.cols())
{
@ -179,7 +179,8 @@ template<typename MatrixType_> class GeneralizedEigenSolver
* \sa eigenvalues()
*/
EigenvectorsType eigenvectors() const {
eigen_assert(m_vectorsOkay && "Eigenvectors for GeneralizedEigenSolver were not calculated.");
eigen_assert(info() == Success && "GeneralizedEigenSolver failed to compute eigenvectors")
eigen_assert(m_computeEigenvectors && "Eigenvectors for GeneralizedEigenSolver were not calculated");
return m_eivec;
}
@ -203,7 +204,7 @@ template<typename MatrixType_> class GeneralizedEigenSolver
*/
EigenvalueType eigenvalues() const
{
eigen_assert(m_valuesOkay && "GeneralizedEigenSolver is not initialized.");
eigen_assert(info() == Success && "GeneralizedEigenSolver failed to compute eigenvalues.");
return EigenvalueType(m_alphas,m_betas);
}
@ -214,7 +215,7 @@ template<typename MatrixType_> class GeneralizedEigenSolver
* \sa betas(), eigenvalues() */
const ComplexVectorType& alphas() const
{
eigen_assert(m_valuesOkay && "GeneralizedEigenSolver is not initialized.");
eigen_assert(info() == Success && "GeneralizedEigenSolver failed to compute alphas.");
return m_alphas;
}
@ -225,7 +226,7 @@ template<typename MatrixType_> class GeneralizedEigenSolver
* \sa alphas(), eigenvalues() */
const VectorType& betas() const
{
eigen_assert(m_valuesOkay && "GeneralizedEigenSolver is not initialized.");
eigen_assert(info() == Success && "GeneralizedEigenSolver failed to compute betas.");
return m_betas;
}
@ -256,7 +257,7 @@ template<typename MatrixType_> class GeneralizedEigenSolver
ComputationInfo info() const
{
eigen_assert(m_valuesOkay && "EigenSolver is not initialized.");
eigen_assert(m_isInitialized && "EigenSolver is not initialized.");
return m_realQZ.info();
}
@ -276,7 +277,8 @@ template<typename MatrixType_> class GeneralizedEigenSolver
EigenvectorsType m_eivec;
ComplexVectorType m_alphas;
VectorType m_betas;
bool m_valuesOkay, m_vectorsOkay;
bool m_computeEigenvectors;
bool m_isInitialized;
RealQZ<MatrixType> m_realQZ;
ComplexVectorType m_tmp;
};
@ -289,8 +291,6 @@ GeneralizedEigenSolver<MatrixType>::compute(const MatrixType& A, const MatrixTyp
using std::abs;
eigen_assert(A.cols() == A.rows() && B.cols() == A.rows() && B.cols() == B.rows());
Index size = A.cols();
m_valuesOkay = false;
m_vectorsOkay = false;
// Reduce to generalized real Schur form:
// A = Q S Z and B = Q T Z
m_realQZ.compute(A, B, computeEigenvectors);
@ -403,10 +403,9 @@ GeneralizedEigenSolver<MatrixType>::compute(const MatrixType& A, const MatrixTyp
i += 2;
}
}
m_valuesOkay = true;
m_vectorsOkay = computeEigenvectors;
}
m_computeEigenvectors = computeEigenvectors;
m_isInitialized = true;
return *this;
}

View File

@ -85,6 +85,42 @@ template<typename MatrixType> void generalized_eigensolver_real(const MatrixType
}
}
template<typename MatrixType>
void generalized_eigensolver_assert() {
GeneralizedEigenSolver<MatrixType> eig;
// all raise assert if uninitialized
VERIFY_RAISES_ASSERT(eig.info());
VERIFY_RAISES_ASSERT(eig.eigenvectors());
VERIFY_RAISES_ASSERT(eig.eigenvalues());
VERIFY_RAISES_ASSERT(eig.alphas());
VERIFY_RAISES_ASSERT(eig.betas());
// none raise assert after compute called
eig.compute(MatrixType::Random(20, 20), MatrixType::Random(20, 20));
VERIFY(eig.info() == Success);
eig.eigenvectors();
eig.eigenvalues();
eig.alphas();
eig.betas();
// eigenvectors() raises assert, if eigenvectors were not requested
eig.compute(MatrixType::Random(20, 20), MatrixType::Random(20, 20), false);
VERIFY(eig.info() == Success);
VERIFY_RAISES_ASSERT(eig.eigenvectors());
eig.eigenvalues();
eig.alphas();
eig.betas();
// all except info raise assert if realQZ did not converge
eig.setMaxIterations(0); // force real QZ to fail.
eig.compute(MatrixType::Random(20, 20), MatrixType::Random(20, 20));
VERIFY(eig.info() == NoConvergence);
VERIFY_RAISES_ASSERT(eig.eigenvectors());
VERIFY_RAISES_ASSERT(eig.eigenvalues());
VERIFY_RAISES_ASSERT(eig.alphas());
VERIFY_RAISES_ASSERT(eig.betas());
}
EIGEN_DECLARE_TEST(eigensolver_generalized_real)
{
for(int i = 0; i < g_repeat; i++) {
@ -98,6 +134,7 @@ EIGEN_DECLARE_TEST(eigensolver_generalized_real)
CALL_SUBTEST_2( generalized_eigensolver_real(MatrixXd(2,2)) );
CALL_SUBTEST_3( generalized_eigensolver_real(Matrix<double,1,1>()) );
CALL_SUBTEST_4( generalized_eigensolver_real(Matrix2d()) );
CALL_SUBTEST_5( generalized_eigensolver_assert<MatrixXd>() );
TEST_SET_BUT_UNUSED_VARIABLE(s)
}
}