11 #ifndef EIGEN_SPARSE_TRIANGULARVIEW_H 12 #define EIGEN_SPARSE_TRIANGULARVIEW_H 25 template<
typename MatrixType,
unsigned int Mode>
class TriangularViewImpl<MatrixType,Mode,
Sparse>
30 SkipLast = !SkipFirst,
32 HasUnitDiag = (Mode&
UnitDiag) ? 1 : 0
44 EIGEN_SPARSE_PUBLIC_INTERFACE(TriangularViewType)
47 class ReverseInnerIterator;
49 typedef typename MatrixType::Nested MatrixTypeNested;
50 typedef typename internal::remove_reference<MatrixTypeNested>::type MatrixTypeNestedNonRef;
51 typedef typename internal::remove_all<MatrixTypeNested>::type MatrixTypeNestedCleaned;
53 template<
typename RhsType,
typename DstType>
55 EIGEN_STRONG_INLINE
void _solve_impl(
const RhsType &rhs, DstType &dst)
const {
56 if(!(internal::is_same<RhsType,DstType>::value && internal::extract_data(dst) == internal::extract_data(rhs)))
58 this->solveInPlace(dst);
66 template<
typename MatrixType,
unsigned int Mode>
67 class TriangularViewImpl<MatrixType,Mode,Sparse>::InnerIterator :
public MatrixTypeNestedCleaned::InnerIterator
69 typedef typename MatrixTypeNestedCleaned::InnerIterator
Base;
72 EIGEN_STRONG_INLINE InnerIterator(
const TriangularViewImpl& view,
Index outer)
73 : Base(view.derived().nestedExpression(), outer), m_returnOne(
false)
77 while((*
this) && ((HasUnitDiag||SkipDiag) ? this->index()<=outer : this->index()<outer))
82 else if(HasUnitDiag && ((!Base::operator
bool()) || Base::index()>=Base::outer()))
84 if((!SkipFirst) && Base::operator
bool())
90 EIGEN_STRONG_INLINE InnerIterator& operator++()
92 if(HasUnitDiag && m_returnOne)
97 if(HasUnitDiag && (!SkipFirst) && ((!Base::operator
bool()) || Base::index()>=Base::outer()))
99 if((!SkipFirst) && Base::operator
bool())
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 StorageIndex index()
const 111 if(HasUnitDiag && m_returnOne)
return Base::outer();
112 else return Base::index();
114 inline Scalar value()
const 116 if(HasUnitDiag && m_returnOne)
return Scalar(1);
117 else return Base::value();
120 EIGEN_STRONG_INLINE
operator bool()
const 122 if(HasUnitDiag && m_returnOne)
124 if(SkipFirst)
return Base::operator bool();
127 if (SkipDiag)
return (Base::operator
bool() && this->index() < this->outer());
128 else return (Base::operator
bool() && this->index() <= this->outer());
135 template<
typename MatrixType,
unsigned int Mode>
136 class TriangularViewImpl<MatrixType,Mode,Sparse>::ReverseInnerIterator :
public MatrixTypeNestedCleaned::ReverseInnerIterator
138 typedef typename MatrixTypeNestedCleaned::ReverseInnerIterator
Base;
144 eigen_assert((!HasUnitDiag) &&
"ReverseInnerIterator does not support yet triangular views with a unit diagonal");
146 while((*
this) && (SkipDiag ? this->index()>=outer : this->index()>outer))
151 EIGEN_STRONG_INLINE ReverseInnerIterator& operator--()
152 { Base::operator--();
return *
this; }
154 inline Index row()
const {
return Base::row(); }
155 inline Index col()
const {
return Base::col(); }
157 EIGEN_STRONG_INLINE
operator bool()
const 159 if (SkipLast)
return Base::operator bool() ;
162 if(SkipDiag)
return (Base::operator
bool() && this->index() > this->outer());
163 else return (Base::operator
bool() && this->index() >= this->outer());
170 template<
typename ArgType,
unsigned int Mode>
171 struct unary_evaluator<TriangularView<ArgType,Mode>, IteratorBased>
172 : evaluator_base<TriangularView<ArgType,Mode> >
178 typedef typename XprType::Scalar Scalar;
179 typedef typename XprType::StorageIndex StorageIndex;
180 typedef typename evaluator<ArgType>::InnerIterator EvalIterator;
184 SkipLast = !SkipFirst,
186 HasUnitDiag = (Mode&
UnitDiag) ? 1 : 0
192 CoeffReadCost = evaluator<ArgType>::CoeffReadCost,
193 Flags = XprType::Flags
196 explicit unary_evaluator(
const XprType &xpr) : m_argImpl(xpr.nestedExpression()) {}
198 inline Index nonZerosEstimate()
const {
199 return m_argImpl.nonZerosEstimate();
202 class InnerIterator :
public EvalIterator
204 typedef EvalIterator
Base;
207 EIGEN_STRONG_INLINE InnerIterator(
const unary_evaluator& xprEval,
Index outer)
208 : Base(xprEval.m_argImpl,outer), m_returnOne(
false)
212 while((*
this) && ((HasUnitDiag||SkipDiag) ? this->index()<=outer : this->index()<outer))
217 else if(HasUnitDiag && ((!Base::operator
bool()) || Base::index()>=Base::outer()))
219 if((!SkipFirst) && Base::operator
bool())
225 EIGEN_STRONG_INLINE InnerIterator& operator++()
227 if(HasUnitDiag && m_returnOne)
232 if(HasUnitDiag && (!SkipFirst) && ((!Base::operator
bool()) || Base::index()>=Base::outer()))
234 if((!SkipFirst) && Base::operator
bool())
242 EIGEN_STRONG_INLINE
operator bool()
const 244 if(HasUnitDiag && m_returnOne)
246 if(SkipFirst)
return Base::operator bool();
249 if (SkipDiag)
return (Base::operator
bool() && this->index() < this->outer());
250 else return (Base::operator
bool() && this->index() <= this->outer());
256 inline StorageIndex index()
const 258 if(HasUnitDiag && m_returnOne)
return internal::convert_index<StorageIndex>(Base::outer());
259 else return Base::index();
261 inline Scalar value()
const 263 if(HasUnitDiag && m_returnOne)
return Scalar(1);
264 else return Base::value();
274 evaluator<ArgType> m_argImpl;
279 template<
typename Derived>
289 #endif // EIGEN_SPARSE_TRIANGULARVIEW_H const NestedExpression & nestedExpression() const
Definition: TriangularMatrix.h:234
Definition: Constants.h:210
Eigen::Index Index
The interface type of indices.
Definition: EigenBase.h:37
const unsigned int RowMajorBit
Definition: Constants.h:61
Definition: Constants.h:204
Base class of any sparse matrices or sparse expressions.
Definition: ForwardDeclarations.h:278
Definition: Constants.h:493
Definition: Eigen_Colamd.h:54
Expression of a triangular part in a matrix.
Definition: TriangularMatrix.h:186
Definition: Constants.h:208
Definition: Constants.h:206
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48