SparseSolve.h
Go to the documentation of this file.
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2010 Gael Guennebaud <gael.guennebaud@inria.fr>
5 //
6 // Eigen is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 3 of the License, or (at your option) any later version.
10 //
11 // Alternatively, you can redistribute it and/or
12 // modify it under the terms of the GNU General Public License as
13 // published by the Free Software Foundation; either version 2 of
14 // the License, or (at your option) any later version.
15 //
16 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
17 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
19 // GNU General Public License for more details.
20 //
21 // You should have received a copy of the GNU Lesser General Public
22 // License and a copy of the GNU General Public License along with
23 // Eigen. If not, see <http://www.gnu.org/licenses/>.
24 
25 #ifndef EIGEN_SPARSE_SOLVE_H
26 #define EIGEN_SPARSE_SOLVE_H
27 
28 namespace Eigen {
29 
30 namespace internal {
31 
32 template<typename _DecompositionType, typename Rhs> struct sparse_solve_retval_base;
33 template<typename _DecompositionType, typename Rhs> struct sparse_solve_retval;
34 
35 template<typename DecompositionType, typename Rhs>
36 struct traits<sparse_solve_retval_base<DecompositionType, Rhs> >
37 {
38  typedef typename DecompositionType::MatrixType MatrixType;
39  typedef SparseMatrix<typename Rhs::Scalar, Rhs::Options, typename Rhs::Index> ReturnType;
40 };
41 
42 template<typename _DecompositionType, typename Rhs> struct sparse_solve_retval_base
43  : public ReturnByValue<sparse_solve_retval_base<_DecompositionType, Rhs> >
44 {
45  typedef typename remove_all<typename Rhs::Nested>::type RhsNestedCleaned;
46  typedef _DecompositionType DecompositionType;
47  typedef ReturnByValue<sparse_solve_retval_base> Base;
48  typedef typename Base::Index Index;
49 
50  sparse_solve_retval_base(const DecompositionType& dec, const Rhs& rhs)
51  : m_dec(dec), m_rhs(rhs)
52  {}
53 
54  inline Index rows() const { return m_dec.cols(); }
55  inline Index cols() const { return m_rhs.cols(); }
56  inline const DecompositionType& dec() const { return m_dec; }
57  inline const RhsNestedCleaned& rhs() const { return m_rhs; }
58 
59  template<typename Dest> inline void evalTo(Dest& dst) const
60  {
61  static_cast<const sparse_solve_retval<DecompositionType,Rhs>*>(this)->evalTo(dst);
62  }
63 
64  protected:
65  const DecompositionType& m_dec;
66  typename Rhs::Nested m_rhs;
67 };
68 
69 #define EIGEN_MAKE_SPARSE_SOLVE_HELPERS(DecompositionType,Rhs) \
70  typedef typename DecompositionType::MatrixType MatrixType; \
71  typedef typename MatrixType::Scalar Scalar; \
72  typedef typename MatrixType::RealScalar RealScalar; \
73  typedef typename MatrixType::Index Index; \
74  typedef Eigen::internal::sparse_solve_retval_base<DecompositionType,Rhs> Base; \
75  using Base::dec; \
76  using Base::rhs; \
77  using Base::rows; \
78  using Base::cols; \
79  sparse_solve_retval(const DecompositionType& dec, const Rhs& rhs) \
80  : Base(dec, rhs) {}
81 
82 
83 
84 template<typename DecompositionType, typename Rhs, typename Guess> struct solve_retval_with_guess;
85 
86 template<typename DecompositionType, typename Rhs, typename Guess>
87 struct traits<solve_retval_with_guess<DecompositionType, Rhs, Guess> >
88 {
89  typedef typename DecompositionType::MatrixType MatrixType;
90  typedef Matrix<typename Rhs::Scalar,
91  MatrixType::ColsAtCompileTime,
92  Rhs::ColsAtCompileTime,
93  Rhs::PlainObject::Options,
94  MatrixType::MaxColsAtCompileTime,
95  Rhs::MaxColsAtCompileTime> ReturnType;
96 };
97 
98 template<typename DecompositionType, typename Rhs, typename Guess> struct solve_retval_with_guess
99  : public ReturnByValue<solve_retval_with_guess<DecompositionType, Rhs, Guess> >
100 {
101  typedef typename DecompositionType::Index Index;
102 
103  solve_retval_with_guess(const DecompositionType& dec, const Rhs& rhs, const Guess& guess)
104  : m_dec(dec), m_rhs(rhs), m_guess(guess)
105  {}
106 
107  inline Index rows() const { return m_dec.cols(); }
108  inline Index cols() const { return m_rhs.cols(); }
109 
110  template<typename Dest> inline void evalTo(Dest& dst) const
111  {
112  dst = m_guess;
113  m_dec._solveWithGuess(m_rhs,dst);
114  }
115 
116  protected:
117  const DecompositionType& m_dec;
118  const typename Rhs::Nested m_rhs;
119  const typename Guess::Nested m_guess;
120 };
121 
122 } // namepsace internal
123 
124 } // end namespace Eigen
125 
126 #endif // EIGEN_SPARSE_SOLVE_H