Eigen  3.3.3
SparseSolverBase.h
00001 // This file is part of Eigen, a lightweight C++ template library
00002 // for linear algebra.
00003 //
00004 // Copyright (C) 2014 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_SPARSESOLVERBASE_H
00011 #define EIGEN_SPARSESOLVERBASE_H
00012 
00013 namespace Eigen { 
00014 
00015 namespace internal {
00016 
00021 template<typename Decomposition, typename Rhs, typename Dest>
00022 typename enable_if<Rhs::ColsAtCompileTime!=1 && Dest::ColsAtCompileTime!=1>::type
00023 solve_sparse_through_dense_panels(const Decomposition &dec, const Rhs& rhs, Dest &dest)
00024 {
00025   EIGEN_STATIC_ASSERT((Dest::Flags&RowMajorBit)==0,THIS_METHOD_IS_ONLY_FOR_COLUMN_MAJOR_MATRICES);
00026   typedef typename Dest::Scalar DestScalar;
00027   // we process the sparse rhs per block of NbColsAtOnce columns temporarily stored into a dense matrix.
00028   static const Index NbColsAtOnce = 4;
00029   Index rhsCols = rhs.cols();
00030   Index size = rhs.rows();
00031   // the temporary matrices do not need more columns than NbColsAtOnce:
00032   Index tmpCols = (std::min)(rhsCols, NbColsAtOnce); 
00033   Eigen::Matrix<DestScalar,Dynamic,Dynamic> tmp(size,tmpCols);
00034   Eigen::Matrix<DestScalar,Dynamic,Dynamic> tmpX(size,tmpCols);
00035   for(Index k=0; k<rhsCols; k+=NbColsAtOnce)
00036   {
00037     Index actualCols = std::min<Index>(rhsCols-k, NbColsAtOnce);
00038     tmp.leftCols(actualCols) = rhs.middleCols(k,actualCols);
00039     tmpX.leftCols(actualCols) = dec.solve(tmp.leftCols(actualCols));
00040     dest.middleCols(k,actualCols) = tmpX.leftCols(actualCols).sparseView();
00041   }
00042 }
00043 
00044 // Overload for vector as rhs
00045 template<typename Decomposition, typename Rhs, typename Dest>
00046 typename enable_if<Rhs::ColsAtCompileTime==1 || Dest::ColsAtCompileTime==1>::type
00047 solve_sparse_through_dense_panels(const Decomposition &dec, const Rhs& rhs, Dest &dest)
00048 {
00049   typedef typename Dest::Scalar DestScalar;
00050   Index size = rhs.rows();
00051   Eigen::Matrix<DestScalar,Dynamic,1> rhs_dense(rhs);
00052   Eigen::Matrix<DestScalar,Dynamic,1> dest_dense(size);
00053   dest_dense = dec.solve(rhs_dense);
00054   dest = dest_dense.sparseView();
00055 }
00056 
00057 } // end namespace internal
00058 
00066 template<typename Derived>
00067 class SparseSolverBase : internal::noncopyable
00068 {
00069   public:
00070 
00072     SparseSolverBase()
00073       : m_isInitialized(false)
00074     {}
00075 
00076     ~SparseSolverBase()
00077     {}
00078 
00079     Derived& derived() { return *static_cast<Derived*>(this); }
00080     const Derived& derived() const { return *static_cast<const Derived*>(this); }
00081     
00086     template<typename Rhs>
00087     inline const Solve<Derived, Rhs>
00088     solve(const MatrixBase<Rhs>& b) const
00089     {
00090       eigen_assert(m_isInitialized && "Solver is not initialized.");
00091       eigen_assert(derived().rows()==b.rows() && "solve(): invalid number of rows of the right hand side matrix b");
00092       return Solve<Derived, Rhs>(derived(), b.derived());
00093     }
00094     
00099     template<typename Rhs>
00100     inline const Solve<Derived, Rhs>
00101     solve(const SparseMatrixBase<Rhs>& b) const
00102     {
00103       eigen_assert(m_isInitialized && "Solver is not initialized.");
00104       eigen_assert(derived().rows()==b.rows() && "solve(): invalid number of rows of the right hand side matrix b");
00105       return Solve<Derived, Rhs>(derived(), b.derived());
00106     }
00107     
00108     #ifndef EIGEN_PARSED_BY_DOXYGEN
00109 
00110     template<typename Rhs,typename Dest>
00111     void _solve_impl(const SparseMatrixBase<Rhs> &b, SparseMatrixBase<Dest> &dest) const
00112     {
00113       internal::solve_sparse_through_dense_panels(derived(), b.derived(), dest.derived());
00114     }
00115     #endif // EIGEN_PARSED_BY_DOXYGEN
00116 
00117   protected:
00118     
00119     mutable bool m_isInitialized;
00120 };
00121 
00122 } // end namespace Eigen
00123 
00124 #endif // EIGEN_SPARSESOLVERBASE_H
 All Classes Functions Variables Typedefs Enumerations Enumerator Friends