IterativeSolverBase.h
Go to the documentation of this file.
00001 // This file is part of Eigen, a lightweight C++ template library
00002 // for linear algebra.
00003 //
00004 // Copyright (C) 2011 Gael Guennebaud <gael.guennebaud@inria.fr>
00005 //
00006 // This Source Code Form is subject to the terms of the Mozilla
00007 // Public License v. 2.0. If a copy of the MPL was not distributed
00008 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
00009 
00010 #ifndef EIGEN_ITERATIVE_SOLVER_BASE_H
00011 #define EIGEN_ITERATIVE_SOLVER_BASE_H
00012 
00013 namespace Eigen { 
00014 
00020 template< typename Derived>
00021 class IterativeSolverBase : internal::noncopyable
00022 {
00023 public:
00024   typedef typename internal::traits<Derived>::MatrixType MatrixType;
00025   typedef typename internal::traits<Derived>::Preconditioner Preconditioner;
00026   typedef typename MatrixType::Scalar Scalar;
00027   typedef typename MatrixType::Index Index;
00028   typedef typename MatrixType::RealScalar RealScalar;
00029 
00030 public:
00031 
00032   Derived& derived() { return *static_cast<Derived*>(this); }
00033   const Derived& derived() const { return *static_cast<const Derived*>(this); }
00034 
00036   IterativeSolverBase()
00037     : mp_matrix(0)
00038   {
00039     init();
00040   }
00041 
00052   IterativeSolverBase(const MatrixType& A)
00053   {
00054     init();
00055     compute(A);
00056   }
00057 
00058   ~IterativeSolverBase() {}
00059   
00065   Derived& analyzePattern(const MatrixType& A)
00066   {
00067     m_preconditioner.analyzePattern(A);
00068     m_isInitialized = true;
00069     m_analysisIsOk = true;
00070     m_info = Success;
00071     return derived();
00072   }
00073   
00083   Derived& factorize(const MatrixType& A)
00084   {
00085     eigen_assert(m_analysisIsOk && "You must first call analyzePattern()"); 
00086     mp_matrix = &A;
00087     m_preconditioner.factorize(A);
00088     m_factorizationIsOk = true;
00089     m_info = Success;
00090     return derived();
00091   }
00092 
00103   Derived& compute(const MatrixType& A)
00104   {
00105     mp_matrix = &A;
00106     m_preconditioner.compute(A);
00107     m_isInitialized = true;
00108     m_analysisIsOk = true;
00109     m_factorizationIsOk = true;
00110     m_info = Success;
00111     return derived();
00112   }
00113 
00115   Index rows() const { return mp_matrix ? mp_matrix->rows() : 0; }
00117   Index cols() const { return mp_matrix ? mp_matrix->cols() : 0; }
00118 
00120   RealScalar tolerance() const { return m_tolerance; }
00121   
00123   Derived& setTolerance(const RealScalar& tolerance)
00124   {
00125     m_tolerance = tolerance;
00126     return derived();
00127   }
00128 
00130   Preconditioner& preconditioner() { return m_preconditioner; }
00131   
00133   const Preconditioner& preconditioner() const { return m_preconditioner; }
00134 
00136   int maxIterations() const
00137   {
00138     return (mp_matrix && m_maxIterations<0) ? mp_matrix->cols() : m_maxIterations;
00139   }
00140   
00142   Derived& setMaxIterations(int maxIters)
00143   {
00144     m_maxIterations = maxIters;
00145     return derived();
00146   }
00147 
00149   int iterations() const
00150   {
00151     eigen_assert(m_isInitialized && "ConjugateGradient is not initialized.");
00152     return m_iterations;
00153   }
00154 
00156   RealScalar error() const
00157   {
00158     eigen_assert(m_isInitialized && "ConjugateGradient is not initialized.");
00159     return m_error;
00160   }
00161 
00166   template<typename Rhs> inline const internal::solve_retval<Derived, Rhs>
00167   solve(const MatrixBase<Rhs>& b) const
00168   {
00169     eigen_assert(m_isInitialized && "IterativeSolverBase is not initialized.");
00170     eigen_assert(rows()==b.rows()
00171               && "IterativeSolverBase::solve(): invalid number of rows of the right hand side matrix b");
00172     return internal::solve_retval<Derived, Rhs>(derived(), b.derived());
00173   }
00174   
00179   template<typename Rhs>
00180   inline const internal::sparse_solve_retval<IterativeSolverBase, Rhs>
00181   solve(const SparseMatrixBase<Rhs>& b) const
00182   {
00183     eigen_assert(m_isInitialized && "IterativeSolverBase is not initialized.");
00184     eigen_assert(rows()==b.rows()
00185               && "IterativeSolverBase::solve(): invalid number of rows of the right hand side matrix b");
00186     return internal::sparse_solve_retval<IterativeSolverBase, Rhs>(*this, b.derived());
00187   }
00188 
00190   ComputationInfo info() const
00191   {
00192     eigen_assert(m_isInitialized && "IterativeSolverBase is not initialized.");
00193     return m_info;
00194   }
00195   
00197   template<typename Rhs, typename DestScalar, int DestOptions, typename DestIndex>
00198   void _solve_sparse(const Rhs& b, SparseMatrix<DestScalar,DestOptions,DestIndex> &dest) const
00199   {
00200     eigen_assert(rows()==b.rows());
00201     
00202     int rhsCols = b.cols();
00203     int size = b.rows();
00204     Eigen::Matrix<DestScalar,Dynamic,1> tb(size);
00205     Eigen::Matrix<DestScalar,Dynamic,1> tx(size);
00206     for(int k=0; k<rhsCols; ++k)
00207     {
00208       tb = b.col(k);
00209       tx = derived().solve(tb);
00210       dest.col(k) = tx.sparseView(0);
00211     }
00212   }
00213 
00214 protected:
00215   void init()
00216   {
00217     m_isInitialized = false;
00218     m_analysisIsOk = false;
00219     m_factorizationIsOk = false;
00220     m_maxIterations = -1;
00221     m_tolerance = NumTraits<Scalar>::epsilon();
00222   }
00223   const MatrixType* mp_matrix;
00224   Preconditioner m_preconditioner;
00225 
00226   int m_maxIterations;
00227   RealScalar m_tolerance;
00228   
00229   mutable RealScalar m_error;
00230   mutable int m_iterations;
00231   mutable ComputationInfo m_info;
00232   mutable bool m_isInitialized, m_analysisIsOk, m_factorizationIsOk;
00233 };
00234 
00235 namespace internal {
00236  
00237 template<typename Derived, typename Rhs>
00238 struct sparse_solve_retval<IterativeSolverBase<Derived>, Rhs>
00239   : sparse_solve_retval_base<IterativeSolverBase<Derived>, Rhs>
00240 {
00241   typedef IterativeSolverBase<Derived> Dec;
00242   EIGEN_MAKE_SPARSE_SOLVE_HELPERS(Dec,Rhs)
00243 
00244   template<typename Dest> void evalTo(Dest& dst) const
00245   {
00246     dec().derived()._solve_sparse(rhs(),dst);
00247   }
00248 };
00249 
00250 } // end namespace internal
00251 
00252 } // end namespace Eigen
00253 
00254 #endif // EIGEN_ITERATIVE_SOLVER_BASE_H


acado
Author(s): Milan Vukov, Rien Quirynen
autogenerated on Thu Aug 27 2015 11:58:42