Fix GeneralizedEigenSolver::info() and Asserts

(cherry picked from commit a7c1cac18bfef26ec61a73c1619ccf0f9b734745)
This commit is contained in:
Arthur 2022-08-25 22:05:04 +00:00 committed by Antonio Sanchez
parent d0e2b3e58d
commit 68f35d76b8
2 changed files with 53 additions and 17 deletions

View File

@ -119,8 +119,8 @@ template<typename _MatrixType> class GeneralizedEigenSolver
: m_eivec(), : m_eivec(),
m_alphas(), m_alphas(),
m_betas(), m_betas(),
m_valuesOkay(false), m_computeEigenvectors(false),
m_vectorsOkay(false), m_isInitialized(false),
m_realQZ() m_realQZ()
{} {}
@ -134,8 +134,8 @@ template<typename _MatrixType> class GeneralizedEigenSolver
: m_eivec(size, size), : m_eivec(size, size),
m_alphas(size), m_alphas(size),
m_betas(size), m_betas(size),
m_valuesOkay(false), m_computeEigenvectors(false),
m_vectorsOkay(false), m_isInitialized(false),
m_realQZ(size), m_realQZ(size),
m_tmp(size) m_tmp(size)
{} {}
@ -156,8 +156,8 @@ template<typename _MatrixType> class GeneralizedEigenSolver
: m_eivec(A.rows(), A.cols()), : m_eivec(A.rows(), A.cols()),
m_alphas(A.cols()), m_alphas(A.cols()),
m_betas(A.cols()), m_betas(A.cols()),
m_valuesOkay(false), m_computeEigenvectors(false),
m_vectorsOkay(false), m_isInitialized(false),
m_realQZ(A.cols()), m_realQZ(A.cols()),
m_tmp(A.cols()) m_tmp(A.cols())
{ {
@ -177,7 +177,8 @@ template<typename _MatrixType> class GeneralizedEigenSolver
* \sa eigenvalues() * \sa eigenvalues()
*/ */
EigenvectorsType eigenvectors() const { 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; return m_eivec;
} }
@ -201,7 +202,7 @@ template<typename _MatrixType> class GeneralizedEigenSolver
*/ */
EigenvalueType eigenvalues() const 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); return EigenvalueType(m_alphas,m_betas);
} }
@ -212,7 +213,7 @@ template<typename _MatrixType> class GeneralizedEigenSolver
* \sa betas(), eigenvalues() */ * \sa betas(), eigenvalues() */
ComplexVectorType alphas() const ComplexVectorType alphas() const
{ {
eigen_assert(m_valuesOkay && "GeneralizedEigenSolver is not initialized."); eigen_assert(info() == Success && "GeneralizedEigenSolver failed to compute alphas.");
return m_alphas; return m_alphas;
} }
@ -223,7 +224,7 @@ template<typename _MatrixType> class GeneralizedEigenSolver
* \sa alphas(), eigenvalues() */ * \sa alphas(), eigenvalues() */
VectorType betas() const VectorType betas() const
{ {
eigen_assert(m_valuesOkay && "GeneralizedEigenSolver is not initialized."); eigen_assert(info() == Success && "GeneralizedEigenSolver failed to compute betas.");
return m_betas; return m_betas;
} }
@ -254,7 +255,7 @@ template<typename _MatrixType> class GeneralizedEigenSolver
ComputationInfo info() const ComputationInfo info() const
{ {
eigen_assert(m_valuesOkay && "EigenSolver is not initialized."); eigen_assert(m_isInitialized && "EigenSolver is not initialized.");
return m_realQZ.info(); return m_realQZ.info();
} }
@ -277,7 +278,8 @@ template<typename _MatrixType> class GeneralizedEigenSolver
EigenvectorsType m_eivec; EigenvectorsType m_eivec;
ComplexVectorType m_alphas; ComplexVectorType m_alphas;
VectorType m_betas; VectorType m_betas;
bool m_valuesOkay, m_vectorsOkay; bool m_computeEigenvectors;
bool m_isInitialized;
RealQZ<MatrixType> m_realQZ; RealQZ<MatrixType> m_realQZ;
ComplexVectorType m_tmp; ComplexVectorType m_tmp;
}; };
@ -292,8 +294,6 @@ GeneralizedEigenSolver<MatrixType>::compute(const MatrixType& A, const MatrixTyp
using std::abs; using std::abs;
eigen_assert(A.cols() == A.rows() && B.cols() == A.rows() && B.cols() == B.rows()); eigen_assert(A.cols() == A.rows() && B.cols() == A.rows() && B.cols() == B.rows());
Index size = A.cols(); Index size = A.cols();
m_valuesOkay = false;
m_vectorsOkay = false;
// Reduce to generalized real Schur form: // Reduce to generalized real Schur form:
// A = Q S Z and B = Q T Z // A = Q S Z and B = Q T Z
m_realQZ.compute(A, B, computeEigenvectors); m_realQZ.compute(A, B, computeEigenvectors);
@ -406,10 +406,9 @@ GeneralizedEigenSolver<MatrixType>::compute(const MatrixType& A, const MatrixTyp
i += 2; i += 2;
} }
} }
m_valuesOkay = true;
m_vectorsOkay = computeEigenvectors;
} }
m_computeEigenvectors = computeEigenvectors;
m_isInitialized = true;
return *this; 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) EIGEN_DECLARE_TEST(eigensolver_generalized_real)
{ {
for(int i = 0; i < g_repeat; i++) { 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_2( generalized_eigensolver_real(MatrixXd(2,2)) );
CALL_SUBTEST_3( generalized_eigensolver_real(Matrix<double,1,1>()) ); CALL_SUBTEST_3( generalized_eigensolver_real(Matrix<double,1,1>()) );
CALL_SUBTEST_4( generalized_eigensolver_real(Matrix2d()) ); CALL_SUBTEST_4( generalized_eigensolver_real(Matrix2d()) );
CALL_SUBTEST_5( generalized_eigensolver_assert<MatrixXd>() );
TEST_SET_BUT_UNUSED_VARIABLE(s) TEST_SET_BUT_UNUSED_VARIABLE(s)
} }
} }