aboutsummaryrefslogtreecommitdiffhomepage
diff options
context:
space:
mode:
authorGravatar Gael Guennebaud <g.gael@free.fr>2015-08-04 16:13:34 +0200
committerGravatar Gael Guennebaud <g.gael@free.fr>2015-08-04 16:13:34 +0200
commit506964fc29dc8ad2f76bed5cdaacbc1d39785912 (patch)
treee217460c3fceb0a5e3f352c273e5d5566dcb117e
parentdb0f5c9d90f76f21f231063d9cbd8ac7f274a011 (diff)
Propagate precondition info to the iterative solver.
-rw-r--r--Eigen/src/IterativeLinearSolvers/BasicPreconditioners.h6
-rw-r--r--Eigen/src/IterativeLinearSolvers/IncompleteLUT.h3
-rw-r--r--Eigen/src/IterativeLinearSolvers/IterativeSolverBase.h6
3 files changed, 10 insertions, 5 deletions
diff --git a/Eigen/src/IterativeLinearSolvers/BasicPreconditioners.h b/Eigen/src/IterativeLinearSolvers/BasicPreconditioners.h
index 3710a8209..ff7f08c1c 100644
--- a/Eigen/src/IterativeLinearSolvers/BasicPreconditioners.h
+++ b/Eigen/src/IterativeLinearSolvers/BasicPreconditioners.h
@@ -95,6 +95,8 @@ class DiagonalPreconditioner
&& "DiagonalPreconditioner::solve(): invalid number of rows of the right hand side matrix b");
return Solve<DiagonalPreconditioner, Rhs>(*this, b.derived());
}
+
+ ComputationInfo info() { return Success; }
protected:
Vector m_invdiag;
@@ -161,6 +163,8 @@ class LeastSquareDiagonalPreconditioner : public DiagonalPreconditioner<_Scalar>
{
return factorize(mat);
}
+
+ ComputationInfo info() { return Success; }
protected:
};
@@ -190,6 +194,8 @@ class IdentityPreconditioner
template<typename Rhs>
inline const Rhs& solve(const Rhs& b) const { return b; }
+
+ ComputationInfo info() { return Success; }
};
} // end namespace Eigen
diff --git a/Eigen/src/IterativeLinearSolvers/IncompleteLUT.h b/Eigen/src/IterativeLinearSolvers/IncompleteLUT.h
index 102e01f76..b644163f1 100644
--- a/Eigen/src/IterativeLinearSolvers/IncompleteLUT.h
+++ b/Eigen/src/IterativeLinearSolvers/IncompleteLUT.h
@@ -232,7 +232,7 @@ void IncompleteLUT<Scalar,StorageIndex>::analyzePattern(const _MatrixType& amat)
m_analysisIsOk = true;
m_factorizationIsOk = false;
- m_isInitialized = false;
+ m_isInitialized = true;
}
template <typename Scalar, typename StorageIndex>
@@ -441,7 +441,6 @@ void IncompleteLUT<Scalar,StorageIndex>::factorize(const _MatrixType& amat)
m_lu.makeCompressed();
m_factorizationIsOk = true;
- m_isInitialized = m_factorizationIsOk;
m_info = Success;
}
diff --git a/Eigen/src/IterativeLinearSolvers/IterativeSolverBase.h b/Eigen/src/IterativeLinearSolvers/IterativeSolverBase.h
index 22c602a91..a7cb30cb9 100644
--- a/Eigen/src/IterativeLinearSolvers/IterativeSolverBase.h
+++ b/Eigen/src/IterativeLinearSolvers/IterativeSolverBase.h
@@ -74,7 +74,7 @@ public:
m_preconditioner.analyzePattern(mp_matrix);
m_isInitialized = true;
m_analysisIsOk = true;
- m_info = Success;
+ m_info = m_preconditioner.info();
return derived();
}
@@ -94,7 +94,7 @@ public:
grab(A.derived());
m_preconditioner.factorize(mp_matrix);
m_factorizationIsOk = true;
- m_info = Success;
+ m_info = m_preconditioner.info();
return derived();
}
@@ -116,7 +116,7 @@ public:
m_isInitialized = true;
m_analysisIsOk = true;
m_factorizationIsOk = true;
- m_info = Success;
+ m_info = m_preconditioner.info();
return derived();
}