00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026 #ifndef EIGEN_TRIANGULARMATRIX_H
00027 #define EIGEN_TRIANGULARMATRIX_H
00028
00029 namespace internal {
00030
00031 template<int Side, typename TriangularType, typename Rhs> struct triangular_solve_retval;
00032
00033 }
00034
00042 template<typename Derived> class TriangularBase : public EigenBase<Derived>
00043 {
00044 public:
00045
00046 enum {
00047 Mode = internal::traits<Derived>::Mode,
00048 CoeffReadCost = internal::traits<Derived>::CoeffReadCost,
00049 RowsAtCompileTime = internal::traits<Derived>::RowsAtCompileTime,
00050 ColsAtCompileTime = internal::traits<Derived>::ColsAtCompileTime,
00051 MaxRowsAtCompileTime = internal::traits<Derived>::MaxRowsAtCompileTime,
00052 MaxColsAtCompileTime = internal::traits<Derived>::MaxColsAtCompileTime
00053 };
00054 typedef typename internal::traits<Derived>::Scalar Scalar;
00055 typedef typename internal::traits<Derived>::StorageKind StorageKind;
00056 typedef typename internal::traits<Derived>::Index Index;
00057 typedef typename internal::traits<Derived>::DenseMatrixType DenseMatrixType;
00058 typedef DenseMatrixType DenseType;
00059
00060 inline TriangularBase() { eigen_assert(!((Mode&UnitDiag) && (Mode&ZeroDiag))); }
00061
00062 inline Index rows() const { return derived().rows(); }
00063 inline Index cols() const { return derived().cols(); }
00064 inline Index outerStride() const { return derived().outerStride(); }
00065 inline Index innerStride() const { return derived().innerStride(); }
00066
00067 inline Scalar coeff(Index row, Index col) const { return derived().coeff(row,col); }
00068 inline Scalar& coeffRef(Index row, Index col) { return derived().coeffRef(row,col); }
00069
00072 template<typename Other>
00073 EIGEN_STRONG_INLINE void copyCoeff(Index row, Index col, Other& other)
00074 {
00075 derived().coeffRef(row, col) = other.coeff(row, col);
00076 }
00077
00078 inline Scalar operator()(Index row, Index col) const
00079 {
00080 check_coordinates(row, col);
00081 return coeff(row,col);
00082 }
00083 inline Scalar& operator()(Index row, Index col)
00084 {
00085 check_coordinates(row, col);
00086 return coeffRef(row,col);
00087 }
00088
00089 #ifndef EIGEN_PARSED_BY_DOXYGEN
00090 inline const Derived& derived() const { return *static_cast<const Derived*>(this); }
00091 inline Derived& derived() { return *static_cast<Derived*>(this); }
00092 #endif // not EIGEN_PARSED_BY_DOXYGEN
00093
00094 template<typename DenseDerived>
00095 void evalTo(MatrixBase<DenseDerived> &other) const;
00096 template<typename DenseDerived>
00097 void evalToLazy(MatrixBase<DenseDerived> &other) const;
00098
00099 DenseMatrixType toDenseMatrix() const
00100 {
00101 DenseMatrixType res(rows(), cols());
00102 evalToLazy(res);
00103 return res;
00104 }
00105
00106 protected:
00107
00108 void check_coordinates(Index row, Index col) const
00109 {
00110 EIGEN_ONLY_USED_FOR_DEBUG(row);
00111 EIGEN_ONLY_USED_FOR_DEBUG(col);
00112 eigen_assert(col>=0 && col<cols() && row>=0 && row<rows());
00113 const int mode = int(Mode) & ~SelfAdjoint;
00114 eigen_assert((mode==Upper && col>=row)
00115 || (mode==Lower && col<=row)
00116 || ((mode==StrictlyUpper || mode==UnitUpper) && col>row)
00117 || ((mode==StrictlyLower || mode==UnitLower) && col<row));
00118 }
00119
00120 #ifdef EIGEN_INTERNAL_DEBUGGING
00121 void check_coordinates_internal(Index row, Index col) const
00122 {
00123 check_coordinates(row, col);
00124 }
00125 #else
00126 void check_coordinates_internal(Index , Index ) const {}
00127 #endif
00128
00129 };
00130
00148 namespace internal {
00149 template<typename MatrixType, unsigned int _Mode>
00150 struct traits<TriangularView<MatrixType, _Mode> > : traits<MatrixType>
00151 {
00152 typedef typename nested<MatrixType>::type MatrixTypeNested;
00153 typedef typename remove_reference<MatrixTypeNested>::type MatrixTypeNestedNonRef;
00154 typedef typename remove_all<MatrixTypeNested>::type MatrixTypeNestedCleaned;
00155 typedef MatrixType ExpressionType;
00156 typedef typename MatrixType::PlainObject DenseMatrixType;
00157 enum {
00158 Mode = _Mode,
00159 Flags = (MatrixTypeNestedCleaned::Flags & (HereditaryBits) & (~(PacketAccessBit | DirectAccessBit | LinearAccessBit))) | Mode,
00160 CoeffReadCost = MatrixTypeNestedCleaned::CoeffReadCost
00161 };
00162 };
00163 }
00164
00165 template<int Mode, bool LhsIsTriangular,
00166 typename Lhs, bool LhsIsVector,
00167 typename Rhs, bool RhsIsVector>
00168 struct TriangularProduct;
00169
00170 template<typename _MatrixType, unsigned int _Mode> class TriangularView
00171 : public TriangularBase<TriangularView<_MatrixType, _Mode> >
00172 {
00173 public:
00174
00175 typedef TriangularBase<TriangularView> Base;
00176 typedef typename internal::traits<TriangularView>::Scalar Scalar;
00177
00178 typedef _MatrixType MatrixType;
00179 typedef typename internal::traits<TriangularView>::DenseMatrixType DenseMatrixType;
00180 typedef DenseMatrixType PlainObject;
00181
00182 protected:
00183 typedef typename internal::traits<TriangularView>::MatrixTypeNested MatrixTypeNested;
00184 typedef typename internal::traits<TriangularView>::MatrixTypeNestedNonRef MatrixTypeNestedNonRef;
00185 typedef typename internal::traits<TriangularView>::MatrixTypeNestedCleaned MatrixTypeNestedCleaned;
00186
00187 typedef typename internal::remove_all<typename MatrixType::ConjugateReturnType>::type MatrixConjugateReturnType;
00188
00189 public:
00190 using Base::evalToLazy;
00191
00192
00193 typedef typename internal::traits<TriangularView>::StorageKind StorageKind;
00194 typedef typename internal::traits<TriangularView>::Index Index;
00195
00196 enum {
00197 Mode = _Mode,
00198 TransposeMode = (Mode & Upper ? Lower : 0)
00199 | (Mode & Lower ? Upper : 0)
00200 | (Mode & (UnitDiag))
00201 | (Mode & (ZeroDiag))
00202 };
00203
00204 inline TriangularView(const MatrixType& matrix) : m_matrix(matrix)
00205 {}
00206
00207 inline Index rows() const { return m_matrix.rows(); }
00208 inline Index cols() const { return m_matrix.cols(); }
00209 inline Index outerStride() const { return m_matrix.outerStride(); }
00210 inline Index innerStride() const { return m_matrix.innerStride(); }
00211
00213 template<typename Other> TriangularView& operator+=(const DenseBase<Other>& other) { return *this = m_matrix + other.derived(); }
00215 template<typename Other> TriangularView& operator-=(const DenseBase<Other>& other) { return *this = m_matrix - other.derived(); }
00217 TriangularView& operator*=(const typename internal::traits<MatrixType>::Scalar& other) { return *this = m_matrix * other; }
00219 TriangularView& operator/=(const typename internal::traits<MatrixType>::Scalar& other) { return *this = m_matrix / other; }
00220
00222 void fill(const Scalar& value) { setConstant(value); }
00224 TriangularView& setConstant(const Scalar& value)
00225 { return *this = MatrixType::Constant(rows(), cols(), value); }
00227 TriangularView& setZero() { return setConstant(Scalar(0)); }
00229 TriangularView& setOnes() { return setConstant(Scalar(1)); }
00230
00234 inline Scalar coeff(Index row, Index col) const
00235 {
00236 Base::check_coordinates_internal(row, col);
00237 return m_matrix.coeff(row, col);
00238 }
00239
00243 inline Scalar& coeffRef(Index row, Index col)
00244 {
00245 Base::check_coordinates_internal(row, col);
00246 return m_matrix.const_cast_derived().coeffRef(row, col);
00247 }
00248
00249 const MatrixTypeNestedCleaned& nestedExpression() const { return m_matrix; }
00250 MatrixTypeNestedCleaned& nestedExpression() { return *const_cast<MatrixTypeNestedCleaned*>(&m_matrix); }
00251
00253 template<typename OtherDerived>
00254 TriangularView& operator=(const TriangularBase<OtherDerived>& other);
00255
00256 template<typename OtherDerived>
00257 TriangularView& operator=(const MatrixBase<OtherDerived>& other);
00258
00259 TriangularView& operator=(const TriangularView& other)
00260 { return *this = other.nestedExpression(); }
00261
00262 template<typename OtherDerived>
00263 void lazyAssign(const TriangularBase<OtherDerived>& other);
00264
00265 template<typename OtherDerived>
00266 void lazyAssign(const MatrixBase<OtherDerived>& other);
00267
00269 inline TriangularView<MatrixConjugateReturnType,Mode> conjugate()
00270 { return m_matrix.conjugate(); }
00272 inline const TriangularView<MatrixConjugateReturnType,Mode> conjugate() const
00273 { return m_matrix.conjugate(); }
00274
00276 inline TriangularView<typename MatrixType::AdjointReturnType,TransposeMode> adjoint()
00277 { return m_matrix.adjoint(); }
00279 inline const TriangularView<typename MatrixType::AdjointReturnType,TransposeMode> adjoint() const
00280 { return m_matrix.adjoint(); }
00281
00283 inline TriangularView<Transpose<MatrixType>,TransposeMode> transpose()
00284 {
00285 EIGEN_STATIC_ASSERT_LVALUE(MatrixType)
00286 return m_matrix.const_cast_derived().transpose();
00287 }
00289 inline const TriangularView<Transpose<MatrixType>,TransposeMode> transpose() const
00290 { return m_matrix.transpose(); }
00291
00293 template<typename OtherDerived>
00294 TriangularProduct<Mode,true,MatrixType,false,OtherDerived,OtherDerived::IsVectorAtCompileTime>
00295 operator*(const MatrixBase<OtherDerived>& rhs) const
00296 {
00297 return TriangularProduct
00298 <Mode,true,MatrixType,false,OtherDerived,OtherDerived::IsVectorAtCompileTime>
00299 (m_matrix, rhs.derived());
00300 }
00301
00303 template<typename OtherDerived> friend
00304 TriangularProduct<Mode,false,OtherDerived,OtherDerived::IsVectorAtCompileTime,MatrixType,false>
00305 operator*(const MatrixBase<OtherDerived>& lhs, const TriangularView& rhs)
00306 {
00307 return TriangularProduct
00308 <Mode,false,OtherDerived,OtherDerived::IsVectorAtCompileTime,MatrixType,false>
00309 (lhs.derived(),rhs.m_matrix);
00310 }
00311
00312 #ifdef EIGEN2_SUPPORT
00313 template<typename OtherDerived>
00314 struct eigen2_product_return_type
00315 {
00316 typedef typename TriangularView<MatrixType,Mode>::DenseMatrixType DenseMatrixType;
00317 typedef typename OtherDerived::PlainObject::DenseType OtherPlainObject;
00318 typedef typename ProductReturnType<DenseMatrixType, OtherPlainObject>::Type ProdRetType;
00319 typedef typename ProdRetType::PlainObject type;
00320 };
00321 template<typename OtherDerived>
00322 const typename eigen2_product_return_type<OtherDerived>::type
00323 operator*(const EigenBase<OtherDerived>& rhs) const
00324 {
00325 typename OtherDerived::PlainObject::DenseType rhsPlainObject;
00326 rhs.evalTo(rhsPlainObject);
00327 return this->toDenseMatrix() * rhsPlainObject;
00328 }
00329 template<typename OtherMatrixType>
00330 bool isApprox(const TriangularView<OtherMatrixType, Mode>& other, typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision()) const
00331 {
00332 return this->toDenseMatrix().isApprox(other.toDenseMatrix(), precision);
00333 }
00334 template<typename OtherDerived>
00335 bool isApprox(const MatrixBase<OtherDerived>& other, typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision()) const
00336 {
00337 return this->toDenseMatrix().isApprox(other, precision);
00338 }
00339 #endif // EIGEN2_SUPPORT
00340
00341 template<int Side, typename Other>
00342 inline const internal::triangular_solve_retval<Side,TriangularView, Other>
00343 solve(const MatrixBase<Other>& other) const;
00344
00345 template<int Side, typename OtherDerived>
00346 void solveInPlace(const MatrixBase<OtherDerived>& other) const;
00347
00348 template<typename Other>
00349 inline const internal::triangular_solve_retval<OnTheLeft,TriangularView, Other>
00350 solve(const MatrixBase<Other>& other) const
00351 { return solve<OnTheLeft>(other); }
00352
00353 template<typename OtherDerived>
00354 void solveInPlace(const MatrixBase<OtherDerived>& other) const
00355 { return solveInPlace<OnTheLeft>(other); }
00356
00357 const SelfAdjointView<MatrixTypeNestedNonRef,Mode> selfadjointView() const
00358 {
00359 EIGEN_STATIC_ASSERT((Mode&UnitDiag)==0,PROGRAMMING_ERROR);
00360 return SelfAdjointView<MatrixTypeNestedNonRef,Mode>(m_matrix);
00361 }
00362 SelfAdjointView<MatrixTypeNestedNonRef,Mode> selfadjointView()
00363 {
00364 EIGEN_STATIC_ASSERT((Mode&UnitDiag)==0,PROGRAMMING_ERROR);
00365 return SelfAdjointView<MatrixTypeNestedNonRef,Mode>(m_matrix);
00366 }
00367
00368 template<typename OtherDerived>
00369 void swap(TriangularBase<OtherDerived> const & other)
00370 {
00371 TriangularView<SwapWrapper<MatrixType>,Mode>(const_cast<MatrixType&>(m_matrix)).lazyAssign(other.derived());
00372 }
00373
00374 template<typename OtherDerived>
00375 void swap(MatrixBase<OtherDerived> const & other)
00376 {
00377 TriangularView<SwapWrapper<MatrixType>,Mode>(const_cast<MatrixType&>(m_matrix)).lazyAssign(other.derived());
00378 }
00379
00380 Scalar determinant() const
00381 {
00382 if (Mode & UnitDiag)
00383 return 1;
00384 else if (Mode & ZeroDiag)
00385 return 0;
00386 else
00387 return m_matrix.diagonal().prod();
00388 }
00389
00390
00391 template<typename ProductDerived, typename Lhs, typename Rhs>
00392 EIGEN_STRONG_INLINE TriangularView& operator=(const ProductBase<ProductDerived, Lhs,Rhs>& other)
00393 {
00394 setZero();
00395 return assignProduct(other,1);
00396 }
00397
00398 template<typename ProductDerived, typename Lhs, typename Rhs>
00399 EIGEN_STRONG_INLINE TriangularView& operator+=(const ProductBase<ProductDerived, Lhs,Rhs>& other)
00400 {
00401 return assignProduct(other,1);
00402 }
00403
00404 template<typename ProductDerived, typename Lhs, typename Rhs>
00405 EIGEN_STRONG_INLINE TriangularView& operator-=(const ProductBase<ProductDerived, Lhs,Rhs>& other)
00406 {
00407 return assignProduct(other,-1);
00408 }
00409
00410
00411 template<typename ProductDerived>
00412 EIGEN_STRONG_INLINE TriangularView& operator=(const ScaledProduct<ProductDerived>& other)
00413 {
00414 setZero();
00415 return assignProduct(other,other.alpha());
00416 }
00417
00418 template<typename ProductDerived>
00419 EIGEN_STRONG_INLINE TriangularView& operator+=(const ScaledProduct<ProductDerived>& other)
00420 {
00421 return assignProduct(other,other.alpha());
00422 }
00423
00424 template<typename ProductDerived>
00425 EIGEN_STRONG_INLINE TriangularView& operator-=(const ScaledProduct<ProductDerived>& other)
00426 {
00427 return assignProduct(other,-other.alpha());
00428 }
00429
00430 protected:
00431
00432 template<typename ProductDerived, typename Lhs, typename Rhs>
00433 EIGEN_STRONG_INLINE TriangularView& assignProduct(const ProductBase<ProductDerived, Lhs,Rhs>& prod, const Scalar& alpha);
00434
00435 const MatrixTypeNested m_matrix;
00436 };
00437
00438
00439
00440
00441
00442 namespace internal {
00443
00444 template<typename Derived1, typename Derived2, unsigned int Mode, int UnrollCount, bool ClearOpposite>
00445 struct triangular_assignment_selector
00446 {
00447 enum {
00448 col = (UnrollCount-1) / Derived1::RowsAtCompileTime,
00449 row = (UnrollCount-1) % Derived1::RowsAtCompileTime
00450 };
00451
00452 typedef typename Derived1::Scalar Scalar;
00453
00454 inline static void run(Derived1 &dst, const Derived2 &src)
00455 {
00456 triangular_assignment_selector<Derived1, Derived2, Mode, UnrollCount-1, ClearOpposite>::run(dst, src);
00457
00458 eigen_assert( Mode == Upper || Mode == Lower
00459 || Mode == StrictlyUpper || Mode == StrictlyLower
00460 || Mode == UnitUpper || Mode == UnitLower);
00461 if((Mode == Upper && row <= col)
00462 || (Mode == Lower && row >= col)
00463 || (Mode == StrictlyUpper && row < col)
00464 || (Mode == StrictlyLower && row > col)
00465 || (Mode == UnitUpper && row < col)
00466 || (Mode == UnitLower && row > col))
00467 dst.copyCoeff(row, col, src);
00468 else if(ClearOpposite)
00469 {
00470 if (Mode&UnitDiag && row==col)
00471 dst.coeffRef(row, col) = Scalar(1);
00472 else
00473 dst.coeffRef(row, col) = Scalar(0);
00474 }
00475 }
00476 };
00477
00478
00479 template<typename Derived1, typename Derived2, unsigned int Mode, bool ClearOpposite>
00480 struct triangular_assignment_selector<Derived1, Derived2, Mode, 0, ClearOpposite>
00481 {
00482 inline static void run(Derived1 &, const Derived2 &) {}
00483 };
00484
00485 template<typename Derived1, typename Derived2, bool ClearOpposite>
00486 struct triangular_assignment_selector<Derived1, Derived2, Upper, Dynamic, ClearOpposite>
00487 {
00488 typedef typename Derived1::Index Index;
00489 typedef typename Derived1::Scalar Scalar;
00490 inline static void run(Derived1 &dst, const Derived2 &src)
00491 {
00492 for(Index j = 0; j < dst.cols(); ++j)
00493 {
00494 Index maxi = std::min(j, dst.rows()-1);
00495 for(Index i = 0; i <= maxi; ++i)
00496 dst.copyCoeff(i, j, src);
00497 if (ClearOpposite)
00498 for(Index i = maxi+1; i < dst.rows(); ++i)
00499 dst.coeffRef(i, j) = Scalar(0);
00500 }
00501 }
00502 };
00503
00504 template<typename Derived1, typename Derived2, bool ClearOpposite>
00505 struct triangular_assignment_selector<Derived1, Derived2, Lower, Dynamic, ClearOpposite>
00506 {
00507 typedef typename Derived1::Index Index;
00508 inline static void run(Derived1 &dst, const Derived2 &src)
00509 {
00510 for(Index j = 0; j < dst.cols(); ++j)
00511 {
00512 for(Index i = j; i < dst.rows(); ++i)
00513 dst.copyCoeff(i, j, src);
00514 Index maxi = std::min(j, dst.rows());
00515 if (ClearOpposite)
00516 for(Index i = 0; i < maxi; ++i)
00517 dst.coeffRef(i, j) = static_cast<typename Derived1::Scalar>(0);
00518 }
00519 }
00520 };
00521
00522 template<typename Derived1, typename Derived2, bool ClearOpposite>
00523 struct triangular_assignment_selector<Derived1, Derived2, StrictlyUpper, Dynamic, ClearOpposite>
00524 {
00525 typedef typename Derived1::Index Index;
00526 inline static void run(Derived1 &dst, const Derived2 &src)
00527 {
00528 for(Index j = 0; j < dst.cols(); ++j)
00529 {
00530 Index maxi = std::min(j, dst.rows());
00531 for(Index i = 0; i < maxi; ++i)
00532 dst.copyCoeff(i, j, src);
00533 if (ClearOpposite)
00534 for(Index i = maxi; i < dst.rows(); ++i)
00535 dst.coeffRef(i, j) = 0;
00536 }
00537 }
00538 };
00539
00540 template<typename Derived1, typename Derived2, bool ClearOpposite>
00541 struct triangular_assignment_selector<Derived1, Derived2, StrictlyLower, Dynamic, ClearOpposite>
00542 {
00543 typedef typename Derived1::Index Index;
00544 inline static void run(Derived1 &dst, const Derived2 &src)
00545 {
00546 for(Index j = 0; j < dst.cols(); ++j)
00547 {
00548 for(Index i = j+1; i < dst.rows(); ++i)
00549 dst.copyCoeff(i, j, src);
00550 Index maxi = std::min(j, dst.rows()-1);
00551 if (ClearOpposite)
00552 for(Index i = 0; i <= maxi; ++i)
00553 dst.coeffRef(i, j) = static_cast<typename Derived1::Scalar>(0);
00554 }
00555 }
00556 };
00557
00558 template<typename Derived1, typename Derived2, bool ClearOpposite>
00559 struct triangular_assignment_selector<Derived1, Derived2, UnitUpper, Dynamic, ClearOpposite>
00560 {
00561 typedef typename Derived1::Index Index;
00562 inline static void run(Derived1 &dst, const Derived2 &src)
00563 {
00564 for(Index j = 0; j < dst.cols(); ++j)
00565 {
00566 Index maxi = std::min(j, dst.rows());
00567 for(Index i = 0; i < maxi; ++i)
00568 dst.copyCoeff(i, j, src);
00569 if (ClearOpposite)
00570 {
00571 for(Index i = maxi+1; i < dst.rows(); ++i)
00572 dst.coeffRef(i, j) = 0;
00573 }
00574 }
00575 dst.diagonal().setOnes();
00576 }
00577 };
00578 template<typename Derived1, typename Derived2, bool ClearOpposite>
00579 struct triangular_assignment_selector<Derived1, Derived2, UnitLower, Dynamic, ClearOpposite>
00580 {
00581 typedef typename Derived1::Index Index;
00582 inline static void run(Derived1 &dst, const Derived2 &src)
00583 {
00584 for(Index j = 0; j < dst.cols(); ++j)
00585 {
00586 Index maxi = std::min(j, dst.rows());
00587 for(Index i = maxi+1; i < dst.rows(); ++i)
00588 dst.copyCoeff(i, j, src);
00589 if (ClearOpposite)
00590 {
00591 for(Index i = 0; i < maxi; ++i)
00592 dst.coeffRef(i, j) = 0;
00593 }
00594 }
00595 dst.diagonal().setOnes();
00596 }
00597 };
00598
00599 }
00600
00601
00602 template<typename MatrixType, unsigned int Mode>
00603 template<typename OtherDerived>
00604 inline TriangularView<MatrixType, Mode>&
00605 TriangularView<MatrixType, Mode>::operator=(const MatrixBase<OtherDerived>& other)
00606 {
00607 if(OtherDerived::Flags & EvalBeforeAssigningBit)
00608 {
00609 typename internal::plain_matrix_type<OtherDerived>::type other_evaluated(other.rows(), other.cols());
00610 other_evaluated.template triangularView<Mode>().lazyAssign(other.derived());
00611 lazyAssign(other_evaluated);
00612 }
00613 else
00614 lazyAssign(other.derived());
00615 return *this;
00616 }
00617
00618
00619 template<typename MatrixType, unsigned int Mode>
00620 template<typename OtherDerived>
00621 void TriangularView<MatrixType, Mode>::lazyAssign(const MatrixBase<OtherDerived>& other)
00622 {
00623 enum {
00624 unroll = MatrixType::SizeAtCompileTime != Dynamic
00625 && internal::traits<OtherDerived>::CoeffReadCost != Dynamic
00626 && MatrixType::SizeAtCompileTime*internal::traits<OtherDerived>::CoeffReadCost/2 <= EIGEN_UNROLLING_LIMIT
00627 };
00628 eigen_assert(m_matrix.rows() == other.rows() && m_matrix.cols() == other.cols());
00629
00630 internal::triangular_assignment_selector
00631 <MatrixType, OtherDerived, int(Mode),
00632 unroll ? int(MatrixType::SizeAtCompileTime) : Dynamic,
00633 false
00634 >::run(m_matrix.const_cast_derived(), other.derived());
00635 }
00636
00637
00638
00639 template<typename MatrixType, unsigned int Mode>
00640 template<typename OtherDerived>
00641 inline TriangularView<MatrixType, Mode>&
00642 TriangularView<MatrixType, Mode>::operator=(const TriangularBase<OtherDerived>& other)
00643 {
00644 eigen_assert(Mode == int(OtherDerived::Mode));
00645 if(internal::traits<OtherDerived>::Flags & EvalBeforeAssigningBit)
00646 {
00647 typename OtherDerived::DenseMatrixType other_evaluated(other.rows(), other.cols());
00648 other_evaluated.template triangularView<Mode>().lazyAssign(other.derived().nestedExpression());
00649 lazyAssign(other_evaluated);
00650 }
00651 else
00652 lazyAssign(other.derived().nestedExpression());
00653 return *this;
00654 }
00655
00656 template<typename MatrixType, unsigned int Mode>
00657 template<typename OtherDerived>
00658 void TriangularView<MatrixType, Mode>::lazyAssign(const TriangularBase<OtherDerived>& other)
00659 {
00660 enum {
00661 unroll = MatrixType::SizeAtCompileTime != Dynamic
00662 && internal::traits<OtherDerived>::CoeffReadCost != Dynamic
00663 && MatrixType::SizeAtCompileTime * internal::traits<OtherDerived>::CoeffReadCost / 2
00664 <= EIGEN_UNROLLING_LIMIT
00665 };
00666 eigen_assert(m_matrix.rows() == other.rows() && m_matrix.cols() == other.cols());
00667
00668 internal::triangular_assignment_selector
00669 <MatrixType, OtherDerived, int(Mode),
00670 unroll ? int(MatrixType::SizeAtCompileTime) : Dynamic,
00671 false
00672 >::run(m_matrix.const_cast_derived(), other.derived().nestedExpression());
00673 }
00674
00675
00676
00677
00678
00681 template<typename Derived>
00682 template<typename DenseDerived>
00683 void TriangularBase<Derived>::evalTo(MatrixBase<DenseDerived> &other) const
00684 {
00685 if(internal::traits<Derived>::Flags & EvalBeforeAssigningBit)
00686 {
00687 typename internal::plain_matrix_type<Derived>::type other_evaluated(rows(), cols());
00688 evalToLazy(other_evaluated);
00689 other.derived().swap(other_evaluated);
00690 }
00691 else
00692 evalToLazy(other.derived());
00693 }
00694
00697 template<typename Derived>
00698 template<typename DenseDerived>
00699 void TriangularBase<Derived>::evalToLazy(MatrixBase<DenseDerived> &other) const
00700 {
00701 enum {
00702 unroll = DenseDerived::SizeAtCompileTime != Dynamic
00703 && internal::traits<Derived>::CoeffReadCost != Dynamic
00704 && DenseDerived::SizeAtCompileTime * internal::traits<Derived>::CoeffReadCost / 2
00705 <= EIGEN_UNROLLING_LIMIT
00706 };
00707 other.derived().resize(this->rows(), this->cols());
00708
00709 internal::triangular_assignment_selector
00710 <DenseDerived, typename internal::traits<Derived>::MatrixTypeNestedCleaned, Derived::Mode,
00711 unroll ? int(DenseDerived::SizeAtCompileTime) : Dynamic,
00712 true
00713 >::run(other.derived(), derived().nestedExpression());
00714 }
00715
00716
00717
00718
00719
00720
00721
00722
00723
00724 #ifdef EIGEN2_SUPPORT
00725
00726
00727
00728 namespace internal {
00729 template<typename MatrixType, unsigned int Mode>
00730 struct eigen2_part_return_type
00731 {
00732 typedef TriangularView<MatrixType, Mode> type;
00733 };
00734
00735 template<typename MatrixType>
00736 struct eigen2_part_return_type<MatrixType, SelfAdjoint>
00737 {
00738 typedef SelfAdjointView<MatrixType, Upper> type;
00739 };
00740 }
00741
00743 template<typename Derived>
00744 template<unsigned int Mode>
00745 const typename internal::eigen2_part_return_type<Derived, Mode>::type MatrixBase<Derived>::part() const
00746 {
00747 return derived();
00748 }
00749
00751 template<typename Derived>
00752 template<unsigned int Mode>
00753 typename internal::eigen2_part_return_type<Derived, Mode>::type MatrixBase<Derived>::part()
00754 {
00755 return derived();
00756 }
00757 #endif
00758
00770 template<typename Derived>
00771 template<unsigned int Mode>
00772 typename MatrixBase<Derived>::template TriangularViewReturnType<Mode>::Type
00773 MatrixBase<Derived>::triangularView()
00774 {
00775 return derived();
00776 }
00777
00779 template<typename Derived>
00780 template<unsigned int Mode>
00781 typename MatrixBase<Derived>::template ConstTriangularViewReturnType<Mode>::Type
00782 MatrixBase<Derived>::triangularView() const
00783 {
00784 return derived();
00785 }
00786
00792 template<typename Derived>
00793 bool MatrixBase<Derived>::isUpperTriangular(RealScalar prec) const
00794 {
00795 RealScalar maxAbsOnUpperPart = static_cast<RealScalar>(-1);
00796 for(Index j = 0; j < cols(); ++j)
00797 {
00798 Index maxi = std::min(j, rows()-1);
00799 for(Index i = 0; i <= maxi; ++i)
00800 {
00801 RealScalar absValue = internal::abs(coeff(i,j));
00802 if(absValue > maxAbsOnUpperPart) maxAbsOnUpperPart = absValue;
00803 }
00804 }
00805 RealScalar threshold = maxAbsOnUpperPart * prec;
00806 for(Index j = 0; j < cols(); ++j)
00807 for(Index i = j+1; i < rows(); ++i)
00808 if(internal::abs(coeff(i, j)) > threshold) return false;
00809 return true;
00810 }
00811
00817 template<typename Derived>
00818 bool MatrixBase<Derived>::isLowerTriangular(RealScalar prec) const
00819 {
00820 RealScalar maxAbsOnLowerPart = static_cast<RealScalar>(-1);
00821 for(Index j = 0; j < cols(); ++j)
00822 for(Index i = j; i < rows(); ++i)
00823 {
00824 RealScalar absValue = internal::abs(coeff(i,j));
00825 if(absValue > maxAbsOnLowerPart) maxAbsOnLowerPart = absValue;
00826 }
00827 RealScalar threshold = maxAbsOnLowerPart * prec;
00828 for(Index j = 1; j < cols(); ++j)
00829 {
00830 Index maxi = std::min(j, rows()-1);
00831 for(Index i = 0; i < maxi; ++i)
00832 if(internal::abs(coeff(i, j)) > threshold) return false;
00833 }
00834 return true;
00835 }
00836
00837 #endif // EIGEN_TRIANGULARMATRIX_H