diff options
| -rw-r--r-- | Eigen/src/Eigenvalues/Tridiagonalization.h | 11 | ||||
| -rw-r--r-- | test/eigensolver_selfadjoint.cpp | 11 |
2 files changed, 14 insertions, 8 deletions
diff --git a/Eigen/src/Eigenvalues/Tridiagonalization.h b/Eigen/src/Eigenvalues/Tridiagonalization.h index f20528f86..bedd1cb34 100644 --- a/Eigen/src/Eigenvalues/Tridiagonalization.h +++ b/Eigen/src/Eigenvalues/Tridiagonalization.h @@ -91,10 +91,8 @@ template<typename _MatrixType> class Tridiagonalization >::type DiagonalReturnType; typedef typename internal::conditional<NumTraits<Scalar>::IsComplex, - typename internal::add_const_on_value_type<typename Diagonal< - Block<const MatrixType,SizeMinusOne,SizeMinusOne> >::RealReturnType>::type, - const Diagonal< - Block<const MatrixType,SizeMinusOne,SizeMinusOne> > + typename internal::add_const_on_value_type<typename Diagonal<const MatrixType, -1>::RealReturnType>::type, + const Diagonal<const MatrixType, -1> >::type SubDiagonalReturnType; /** \brief Return type of matrixQ() */ @@ -307,7 +305,7 @@ typename Tridiagonalization<MatrixType>::DiagonalReturnType Tridiagonalization<MatrixType>::diagonal() const { eigen_assert(m_isInitialized && "Tridiagonalization is not initialized."); - return m_matrix.diagonal(); + return m_matrix.diagonal().real(); } template<typename MatrixType> @@ -315,8 +313,7 @@ typename Tridiagonalization<MatrixType>::SubDiagonalReturnType Tridiagonalization<MatrixType>::subDiagonal() const { eigen_assert(m_isInitialized && "Tridiagonalization is not initialized."); - Index n = m_matrix.rows(); - return Block<const MatrixType,SizeMinusOne,SizeMinusOne>(m_matrix, 1, 0, n-1,n-1).diagonal(); + return m_matrix.template diagonal<-1>().real(); } namespace internal { diff --git a/test/eigensolver_selfadjoint.cpp b/test/eigensolver_selfadjoint.cpp index 3851f9df2..41b598795 100644 --- a/test/eigensolver_selfadjoint.cpp +++ b/test/eigensolver_selfadjoint.cpp @@ -111,8 +111,17 @@ template<typename MatrixType> void selfadjointeigensolver(const MatrixType& m) // test Tridiagonalization's methods Tridiagonalization<MatrixType> tridiag(symmC); - // FIXME tridiag.matrixQ().adjoint() does not work + VERIFY_IS_APPROX(tridiag.diagonal(), tridiag.matrixT().template diagonal()); + VERIFY_IS_APPROX(tridiag.subDiagonal(), tridiag.matrixT().template diagonal<-1>()); + MatrixType T = tridiag.matrixT(); + if(rows>1 && cols>1) { + // FIXME check that upper and lower part are 0: + //VERIFY(T.topRightCorner(rows-2, cols-2).template triangularView<Upper>().isZero()); + } + VERIFY_IS_APPROX(tridiag.diagonal(), T.diagonal().real()); + VERIFY_IS_APPROX(tridiag.subDiagonal(), T.template diagonal<1>().real()); VERIFY_IS_APPROX(MatrixType(symmC.template selfadjointView<Lower>()), tridiag.matrixQ() * tridiag.matrixT().eval() * MatrixType(tridiag.matrixQ()).adjoint()); + VERIFY_IS_APPROX(MatrixType(symmC.template selfadjointView<Lower>()), tridiag.matrixQ() * tridiag.matrixT() * tridiag.matrixQ().adjoint()); // Test computation of eigenvalues from tridiagonal matrix if(rows > 1) |