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