$treeview $search $mathjax
Eigen
3.2.5
$projectbrief
|
$projectbrief
|
$searchbox |
00001 // This file is part of Eigen, a lightweight C++ template library 00002 // for linear algebra. 00003 // 00004 // Copyright (C) 2010 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_SPARSE_SOLVE_H 00011 #define EIGEN_SPARSE_SOLVE_H 00012 00013 namespace Eigen { 00014 00015 namespace internal { 00016 00017 template<typename _DecompositionType, typename Rhs> struct sparse_solve_retval_base; 00018 template<typename _DecompositionType, typename Rhs> struct sparse_solve_retval; 00019 00020 template<typename DecompositionType, typename Rhs> 00021 struct traits<sparse_solve_retval_base<DecompositionType, Rhs> > 00022 { 00023 typedef typename DecompositionType::MatrixType MatrixType; 00024 typedef SparseMatrix<typename Rhs::Scalar, Rhs::Options, typename Rhs::Index> ReturnType; 00025 }; 00026 00027 template<typename _DecompositionType, typename Rhs> struct sparse_solve_retval_base 00028 : public ReturnByValue<sparse_solve_retval_base<_DecompositionType, Rhs> > 00029 { 00030 typedef typename remove_all<typename Rhs::Nested>::type RhsNestedCleaned; 00031 typedef _DecompositionType DecompositionType; 00032 typedef ReturnByValue<sparse_solve_retval_base> Base; 00033 typedef typename Base::Index Index; 00034 00035 sparse_solve_retval_base(const DecompositionType& dec, const Rhs& rhs) 00036 : m_dec(dec), m_rhs(rhs) 00037 {} 00038 00039 inline Index rows() const { return m_dec.cols(); } 00040 inline Index cols() const { return m_rhs.cols(); } 00041 inline const DecompositionType& dec() const { return m_dec; } 00042 inline const RhsNestedCleaned& rhs() const { return m_rhs; } 00043 00044 template<typename Dest> inline void evalTo(Dest& dst) const 00045 { 00046 static_cast<const sparse_solve_retval<DecompositionType,Rhs>*>(this)->evalTo(dst); 00047 } 00048 00049 protected: 00050 template<typename DestScalar, int DestOptions, typename DestIndex> 00051 inline void defaultEvalTo(SparseMatrix<DestScalar,DestOptions,DestIndex>& dst) const 00052 { 00053 // we process the sparse rhs per block of NbColsAtOnce columns temporarily stored into a dense matrix. 00054 static const int NbColsAtOnce = 4; 00055 int rhsCols = m_rhs.cols(); 00056 int size = m_rhs.rows(); 00057 Eigen::Matrix<DestScalar,Dynamic,Dynamic> tmp(size,rhsCols); 00058 Eigen::Matrix<DestScalar,Dynamic,Dynamic> tmpX(size,rhsCols); 00059 for(int k=0; k<rhsCols; k+=NbColsAtOnce) 00060 { 00061 int actualCols = std::min<int>(rhsCols-k, NbColsAtOnce); 00062 tmp.leftCols(actualCols) = m_rhs.middleCols(k,actualCols); 00063 tmpX.leftCols(actualCols) = m_dec.solve(tmp.leftCols(actualCols)); 00064 dst.middleCols(k,actualCols) = tmpX.leftCols(actualCols).sparseView(); 00065 } 00066 } 00067 const DecompositionType& m_dec; 00068 typename Rhs::Nested m_rhs; 00069 }; 00070 00071 #define EIGEN_MAKE_SPARSE_SOLVE_HELPERS(DecompositionType,Rhs) \ 00072 typedef typename DecompositionType::MatrixType MatrixType; \ 00073 typedef typename MatrixType::Scalar Scalar; \ 00074 typedef typename MatrixType::RealScalar RealScalar; \ 00075 typedef typename MatrixType::Index Index; \ 00076 typedef Eigen::internal::sparse_solve_retval_base<DecompositionType,Rhs> Base; \ 00077 using Base::dec; \ 00078 using Base::rhs; \ 00079 using Base::rows; \ 00080 using Base::cols; \ 00081 sparse_solve_retval(const DecompositionType& dec, const Rhs& rhs) \ 00082 : Base(dec, rhs) {} 00083 00084 00085 00086 template<typename DecompositionType, typename Rhs, typename Guess> struct solve_retval_with_guess; 00087 00088 template<typename DecompositionType, typename Rhs, typename Guess> 00089 struct traits<solve_retval_with_guess<DecompositionType, Rhs, Guess> > 00090 { 00091 typedef typename DecompositionType::MatrixType MatrixType; 00092 typedef Matrix<typename Rhs::Scalar, 00093 MatrixType::ColsAtCompileTime, 00094 Rhs::ColsAtCompileTime, 00095 Rhs::PlainObject::Options, 00096 MatrixType::MaxColsAtCompileTime, 00097 Rhs::MaxColsAtCompileTime> ReturnType; 00098 }; 00099 00100 template<typename DecompositionType, typename Rhs, typename Guess> struct solve_retval_with_guess 00101 : public ReturnByValue<solve_retval_with_guess<DecompositionType, Rhs, Guess> > 00102 { 00103 typedef typename DecompositionType::Index Index; 00104 00105 solve_retval_with_guess(const DecompositionType& dec, const Rhs& rhs, const Guess& guess) 00106 : m_dec(dec), m_rhs(rhs), m_guess(guess) 00107 {} 00108 00109 inline Index rows() const { return m_dec.cols(); } 00110 inline Index cols() const { return m_rhs.cols(); } 00111 00112 template<typename Dest> inline void evalTo(Dest& dst) const 00113 { 00114 dst = m_guess; 00115 m_dec._solveWithGuess(m_rhs,dst); 00116 } 00117 00118 protected: 00119 const DecompositionType& m_dec; 00120 const typename Rhs::Nested m_rhs; 00121 const typename Guess::Nested m_guess; 00122 }; 00123 00124 } // namepsace internal 00125 00126 } // end namespace Eigen 00127 00128 #endif // EIGEN_SPARSE_SOLVE_H