Fix misuse of dummy_precesion in eigenvalues solvers

This commit is contained in:
Gael Guennebaud 2016-07-23 17:52:31 +02:00
parent 72744d93ef
commit 1b0353c659
2 changed files with 6 additions and 3 deletions

View File

@ -324,11 +324,12 @@ template<typename MatrixType>
MatrixType EigenSolver<MatrixType>::pseudoEigenvalueMatrix() const MatrixType EigenSolver<MatrixType>::pseudoEigenvalueMatrix() const
{ {
eigen_assert(m_isInitialized && "EigenSolver is not initialized."); eigen_assert(m_isInitialized && "EigenSolver is not initialized.");
const RealScalar precision = RealScalar(2)*NumTraits<RealScalar>::epsilon();
Index n = m_eivalues.rows(); Index n = m_eivalues.rows();
MatrixType matD = MatrixType::Zero(n,n); MatrixType matD = MatrixType::Zero(n,n);
for (Index i=0; i<n; ++i) for (Index i=0; i<n; ++i)
{ {
if (internal::isMuchSmallerThan(numext::imag(m_eivalues.coeff(i)), numext::real(m_eivalues.coeff(i)))) if (internal::isMuchSmallerThan(numext::imag(m_eivalues.coeff(i)), numext::real(m_eivalues.coeff(i)), precision))
matD.coeffRef(i,i) = numext::real(m_eivalues.coeff(i)); matD.coeffRef(i,i) = numext::real(m_eivalues.coeff(i));
else else
{ {
@ -345,11 +346,12 @@ typename EigenSolver<MatrixType>::EigenvectorsType EigenSolver<MatrixType>::eige
{ {
eigen_assert(m_isInitialized && "EigenSolver is not initialized."); eigen_assert(m_isInitialized && "EigenSolver is not initialized.");
eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues."); eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
const RealScalar precision = RealScalar(2)*NumTraits<RealScalar>::epsilon();
Index n = m_eivec.cols(); Index n = m_eivec.cols();
EigenvectorsType matV(n,n); EigenvectorsType matV(n,n);
for (Index j=0; j<n; ++j) for (Index j=0; j<n; ++j)
{ {
if (internal::isMuchSmallerThan(numext::imag(m_eivalues.coeff(j)), numext::real(m_eivalues.coeff(j))) || j+1==n) if (internal::isMuchSmallerThan(numext::imag(m_eivalues.coeff(j)), numext::real(m_eivalues.coeff(j)), precision) || j+1==n)
{ {
// we have a real eigen value // we have a real eigen value
matV.col(j) = m_eivec.col(j).template cast<ComplexScalar>(); matV.col(j) = m_eivec.col(j).template cast<ComplexScalar>();

View File

@ -492,11 +492,12 @@ ComputationInfo computeFromTridiagonal_impl(DiagType& diag, SubDiagType& subdiag
typedef typename DiagType::RealScalar RealScalar; typedef typename DiagType::RealScalar RealScalar;
const RealScalar considerAsZero = (std::numeric_limits<RealScalar>::min)(); const RealScalar considerAsZero = (std::numeric_limits<RealScalar>::min)();
const RealScalar precision = RealScalar(2)*NumTraits<RealScalar>::epsilon();
while (end>0) while (end>0)
{ {
for (Index i = start; i<end; ++i) for (Index i = start; i<end; ++i)
if (internal::isMuchSmallerThan(abs(subdiag[i]),(abs(diag[i])+abs(diag[i+1]))) || abs(subdiag[i]) <= considerAsZero) if (internal::isMuchSmallerThan(abs(subdiag[i]),(abs(diag[i])+abs(diag[i+1])),precision) || abs(subdiag[i]) <= considerAsZero)
subdiag[i] = 0; subdiag[i] = 0;
// find the largest unreduced block // find the largest unreduced block