SparseTriangularView.h
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9 
10 #ifndef EIGEN_SPARSE_TRIANGULARVIEW_H
11 #define EIGEN_SPARSE_TRIANGULARVIEW_H
12 
13 namespace Eigen {
14 
15 namespace internal {
16 
17 template<typename MatrixType, int Mode>
18 struct traits<SparseTriangularView<MatrixType,Mode> >
19 : public traits<MatrixType>
20 {};
21 
22 } // namespace internal
23 
24 template<typename MatrixType, int Mode> class SparseTriangularView
25  : public SparseMatrixBase<SparseTriangularView<MatrixType,Mode> >
26 {
27  enum { SkipFirst = ((Mode&Lower) && !(MatrixType::Flags&RowMajorBit))
28  || ((Mode&Upper) && (MatrixType::Flags&RowMajorBit)),
29  SkipLast = !SkipFirst,
30  HasUnitDiag = (Mode&UnitDiag) ? 1 : 0
31  };
32 
33  public:
34 
35  EIGEN_SPARSE_PUBLIC_INTERFACE(SparseTriangularView)
36 
37  class InnerIterator;
38  class ReverseInnerIterator;
39 
40  inline Index rows() const { return m_matrix.rows(); }
41  inline Index cols() const { return m_matrix.cols(); }
42 
43  typedef typename MatrixType::Nested MatrixTypeNested;
44  typedef typename internal::remove_reference<MatrixTypeNested>::type MatrixTypeNestedNonRef;
45  typedef typename internal::remove_all<MatrixTypeNested>::type MatrixTypeNestedCleaned;
46 
47  inline SparseTriangularView(const MatrixType& matrix) : m_matrix(matrix) {}
48 
50  inline const MatrixTypeNestedCleaned& nestedExpression() const { return m_matrix; }
51 
52  template<typename OtherDerived>
53  typename internal::plain_matrix_type_column_major<OtherDerived>::type
54  solve(const MatrixBase<OtherDerived>& other) const;
55 
56  template<typename OtherDerived> void solveInPlace(MatrixBase<OtherDerived>& other) const;
57  template<typename OtherDerived> void solveInPlace(SparseMatrixBase<OtherDerived>& other) const;
58 
59  protected:
60  MatrixTypeNested m_matrix;
61 };
62 
63 template<typename MatrixType, int Mode>
64 class SparseTriangularView<MatrixType,Mode>::InnerIterator : public MatrixTypeNestedCleaned::InnerIterator
65 {
66  typedef typename MatrixTypeNestedCleaned::InnerIterator Base;
67  public:
68 
69  EIGEN_STRONG_INLINE InnerIterator(const SparseTriangularView& view, Index outer)
70  : Base(view.nestedExpression(), outer), m_returnOne(false)
71  {
72  if(SkipFirst)
73  {
74  while((*this) && (HasUnitDiag ? this->index()<=outer : this->index()<outer))
75  Base::operator++();
76  if(HasUnitDiag)
77  m_returnOne = true;
78  }
79  else if(HasUnitDiag && ((!Base::operator bool()) || Base::index()>=Base::outer()))
80  {
81  if((!SkipFirst) && Base::operator bool())
82  Base::operator++();
83  m_returnOne = true;
84  }
85  }
86 
87  EIGEN_STRONG_INLINE InnerIterator& operator++()
88  {
89  if(HasUnitDiag && m_returnOne)
90  m_returnOne = false;
91  else
92  {
93  Base::operator++();
94  if(HasUnitDiag && (!SkipFirst) && ((!Base::operator bool()) || Base::index()>=Base::outer()))
95  {
96  if((!SkipFirst) && Base::operator bool())
97  Base::operator++();
98  m_returnOne = true;
99  }
100  }
101  return *this;
102  }
103 
104  inline Index row() const { return Base::row(); }
105  inline Index col() const { return Base::col(); }
106  inline Index index() const
107  {
108  if(HasUnitDiag && m_returnOne) return Base::outer();
109  else return Base::index();
110  }
111  inline Scalar value() const
112  {
113  if(HasUnitDiag && m_returnOne) return Scalar(1);
114  else return Base::value();
115  }
116 
117  EIGEN_STRONG_INLINE operator bool() const
118  {
119  if(HasUnitDiag && m_returnOne)
120  return true;
121  return (SkipFirst ? Base::operator bool() : (Base::operator bool() && this->index() <= this->outer()));
122  }
123  protected:
124  bool m_returnOne;
125 };
126 
127 template<typename MatrixType, int Mode>
128 class SparseTriangularView<MatrixType,Mode>::ReverseInnerIterator : public MatrixTypeNestedCleaned::ReverseInnerIterator
129 {
130  typedef typename MatrixTypeNestedCleaned::ReverseInnerIterator Base;
131  public:
132 
133  EIGEN_STRONG_INLINE ReverseInnerIterator(const SparseTriangularView& view, Index outer)
134  : Base(view.nestedExpression(), outer)
135  {
136  eigen_assert((!HasUnitDiag) && "ReverseInnerIterator does not support yet triangular views with a unit diagonal");
137  if(SkipLast)
138  while((*this) && this->index()>outer)
139  --(*this);
140  }
141 
142  EIGEN_STRONG_INLINE InnerIterator& operator--()
143  { Base::operator--(); return *this; }
144 
145  inline Index row() const { return Base::row(); }
146  inline Index col() const { return Base::col(); }
147 
148  EIGEN_STRONG_INLINE operator bool() const
149  {
150  return SkipLast ? Base::operator bool() : (Base::operator bool() && this->index() >= this->outer());
151  }
152 };
153 
154 template<typename Derived>
155 template<int Mode>
156 inline const SparseTriangularView<Derived, Mode>
157 SparseMatrixBase<Derived>::triangularView() const
158 {
159  return derived();
160 }
161 
162 } // end namespace Eigen
163 
164 #endif // EIGEN_SPARSE_TRIANGULARVIEW_H