10 #ifndef EIGEN_SPARSE_TRIANGULARVIEW_H
11 #define EIGEN_SPARSE_TRIANGULARVIEW_H
17 template<
typename MatrixType,
int Mode>
18 struct traits<SparseTriangularView<MatrixType,Mode> >
19 :
public traits<MatrixType>
24 template<
typename MatrixType,
int Mode>
class SparseTriangularView
25 :
public SparseMatrixBase<SparseTriangularView<MatrixType,Mode> >
29 SkipLast = !SkipFirst,
30 HasUnitDiag = (Mode&
UnitDiag) ? 1 : 0
35 EIGEN_SPARSE_PUBLIC_INTERFACE(SparseTriangularView)
38 class ReverseInnerIterator;
40 inline Index rows()
const {
return m_matrix.rows(); }
41 inline Index cols()
const {
return m_matrix.cols(); }
43 typedef typename MatrixType::Nested MatrixTypeNested;
44 typedef typename internal::remove_reference<MatrixTypeNested>::type MatrixTypeNestedNonRef;
45 typedef typename internal::remove_all<MatrixTypeNested>::type MatrixTypeNestedCleaned;
47 inline SparseTriangularView(
const MatrixType& matrix) : m_matrix(matrix) {}
50 inline const MatrixTypeNestedCleaned& nestedExpression()
const {
return m_matrix; }
52 template<
typename OtherDerived>
53 typename internal::plain_matrix_type_column_major<OtherDerived>::type
54 solve(
const MatrixBase<OtherDerived>& other)
const;
56 template<
typename OtherDerived>
void solveInPlace(MatrixBase<OtherDerived>& other)
const;
57 template<
typename OtherDerived>
void solveInPlace(SparseMatrixBase<OtherDerived>& other)
const;
60 MatrixTypeNested m_matrix;
63 template<
typename MatrixType,
int Mode>
64 class SparseTriangularView<MatrixType,Mode>::InnerIterator :
public MatrixTypeNestedCleaned::InnerIterator
66 typedef typename MatrixTypeNestedCleaned::InnerIterator Base;
69 EIGEN_STRONG_INLINE InnerIterator(
const SparseTriangularView& view, Index outer)
70 : Base(view.nestedExpression(), outer), m_returnOne(false)
74 while((*
this) && (HasUnitDiag ? this->index()<=outer : this->index()<outer))
79 else if(HasUnitDiag && ((!Base::operator
bool()) || Base::index()>=Base::outer()))
81 if((!SkipFirst) && Base::operator
bool())
87 EIGEN_STRONG_INLINE InnerIterator& operator++()
89 if(HasUnitDiag && m_returnOne)
94 if(HasUnitDiag && (!SkipFirst) && ((!Base::operator
bool()) || Base::index()>=Base::outer()))
96 if((!SkipFirst) && Base::operator
bool())
104 inline Index row()
const {
return Base::row(); }
105 inline Index col()
const {
return Base::col(); }
106 inline Index index()
const
108 if(HasUnitDiag && m_returnOne)
return Base::outer();
109 else return Base::index();
111 inline Scalar value()
const
113 if(HasUnitDiag && m_returnOne)
return Scalar(1);
114 else return Base::value();
117 EIGEN_STRONG_INLINE
operator bool()
const
119 if(HasUnitDiag && m_returnOne)
121 return (SkipFirst ? Base::operator
bool() : (Base::operator
bool() && this->index() <= this->outer()));
127 template<
typename MatrixType,
int Mode>
128 class SparseTriangularView<MatrixType,Mode>::ReverseInnerIterator :
public MatrixTypeNestedCleaned::ReverseInnerIterator
130 typedef typename MatrixTypeNestedCleaned::ReverseInnerIterator Base;
133 EIGEN_STRONG_INLINE ReverseInnerIterator(
const SparseTriangularView& view, Index outer)
134 : Base(view.nestedExpression(), outer)
136 eigen_assert((!HasUnitDiag) &&
"ReverseInnerIterator does not support yet triangular views with a unit diagonal");
138 while((*
this) && this->index()>outer)
142 EIGEN_STRONG_INLINE InnerIterator& operator--()
143 { Base::operator--();
return *
this; }
145 inline Index row()
const {
return Base::row(); }
146 inline Index col()
const {
return Base::col(); }
148 EIGEN_STRONG_INLINE
operator bool()
const
150 return SkipLast ? Base::operator bool() : (Base::operator bool() && this->index() >= this->outer());
154 template<
typename Derived>
156 inline const SparseTriangularView<Derived, Mode>
157 SparseMatrixBase<Derived>::triangularView()
const
164 #endif // EIGEN_SPARSE_TRIANGULARVIEW_H