|
| | LeastSquaresConjugateGradient () |
| |
| template<typename MatrixDerived > |
| | LeastSquaresConjugateGradient (const EigenBase< MatrixDerived > &A) |
| |
| | ~LeastSquaresConjugateGradient () |
| |
| template<typename Rhs , typename Dest > |
| void | _solve_with_guess_impl (const Rhs &b, Dest &x) const |
| |
| template<typename Rhs , typename Dest > |
| void | _solve_impl (const MatrixBase< Rhs > &b, Dest &x) const |
| |
| template<typename Rhs , typename DestDerived > |
| void | _solve_impl (const Rhs &b, SparseMatrixBase< DestDerived > &aDest) const |
| |
| LeastSquaresConjugateGradient< _MatrixType, _Preconditioner > & | analyzePattern (const EigenBase< MatrixDerived > &A) |
| |
| LeastSquaresConjugateGradient< _MatrixType, _Preconditioner > & | factorize (const EigenBase< MatrixDerived > &A) |
| |
| LeastSquaresConjugateGradient< _MatrixType, _Preconditioner > & | compute (const EigenBase< MatrixDerived > &A) |
| |
| Index | rows () const |
| |
| Index | cols () const |
| |
| RealScalar | tolerance () const |
| |
| LeastSquaresConjugateGradient< _MatrixType, _Preconditioner > & | setTolerance (const RealScalar &tolerance) |
| |
| Preconditioner & | preconditioner () |
| |
| const Preconditioner & | preconditioner () const |
| |
| Index | maxIterations () const |
| |
| LeastSquaresConjugateGradient< _MatrixType, _Preconditioner > & | setMaxIterations (Index maxIters) |
| |
| Index | iterations () const |
| |
| RealScalar | error () const |
| |
| const SolveWithGuess< LeastSquaresConjugateGradient< _MatrixType, _Preconditioner >, Rhs, Guess > | solveWithGuess (const MatrixBase< Rhs > &b, const Guess &x0) const |
| |
| ComputationInfo | info () const |
| |
| template<typename Rhs , typename Dest > |
| void | _solve_impl (const SparseMatrixBase< Rhs > &b, SparseMatrixBase< Dest > &dest) const |
| |
| LeastSquaresConjugateGradient< _MatrixType, _Preconditioner > & | derived () |
| |
| const LeastSquaresConjugateGradient< _MatrixType, _Preconditioner > & | derived () const |
| |
| template<typename Rhs > |
| const Solve< Derived, Rhs > | solve (const MatrixBase< Rhs > &b) const |
| |
| template<typename Rhs > |
| const Solve< Derived, Rhs > | solve (const SparseMatrixBase< Rhs > &b) const |
| |
template<typename _MatrixType, typename _Preconditioner>
class Eigen::LeastSquaresConjugateGradient< _MatrixType, _Preconditioner >
A conjugate gradient solver for sparse (or dense) least-square problems.
This class allows to solve for A x = b linear problems using an iterative conjugate gradient algorithm. The matrix A can be non symmetric and rectangular, but the matrix A' A should be positive-definite to guaranty stability. Otherwise, the SparseLU or SparseQR classes might be preferable. The matrix A and the vectors x and b can be either dense or sparse.
- Template Parameters
-
\implsparsesolverconcept
The maximal number of iterations and tolerance value can be controlled via the setMaxIterations() and setTolerance() methods. The defaults are the size of the problem for the maximal number of iterations and NumTraits<Scalar>::epsilon() for the tolerance.
This class can be used as the direct solver classes. Here is a typical usage example:
int m=1000000, n = 10000;
VectorXd x(n), b(m);
std::cout <<
"#iterations: " << lscg.
iterations() << std::endl;
std::cout <<
"estimated error: " << lscg.
error() << std::endl;
RealScalar error() const
Definition IterativeSolverBase.h:305
Derived & compute(const EigenBase< MatrixDerived > &A)
Definition IterativeSolverBase.h:238
Index iterations() const
Definition IterativeSolverBase.h:296
A conjugate gradient solver for sparse (or dense) least-square problems.
Definition LeastSquareConjugateGradient.h:150
A versatible sparse matrix representation.
Definition SparseMatrix.h:98
const Solve< Derived, Rhs > solve(const MatrixBase< Rhs > &b) const
Definition SparseSolverBase.h:88
By default the iterations start with x=0 as an initial guess of the solution. One can control the start using the solveWithGuess() method.
- See also
- class ConjugateGradient, SparseLU, SparseQR
template<typename _MatrixType , typename _Preconditioner >
template<typename Rhs , typename Dest >
186 {
189
190 for(
Index j=0; j<
b.cols(); ++j)
191 {
194
195 typename Dest::ColXpr xj(x,j);
197 }
198
201 }
Index maxIterations() const
Definition IterativeSolverBase.h:281
Preconditioner m_preconditioner
Definition IterativeSolverBase.h:381
RealScalar m_tolerance
Definition IterativeSolverBase.h:384
ComputationInfo m_info
Definition IterativeSolverBase.h:388
RealScalar m_error
Definition IterativeSolverBase.h:386
Index m_iterations
Definition IterativeSolverBase.h:387
bool m_isInitialized
Definition SparseSolverBase.h:119
const ActualMatrixType & matrix() const
Definition IterativeSolverBase.h:369
@ Success
Definition Constants.h:432
@ NoConvergence
Definition Constants.h:436
EIGEN_DONT_INLINE void least_square_conjugate_gradient(const MatrixType &mat, const Rhs &rhs, Dest &x, const Preconditioner &precond, Index &iters, typename Dest::RealScalar &tol_error)
Definition LeastSquareConjugateGradient.h:28
References Eigen::internal::least_square_conjugate_gradient(), Eigen::LeastSquaresConjugateGradient< _MatrixType, _Preconditioner >::m_error, Eigen::LeastSquaresConjugateGradient< _MatrixType, _Preconditioner >::m_info, Eigen::LeastSquaresConjugateGradient< _MatrixType, _Preconditioner >::m_isInitialized, Eigen::LeastSquaresConjugateGradient< _MatrixType, _Preconditioner >::m_iterations, Eigen::IterativeSolverBase< Derived >::m_preconditioner, Eigen::IterativeSolverBase< Derived >::m_tolerance, Eigen::LeastSquaresConjugateGradient< _MatrixType, _Preconditioner >::matrix(), Eigen::IterativeSolverBase< Derived >::maxIterations(), Eigen::NoConvergence, and Eigen::Success.
Referenced by Eigen::LeastSquaresConjugateGradient< _MatrixType, _Preconditioner >::_solve_impl().