11 #ifndef EIGEN_TRIANGULARMATRIX_H
12 #define EIGEN_TRIANGULARMATRIX_H
14 #include "./InternalHeaderCheck.h"
20 template<
int S
ide,
typename TriangularType,
typename Rhs>
struct triangular_solve_retval;
34 Mode = internal::traits<Derived>::Mode,
35 RowsAtCompileTime = internal::traits<Derived>::RowsAtCompileTime,
36 ColsAtCompileTime = internal::traits<Derived>::ColsAtCompileTime,
37 MaxRowsAtCompileTime = internal::traits<Derived>::MaxRowsAtCompileTime,
38 MaxColsAtCompileTime = internal::traits<Derived>::MaxColsAtCompileTime,
45 MaxSizeAtCompileTime = internal::size_at_compile_time(internal::traits<Derived>::MaxRowsAtCompileTime,
46 internal::traits<Derived>::MaxColsAtCompileTime)
49 typedef typename internal::traits<Derived>::Scalar Scalar;
50 typedef typename internal::traits<Derived>::StorageKind StorageKind;
51 typedef typename internal::traits<Derived>::StorageIndex StorageIndex;
52 typedef typename internal::traits<Derived>::FullMatrixType DenseMatrixType;
53 typedef DenseMatrixType DenseType;
54 typedef Derived
const& Nested;
59 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR
60 inline Index rows()
const EIGEN_NOEXCEPT {
return derived().rows(); }
61 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR
62 inline Index cols()
const EIGEN_NOEXCEPT {
return derived().cols(); }
63 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR
64 inline Index outerStride()
const EIGEN_NOEXCEPT {
return derived().outerStride(); }
65 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR
66 inline Index innerStride()
const EIGEN_NOEXCEPT {
return derived().innerStride(); }
72 EIGEN_UNUSED_VARIABLE(rows);
73 EIGEN_UNUSED_VARIABLE(cols);
74 eigen_assert(rows==this->rows() && cols==this->cols());
78 inline Scalar coeff(
Index row,
Index col)
const {
return derived().coeff(row,col); }
80 inline Scalar& coeffRef(
Index row,
Index col) {
return derived().coeffRef(row,col); }
84 template<
typename Other>
88 derived().coeffRef(row, col) = other.coeff(row, col);
92 inline Scalar operator()(
Index row,
Index col)
const
94 check_coordinates(row, col);
95 return coeff(row,col);
98 inline Scalar& operator()(
Index row,
Index col)
100 check_coordinates(row, col);
101 return coeffRef(row,col);
104 #ifndef EIGEN_PARSED_BY_DOXYGEN
106 inline const Derived&
derived()
const {
return *
static_cast<const Derived*
>(
this); }
108 inline Derived&
derived() {
return *
static_cast<Derived*
>(
this); }
111 template<
typename DenseDerived>
114 template<
typename DenseDerived>
119 DenseMatrixType toDenseMatrix()
const
121 DenseMatrixType res(rows(), cols());
128 void check_coordinates(
Index row,
Index col)
const
130 EIGEN_ONLY_USED_FOR_DEBUG(row);
131 EIGEN_ONLY_USED_FOR_DEBUG(col);
132 eigen_assert(col>=0 && col<cols() && row>=0 && row<rows());
134 EIGEN_ONLY_USED_FOR_DEBUG(mode);
135 eigen_assert((mode==
Upper && col>=row)
136 || (mode==
Lower && col<=row)
141 #ifdef EIGEN_INTERNAL_DEBUGGING
142 void check_coordinates_internal(
Index row,
Index col)
const
144 check_coordinates(row, col);
147 void check_coordinates_internal(
Index ,
Index )
const {}
170 template<
typename MatrixType,
unsigned int Mode_>
171 struct traits<TriangularView<MatrixType, Mode_> > : traits<MatrixType>
173 typedef typename ref_selector<MatrixType>::non_const_type MatrixTypeNested;
174 typedef std::remove_reference_t<MatrixTypeNested> MatrixTypeNestedNonRef;
175 typedef remove_all_t<MatrixTypeNested> MatrixTypeNestedCleaned;
176 typedef typename MatrixType::PlainObject FullMatrixType;
177 typedef MatrixType ExpressionType;
180 FlagsLvalueBit = is_lvalue<MatrixType>::value ?
LvalueBit : 0,
186 template<
typename MatrixType_,
unsigned int Mode_,
typename StorageKind>
class TriangularViewImpl;
189 :
public TriangularViewImpl<MatrixType_, Mode_, typename internal::traits<MatrixType_>::StorageKind >
193 typedef TriangularViewImpl<MatrixType_, Mode_, typename internal::traits<MatrixType_>::StorageKind > Base;
194 typedef typename internal::traits<TriangularView>::Scalar Scalar;
195 typedef MatrixType_ MatrixType;
198 typedef typename internal::traits<TriangularView>::MatrixTypeNested MatrixTypeNested;
199 typedef typename internal::traits<TriangularView>::MatrixTypeNestedNonRef MatrixTypeNestedNonRef;
201 typedef internal::remove_all_t<typename MatrixType::ConjugateReturnType> MatrixConjugateReturnType;
206 typedef typename internal::traits<TriangularView>::StorageKind StorageKind;
207 typedef typename internal::traits<TriangularView>::MatrixTypeNestedCleaned NestedExpression;
211 Flags = internal::traits<TriangularView>::Flags,
216 IsVectorAtCompileTime =
false
220 explicit inline TriangularView(MatrixType& matrix) : m_matrix(matrix)
226 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR
227 inline Index rows() const EIGEN_NOEXCEPT {
return m_matrix.rows(); }
229 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR
230 inline Index cols() const EIGEN_NOEXCEPT {
return m_matrix.cols(); }
251 inline std::conditional_t<Cond,ConjugateReturnType,ConstTriangularView>
254 typedef std::conditional_t<Cond,ConjugateReturnType,ConstTriangularView> ReturnType;
255 return ReturnType(m_matrix.template conjugateIf<Cond>());
266 template<
class Dummy=
int>
270 typename MatrixType::TransposeReturnType tmp(m_matrix);
282 template<
typename Other>
290 template<
int S
ide,
typename Other>
292 inline const internal::triangular_solve_retval<Side,TriangularView, Other>
293 solve(
const MatrixBase<Other>& other)
const
294 {
return Base::template solve<Side>(other); }
329 return m_matrix.diagonal().prod();
334 MatrixTypeNested m_matrix;
346 template<
typename MatrixType_,
unsigned int Mode_>
class TriangularViewImpl<MatrixType_,Mode_,
Dense>
354 typedef typename internal::traits<TriangularViewType>::Scalar Scalar;
356 typedef MatrixType_ MatrixType;
357 typedef typename MatrixType::PlainObject DenseMatrixType;
358 typedef DenseMatrixType PlainObject;
361 using Base::evalToLazy;
364 typedef typename internal::traits<TriangularViewType>::StorageKind StorageKind;
368 Flags = internal::traits<TriangularViewType>::Flags
381 template<
typename Other>
384 internal::call_assignment_no_alias(derived(), other.
derived(), internal::add_assign_op<Scalar,typename Other::Scalar>());
388 template<
typename Other>
391 internal::call_assignment_no_alias(derived(), other.
derived(), internal::sub_assign_op<Scalar,typename Other::Scalar>());
404 void fill(
const Scalar& value) { setConstant(value); }
408 {
return *
this = MatrixType::Constant(derived().rows(), derived().cols(), value); }
422 Base::check_coordinates_internal(row, col);
423 return derived().nestedExpression().coeff(row, col);
433 Base::check_coordinates_internal(row, col);
434 return derived().nestedExpression().coeffRef(row, col);
438 template<
typename OtherDerived>
443 template<
typename OtherDerived>
447 #ifndef EIGEN_PARSED_BY_DOXYGEN
450 {
return *
this = other.derived().nestedExpression(); }
452 template<
typename OtherDerived>
454 EIGEN_DEPRECATED EIGEN_DEVICE_FUNC
457 template<
typename OtherDerived>
459 EIGEN_DEPRECATED EIGEN_DEVICE_FUNC
464 template<
typename OtherDerived>
473 template<
typename OtherDerived>
friend
504 template<
int S
ide,
typename Other>
505 inline const internal::triangular_solve_retval<Side,TriangularViewType, Other>
517 template<
int S
ide,
typename OtherDerived>
521 template<
typename OtherDerived>
524 {
return solveInPlace<OnTheLeft>(other); }
527 template<
typename OtherDerived>
529 #ifdef EIGEN_PARSED_BY_DOXYGEN
535 EIGEN_STATIC_ASSERT_LVALUE(OtherDerived);
536 call_assignment(derived(), other.const_cast_derived(), internal::swap_assign_op<Scalar>());
540 template<
typename OtherDerived>
542 EIGEN_DEPRECATED EIGEN_DEVICE_FUNC
545 EIGEN_STATIC_ASSERT_LVALUE(OtherDerived);
546 call_assignment(derived(), other.const_cast_derived(), internal::swap_assign_op<Scalar>());
549 template<
typename RhsType,
typename DstType>
551 EIGEN_STRONG_INLINE
void _solve_impl(
const RhsType &rhs, DstType &dst)
const {
552 if(!internal::is_same_dense(dst,rhs))
554 this->solveInPlace(dst);
557 template<
typename ProductType>
559 EIGEN_STRONG_INLINE TriangularViewType& _assignProduct(
const ProductType& prod,
const Scalar& alpha,
bool beta);
561 EIGEN_DEFAULT_COPY_CONSTRUCTOR(TriangularViewImpl)
562 EIGEN_DEFAULT_EMPTY_CONSTRUCTOR_AND_DESTRUCTOR(TriangularViewImpl)
570 #ifndef EIGEN_PARSED_BY_DOXYGEN
572 template<
typename MatrixType,
unsigned int Mode>
573 template<
typename OtherDerived>
574 EIGEN_DEVICE_FUNC
inline TriangularView<MatrixType, Mode>&
575 TriangularViewImpl<MatrixType, Mode, Dense>::operator=(
const MatrixBase<OtherDerived>& other)
577 internal::call_assignment_no_alias(derived(), other.derived(), internal::assign_op<Scalar,typename OtherDerived::Scalar>());
582 template<
typename MatrixType,
unsigned int Mode>
583 template<
typename OtherDerived>
584 EIGEN_DEVICE_FUNC
void TriangularViewImpl<MatrixType, Mode, Dense>::lazyAssign(
const MatrixBase<OtherDerived>& other)
586 internal::call_assignment_no_alias(derived(), other.template triangularView<Mode>());
591 template<
typename MatrixType,
unsigned int Mode>
592 template<
typename OtherDerived>
593 EIGEN_DEVICE_FUNC
inline TriangularView<MatrixType, Mode>&
594 TriangularViewImpl<MatrixType, Mode, Dense>::operator=(
const TriangularBase<OtherDerived>& other)
596 eigen_assert(Mode ==
int(OtherDerived::Mode));
597 internal::call_assignment(derived(), other.derived());
601 template<
typename MatrixType,
unsigned int Mode>
602 template<
typename OtherDerived>
603 EIGEN_DEVICE_FUNC
void TriangularViewImpl<MatrixType, Mode, Dense>::lazyAssign(
const TriangularBase<OtherDerived>& other)
605 eigen_assert(Mode ==
int(OtherDerived::Mode));
606 internal::call_assignment_no_alias(derived(), other.derived());
616 template<
typename Derived>
617 template<
typename DenseDerived>
642 template<
typename Derived>
643 template<
unsigned int Mode>
652 template<
typename Derived>
653 template<
unsigned int Mode>
666 template<
typename Derived>
669 RealScalar maxAbsOnUpperPart =
static_cast<RealScalar
>(-1);
670 for(
Index j = 0; j < cols(); ++j)
672 Index maxi = numext::mini(j, rows()-1);
673 for(
Index i = 0; i <= maxi; ++i)
675 RealScalar absValue = numext::abs(coeff(i,j));
676 if(absValue > maxAbsOnUpperPart) maxAbsOnUpperPart = absValue;
679 RealScalar threshold = maxAbsOnUpperPart * prec;
680 for(
Index j = 0; j < cols(); ++j)
681 for(
Index i = j+1; i < rows(); ++i)
682 if(numext::abs(coeff(i, j)) > threshold)
return false;
691 template<
typename Derived>
694 RealScalar maxAbsOnLowerPart =
static_cast<RealScalar
>(-1);
695 for(
Index j = 0; j < cols(); ++j)
696 for(
Index i = j; i < rows(); ++i)
698 RealScalar absValue = numext::abs(coeff(i,j));
699 if(absValue > maxAbsOnLowerPart) maxAbsOnLowerPart = absValue;
701 RealScalar threshold = maxAbsOnLowerPart * prec;
702 for(
Index j = 1; j < cols(); ++j)
704 Index maxi = numext::mini(j, rows()-1);
705 for(
Index i = 0; i < maxi; ++i)
706 if(numext::abs(coeff(i, j)) > threshold)
return false;
724 template<
typename MatrixType,
unsigned int Mode>
727 typedef typename storage_kind_to_evaluator_kind<typename MatrixType::StorageKind>::Kind Kind;
728 typedef typename glue_shapes<typename evaluator_traits<MatrixType>::Shape, TriangularShape>::type Shape;
731 template<
typename MatrixType,
unsigned int Mode>
732 struct unary_evaluator<TriangularView<MatrixType,Mode>, IndexBased>
733 : evaluator<internal::remove_all_t<MatrixType>>
735 typedef TriangularView<MatrixType,Mode> XprType;
736 typedef evaluator<internal::remove_all_t<MatrixType>> Base;
738 unary_evaluator(
const XprType &xpr) : Base(xpr.nestedExpression()) {}
742 struct Triangular2Triangular {};
743 struct Triangular2Dense {};
744 struct Dense2Triangular {};
747 template<
typename Kernel,
unsigned int Mode,
int UnrollCount,
bool ClearOpposite>
struct triangular_assignment_loop;
755 template<
int UpLo,
int Mode,
int SetOpposite,
typename DstEvaluatorTypeT,
typename SrcEvaluatorTypeT,
typename Functor,
int Version = Specialized>
756 class triangular_dense_assignment_kernel :
public generic_dense_assignment_kernel<DstEvaluatorTypeT, SrcEvaluatorTypeT, Functor, Version>
759 typedef generic_dense_assignment_kernel<DstEvaluatorTypeT, SrcEvaluatorTypeT, Functor, Version> Base;
760 typedef typename Base::DstXprType DstXprType;
761 typedef typename Base::SrcXprType SrcXprType;
764 using Base::m_functor;
767 typedef typename Base::DstEvaluatorType DstEvaluatorType;
768 typedef typename Base::SrcEvaluatorType SrcEvaluatorType;
769 typedef typename Base::Scalar Scalar;
770 typedef typename Base::AssignmentTraits AssignmentTraits;
773 EIGEN_DEVICE_FUNC triangular_dense_assignment_kernel(DstEvaluatorType &dst,
const SrcEvaluatorType &src,
const Functor &func, DstXprType& dstExpr)
774 : Base(dst, src, func, dstExpr)
777 #ifdef EIGEN_INTERNAL_DEBUGGING
778 EIGEN_DEVICE_FUNC
void assignCoeff(
Index row,
Index col)
780 eigen_internal_assert(row!=col);
781 Base::assignCoeff(row,col);
784 using Base::assignCoeff;
787 EIGEN_DEVICE_FUNC
void assignDiagonalCoeff(
Index id)
789 if(Mode==
UnitDiag && SetOpposite) m_functor.assignCoeff(m_dst.coeffRef(
id,
id), Scalar(1));
790 else if(Mode==
ZeroDiag && SetOpposite) m_functor.assignCoeff(m_dst.coeffRef(
id,
id), Scalar(0));
791 else if(Mode==0) Base::assignCoeff(
id,
id);
794 EIGEN_DEVICE_FUNC
void assignOppositeCoeff(
Index row,
Index col)
796 eigen_internal_assert(row!=col);
798 m_functor.assignCoeff(m_dst.coeffRef(row,col), Scalar(0));
802 template<
int Mode,
bool SetOpposite,
typename DstXprType,
typename SrcXprType,
typename Functor>
803 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
804 void call_triangular_assignment_loop(DstXprType& dst,
const SrcXprType& src,
const Functor &func)
806 typedef evaluator<DstXprType> DstEvaluatorType;
807 typedef evaluator<SrcXprType> SrcEvaluatorType;
809 SrcEvaluatorType srcEvaluator(src);
811 Index dstRows = src.rows();
812 Index dstCols = src.cols();
813 if((dst.rows()!=dstRows) || (dst.cols()!=dstCols))
814 dst.resize(dstRows, dstCols);
815 DstEvaluatorType dstEvaluator(dst);
818 DstEvaluatorType,SrcEvaluatorType,Functor> Kernel;
819 Kernel kernel(dstEvaluator, srcEvaluator, func, dst.const_cast_derived());
822 unroll = DstXprType::SizeAtCompileTime !=
Dynamic
823 && SrcEvaluatorType::CoeffReadCost <
HugeCost
824 && DstXprType::SizeAtCompileTime * (int(DstEvaluatorType::CoeffReadCost) + int(SrcEvaluatorType::CoeffReadCost)) / 2 <= EIGEN_UNROLLING_LIMIT
827 triangular_assignment_loop<Kernel, Mode, unroll ? int(DstXprType::SizeAtCompileTime) :
Dynamic, SetOpposite>::run(kernel);
830 template<
int Mode,
bool SetOpposite,
typename DstXprType,
typename SrcXprType>
831 EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE
832 void call_triangular_assignment_loop(DstXprType& dst,
const SrcXprType& src)
834 call_triangular_assignment_loop<Mode,SetOpposite>(dst, src, internal::assign_op<typename DstXprType::Scalar,typename SrcXprType::Scalar>());
837 template<>
struct AssignmentKind<TriangularShape,TriangularShape> {
typedef Triangular2Triangular Kind; };
838 template<>
struct AssignmentKind<DenseShape,TriangularShape> {
typedef Triangular2Dense Kind; };
839 template<>
struct AssignmentKind<TriangularShape,DenseShape> {
typedef Dense2Triangular Kind; };
842 template<
typename DstXprType,
typename SrcXprType,
typename Functor>
843 struct Assignment<DstXprType, SrcXprType, Functor, Triangular2Triangular>
845 EIGEN_DEVICE_FUNC
static void run(DstXprType &dst,
const SrcXprType &src,
const Functor &func)
847 eigen_assert(
int(DstXprType::Mode) ==
int(SrcXprType::Mode));
849 call_triangular_assignment_loop<DstXprType::Mode, false>(dst, src, func);
853 template<
typename DstXprType,
typename SrcXprType,
typename Functor>
854 struct Assignment<DstXprType, SrcXprType, Functor, Triangular2Dense>
856 EIGEN_DEVICE_FUNC
static void run(DstXprType &dst,
const SrcXprType &src,
const Functor &func)
858 call_triangular_assignment_loop<SrcXprType::Mode, (int(SrcXprType::Mode) & int(SelfAdjoint)) == 0>(dst, src, func);
862 template<
typename DstXprType,
typename SrcXprType,
typename Functor>
863 struct Assignment<DstXprType, SrcXprType, Functor, Dense2Triangular>
865 EIGEN_DEVICE_FUNC
static void run(DstXprType &dst,
const SrcXprType &src,
const Functor &func)
867 call_triangular_assignment_loop<DstXprType::Mode, false>(dst, src, func);
872 template<
typename Kernel,
unsigned int Mode,
int UnrollCount,
bool SetOpposite>
873 struct triangular_assignment_loop
876 typedef typename Kernel::DstEvaluatorType DstEvaluatorType;
877 typedef typename DstEvaluatorType::XprType DstXprType;
880 col = (UnrollCount-1) / DstXprType::RowsAtCompileTime,
881 row = (UnrollCount-1) % DstXprType::RowsAtCompileTime
884 typedef typename Kernel::Scalar Scalar;
887 static inline void run(Kernel &kernel)
889 triangular_assignment_loop<Kernel, Mode, UnrollCount-1, SetOpposite>::run(kernel);
892 kernel.assignDiagonalCoeff(row);
893 else if( ((Mode&
Lower) && row>col) || ((Mode&
Upper) && row<col) )
894 kernel.assignCoeff(row,col);
896 kernel.assignOppositeCoeff(row,col);
901 template<
typename Kernel,
unsigned int Mode,
bool SetOpposite>
902 struct triangular_assignment_loop<Kernel, Mode, 0, SetOpposite>
905 static inline void run(Kernel &) {}
914 template<
typename Kernel,
unsigned int Mode,
bool SetOpposite>
915 struct triangular_assignment_loop<Kernel, Mode,
Dynamic, SetOpposite>
917 typedef typename Kernel::Scalar Scalar;
919 static inline void run(Kernel &kernel)
921 for(
Index j = 0; j < kernel.cols(); ++j)
923 Index maxi = numext::mini(j, kernel.rows());
925 if (((Mode&
Lower) && SetOpposite) || (Mode&
Upper))
928 if(Mode&
Upper) kernel.assignCoeff(i, j);
929 else kernel.assignOppositeCoeff(i, j);
935 kernel.assignDiagonalCoeff(i++);
937 if (((Mode&
Upper) && SetOpposite) || (Mode&
Lower))
939 for(; i < kernel.rows(); ++i)
940 if(Mode&
Lower) kernel.assignCoeff(i, j);
941 else kernel.assignOppositeCoeff(i, j);
951 template<
typename Derived>
952 template<
typename DenseDerived>
955 other.
derived().resize(this->rows(), this->cols());
956 internal::call_triangular_assignment_loop<Derived::Mode, (int(Derived::Mode) & int(
SelfAdjoint)) == 0 >(other.
derived(), derived().nestedExpression());
962 template<
typename DstXprType,
typename Lhs,
typename Rhs,
typename Scalar>
963 struct Assignment<DstXprType,
Product<Lhs,Rhs,DefaultProduct>, internal::assign_op<Scalar,typename Product<Lhs,Rhs,DefaultProduct>::Scalar>, Dense2Triangular>
966 static void run(DstXprType &dst,
const SrcXprType &src,
const internal::assign_op<Scalar,typename SrcXprType::Scalar> &)
968 Index dstRows = src.rows();
969 Index dstCols = src.cols();
970 if((dst.rows()!=dstRows) || (dst.cols()!=dstCols))
971 dst.resize(dstRows, dstCols);
973 dst._assignProduct(src, Scalar(1),
false);
978 template<
typename DstXprType,
typename Lhs,
typename Rhs,
typename Scalar>
979 struct Assignment<DstXprType, Product<Lhs,Rhs,DefaultProduct>, internal::add_assign_op<Scalar,typename Product<Lhs,Rhs,DefaultProduct>::Scalar>, Dense2Triangular>
981 typedef Product<Lhs,Rhs,DefaultProduct> SrcXprType;
982 static void run(DstXprType &dst,
const SrcXprType &src,
const internal::add_assign_op<Scalar,typename SrcXprType::Scalar> &)
984 dst._assignProduct(src, Scalar(1),
true);
989 template<
typename DstXprType,
typename Lhs,
typename Rhs,
typename Scalar>
990 struct Assignment<DstXprType, Product<Lhs,Rhs,DefaultProduct>, internal::sub_assign_op<Scalar,typename Product<Lhs,Rhs,DefaultProduct>::Scalar>, Dense2Triangular>
992 typedef Product<Lhs,Rhs,DefaultProduct> SrcXprType;
993 static void run(DstXprType &dst,
const SrcXprType &src,
const internal::sub_assign_op<Scalar,typename SrcXprType::Scalar> &)
995 dst._assignProduct(src, Scalar(-1),
true);
Base class for all dense matrices, vectors, and arrays.
Definition: DenseBase.h:42
Derived & derived()
Definition: EigenBase.h:48
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:52
bool isLowerTriangular(const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: TriangularMatrix.h:692
bool isUpperTriangular(const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: TriangularMatrix.h:667
Expression of the product of two arbitrary matrices or vectors.
Definition: Product.h:77
Expression of a selfadjoint matrix from a triangular part of a dense matrix.
Definition: SelfAdjointView.h:53
Pseudo expression representing a solving operation.
Definition: Solve.h:65
Expression of the transpose of a matrix.
Definition: Transpose.h:56
Base class for triangular part in a matrix.
Definition: TriangularMatrix.h:30
void copyCoeff(Index row, Index col, Other &other)
Definition: TriangularMatrix.h:86
@ SizeAtCompileTime
Definition: TriangularMatrix.h:40
void evalTo(MatrixBase< DenseDerived > &other) const
Definition: TriangularMatrix.h:618
void evalToLazy(MatrixBase< DenseDerived > &other) const
Definition: TriangularMatrix.h:953
TriangularViewType & setZero()
Definition: TriangularMatrix.h:411
Scalar coeff(Index row, Index col) const
Definition: TriangularMatrix.h:420
TriangularViewType & operator/=(const typename internal::traits< MatrixType >::Scalar &other)
Definition: TriangularMatrix.h:400
void fill(const Scalar &value)
Definition: TriangularMatrix.h:404
TriangularViewType & operator-=(const DenseBase< Other > &other)
Definition: TriangularMatrix.h:390
Index outerStride() const
Definition: TriangularMatrix.h:374
const Product< TriangularViewType, OtherDerived > operator*(const MatrixBase< OtherDerived > &rhs) const
Definition: TriangularMatrix.h:467
EIGEN_DEPRECATED void swap(MatrixBase< OtherDerived > const &other)
Definition: TriangularMatrix.h:543
TriangularViewType & operator=(const MatrixBase< OtherDerived > &other)
Index innerStride() const
Definition: TriangularMatrix.h:378
TriangularViewType & setOnes()
Definition: TriangularMatrix.h:414
TriangularViewType & operator*=(const typename internal::traits< MatrixType >::Scalar &other)
Definition: TriangularMatrix.h:397
TriangularViewType & operator=(const TriangularBase< OtherDerived > &other)
TriangularViewType & setConstant(const Scalar &value)
Definition: TriangularMatrix.h:407
const internal::triangular_solve_retval< Side, TriangularViewType, Other > solve(const MatrixBase< Other > &other) const
friend const Product< OtherDerived, TriangularViewType > operator*(const MatrixBase< OtherDerived > &lhs, const TriangularViewImpl &rhs)
Definition: TriangularMatrix.h:476
void solveInPlace(const MatrixBase< OtherDerived > &other) const
TriangularViewType & operator+=(const DenseBase< Other > &other)
Definition: TriangularMatrix.h:383
Scalar & coeffRef(Index row, Index col)
Definition: TriangularMatrix.h:430
void swap(TriangularBase< OtherDerived > &other)
Definition: TriangularMatrix.h:530
Expression of a triangular part in a matrix.
Definition: TriangularMatrix.h:190
const ConjugateReturnType conjugate() const
Definition: TriangularMatrix.h:243
EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
Definition: TriangularMatrix.h:227
Scalar determinant() const
Definition: TriangularMatrix.h:322
TransposeReturnType transpose(std::enable_if_t< Eigen::internal::is_lvalue< MatrixType >::value, Dummy * >=nullptr)
Definition: TriangularMatrix.h:268
NestedExpression & nestedExpression()
Definition: TriangularMatrix.h:238
const SelfAdjointView< MatrixTypeNestedNonRef, Mode > selfadjointView() const
Definition: TriangularMatrix.h:312
const ConstTransposeReturnType transpose() const
Definition: TriangularMatrix.h:277
std::conditional_t< Cond, ConjugateReturnType, ConstTriangularView > conjugateIf() const
Definition: TriangularMatrix.h:252
const AdjointReturnType adjoint() const
Definition: TriangularMatrix.h:261
const NestedExpression & nestedExpression() const
Definition: TriangularMatrix.h:234
SelfAdjointView< MatrixTypeNestedNonRef, Mode > selfadjointView()
Definition: TriangularMatrix.h:304
EIGEN_CONSTEXPR Index cols() const EIGEN_NOEXCEPT
Definition: TriangularMatrix.h:230
@ StrictlyLower
Definition: Constants.h:223
@ UnitDiag
Definition: Constants.h:215
@ StrictlyUpper
Definition: Constants.h:225
@ UnitLower
Definition: Constants.h:219
@ ZeroDiag
Definition: Constants.h:217
@ SelfAdjoint
Definition: Constants.h:227
@ UnitUpper
Definition: Constants.h:221
@ Lower
Definition: Constants.h:211
@ Upper
Definition: Constants.h:213
const unsigned int PacketAccessBit
Definition: Constants.h:96
const unsigned int LinearAccessBit
Definition: Constants.h:132
const unsigned int DirectAccessBit
Definition: Constants.h:157
const unsigned int LvalueBit
Definition: Constants.h:146
Namespace containing all symbols from the Eigen library.
Definition: Core:139
const int HugeCost
Definition: Constants.h:46
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:59
const int Dynamic
Definition: Constants.h:24
Definition: Constants.h:509
Definition: EigenBase.h:32
Derived & derived()
Definition: EigenBase.h:48
Eigen::Index Index
The interface type of indices.
Definition: EigenBase.h:41