11 #ifndef EIGEN_SPARSE_TRIANGULARVIEW_H
12 #define EIGEN_SPARSE_TRIANGULARVIEW_H
14 #include "./InternalHeaderCheck.h"
27 template<
typename MatrixType,
unsigned int Mode>
class TriangularViewImpl<MatrixType,Mode,
Sparse>
32 SkipLast = !SkipFirst,
34 HasUnitDiag = (Mode&
UnitDiag) ? 1 : 0
48 typedef typename MatrixType::Nested MatrixTypeNested;
49 typedef std::remove_reference_t<MatrixTypeNested> MatrixTypeNestedNonRef;
50 typedef internal::remove_all_t<MatrixTypeNested> MatrixTypeNestedCleaned;
52 template<
typename RhsType,
typename DstType>
54 EIGEN_STRONG_INLINE
void _solve_impl(
const RhsType &rhs, DstType &dst)
const {
55 if(!(internal::is_same<RhsType,DstType>::value && internal::extract_data(dst) == internal::extract_data(rhs)))
57 this->solveInPlace(dst);
70 template<
typename ArgType,
unsigned int Mode>
71 struct unary_evaluator<
TriangularView<ArgType,Mode>, IteratorBased>
72 : evaluator_base<TriangularView<ArgType,Mode> >
78 typedef typename XprType::Scalar Scalar;
79 typedef typename XprType::StorageIndex StorageIndex;
80 typedef typename evaluator<ArgType>::InnerIterator EvalIterator;
84 SkipLast = !SkipFirst,
86 HasUnitDiag = (Mode&
UnitDiag) ? 1 : 0
92 CoeffReadCost = evaluator<ArgType>::CoeffReadCost,
93 Flags = XprType::Flags
96 explicit unary_evaluator(
const XprType &xpr) : m_argImpl(xpr.nestedExpression()), m_arg(xpr.nestedExpression()) {}
98 inline Index nonZerosEstimate()
const {
99 return m_argImpl.nonZerosEstimate();
102 class InnerIterator :
public EvalIterator
104 typedef EvalIterator Base;
107 EIGEN_STRONG_INLINE InnerIterator(
const unary_evaluator& xprEval,
Index outer)
108 : Base(xprEval.m_argImpl,outer), m_returnOne(false), m_containsDiag(Base::outer()<xprEval.m_arg.innerSize())
112 while((*
this) && ((HasUnitDiag||SkipDiag) ? this->index()<=outer : this->index()<outer))
115 m_returnOne = m_containsDiag;
117 else if(HasUnitDiag && ((!Base::operator
bool()) || Base::index()>=Base::outer()))
119 if((!SkipFirst) && Base::operator
bool())
121 m_returnOne = m_containsDiag;
125 EIGEN_STRONG_INLINE InnerIterator& operator++()
127 if(HasUnitDiag && m_returnOne)
132 if(HasUnitDiag && (!SkipFirst) && ((!Base::operator
bool()) || Base::index()>=Base::outer()))
134 if((!SkipFirst) && Base::operator
bool())
136 m_returnOne = m_containsDiag;
142 EIGEN_STRONG_INLINE
operator bool()
const
144 if(HasUnitDiag && m_returnOne)
146 if(SkipFirst)
return Base::operator bool();
149 if (SkipDiag)
return (Base::operator
bool() && this->index() < this->outer());
150 else return (Base::operator
bool() && this->index() <= this->outer());
156 inline StorageIndex index()
const
158 if(HasUnitDiag && m_returnOne)
return internal::convert_index<StorageIndex>(Base::outer());
159 else return Base::index();
161 inline Scalar value()
const
163 if(HasUnitDiag && m_returnOne)
return Scalar(1);
164 else return Base::value();
175 evaluator<ArgType> m_argImpl;
176 const ArgType& m_arg;
181 template<
typename Derived>
183 inline const TriangularView<const Derived, Mode>
184 SparseMatrixBase<Derived>::triangularView()
const
186 return TriangularView<const Derived, Mode>(derived());
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:52
Base class of any sparse matrices or sparse expressions.
Definition: SparseMatrixBase.h:30
void solveInPlace(MatrixBase< OtherDerived > &other) const
void solveInPlace(SparseMatrixBase< OtherDerived > &other) const
Expression of a triangular part in a matrix.
Definition: TriangularMatrix.h:190
@ UnitDiag
Definition: Constants.h:215
@ ZeroDiag
Definition: Constants.h:217
@ Lower
Definition: Constants.h:211
@ Upper
Definition: Constants.h:213
const unsigned int RowMajorBit
Definition: Constants.h:68
Namespace containing all symbols from the Eigen library.
Definition: Core:139
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:59
Definition: Constants.h:512