TriangularMatrix.h
Go to the documentation of this file.
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2008 Benoit Jacob <jacob.benoit.1@gmail.com>
5 // Copyright (C) 2008-2009 Gael Guennebaud <gael.guennebaud@inria.fr>
6 //
7 // Eigen is free software; you can redistribute it and/or
8 // modify it under the terms of the GNU Lesser General Public
9 // License as published by the Free Software Foundation; either
10 // version 3 of the License, or (at your option) any later version.
11 //
12 // Alternatively, you can redistribute it and/or
13 // modify it under the terms of the GNU General Public License as
14 // published by the Free Software Foundation; either version 2 of
15 // the License, or (at your option) any later version.
16 //
17 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
18 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
19 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
20 // GNU General Public License for more details.
21 //
22 // You should have received a copy of the GNU Lesser General Public
23 // License and a copy of the GNU General Public License along with
24 // Eigen. If not, see <http://www.gnu.org/licenses/>.
25 
26 #ifndef EIGEN_TRIANGULARMATRIX_H
27 #define EIGEN_TRIANGULARMATRIX_H
28 
29 namespace Eigen {
30 
31 namespace internal {
32 
33 template<int Side, typename TriangularType, typename Rhs> struct triangular_solve_retval;
34 
35 }
36 
44 template<typename Derived> class TriangularBase : public EigenBase<Derived>
45 {
46  public:
47 
48  enum {
49  Mode = internal::traits<Derived>::Mode,
50  CoeffReadCost = internal::traits<Derived>::CoeffReadCost,
51  RowsAtCompileTime = internal::traits<Derived>::RowsAtCompileTime,
52  ColsAtCompileTime = internal::traits<Derived>::ColsAtCompileTime,
53  MaxRowsAtCompileTime = internal::traits<Derived>::MaxRowsAtCompileTime,
54  MaxColsAtCompileTime = internal::traits<Derived>::MaxColsAtCompileTime
55  };
56  typedef typename internal::traits<Derived>::Scalar Scalar;
57  typedef typename internal::traits<Derived>::StorageKind StorageKind;
58  typedef typename internal::traits<Derived>::Index Index;
59  typedef typename internal::traits<Derived>::DenseMatrixType DenseMatrixType;
61 
63 
64  inline Index rows() const { return derived().rows(); }
65  inline Index cols() const { return derived().cols(); }
66  inline Index outerStride() const { return derived().outerStride(); }
67  inline Index innerStride() const { return derived().innerStride(); }
68 
69  inline Scalar coeff(Index row, Index col) const { return derived().coeff(row,col); }
70  inline Scalar& coeffRef(Index row, Index col) { return derived().coeffRef(row,col); }
71 
74  template<typename Other>
76  {
77  derived().coeffRef(row, col) = other.coeff(row, col);
78  }
79 
81  {
82  check_coordinates(row, col);
83  return coeff(row,col);
84  }
86  {
87  check_coordinates(row, col);
88  return coeffRef(row,col);
89  }
90 
91  #ifndef EIGEN_PARSED_BY_DOXYGEN
92  inline const Derived& derived() const { return *static_cast<const Derived*>(this); }
93  inline Derived& derived() { return *static_cast<Derived*>(this); }
94  #endif // not EIGEN_PARSED_BY_DOXYGEN
95 
96  template<typename DenseDerived>
97  void evalTo(MatrixBase<DenseDerived> &other) const;
98  template<typename DenseDerived>
99  void evalToLazy(MatrixBase<DenseDerived> &other) const;
100 
102  {
103  DenseMatrixType res(rows(), cols());
104  evalToLazy(res);
105  return res;
106  }
107 
108  protected:
109 
111  {
114  eigen_assert(col>=0 && col<cols() && row>=0 && row<rows());
115  const int mode = int(Mode) & ~SelfAdjoint;
117  eigen_assert((mode==Upper && col>=row)
118  || (mode==Lower && col<=row)
119  || ((mode==StrictlyUpper || mode==UnitUpper) && col>row)
120  || ((mode==StrictlyLower || mode==UnitLower) && col<row));
121  }
122 
123  #ifdef EIGEN_INTERNAL_DEBUGGING
125  {
126  check_coordinates(row, col);
127  }
128  #else
130  #endif
131 
132 };
133 
151 namespace internal {
152 template<typename MatrixType, unsigned int _Mode>
153 struct traits<TriangularView<MatrixType, _Mode> > : traits<MatrixType>
154 {
155  typedef typename nested<MatrixType>::type MatrixTypeNested;
156  typedef typename remove_reference<MatrixTypeNested>::type MatrixTypeNestedNonRef;
157  typedef typename remove_all<MatrixTypeNested>::type MatrixTypeNestedCleaned;
158  typedef MatrixType ExpressionType;
159  typedef typename MatrixType::PlainObject DenseMatrixType;
160  enum {
161  Mode = _Mode,
162  Flags = (MatrixTypeNestedCleaned::Flags & (HereditaryBits) & (~(PacketAccessBit | DirectAccessBit | LinearAccessBit))) | Mode,
163  CoeffReadCost = MatrixTypeNestedCleaned::CoeffReadCost
164  };
165 };
166 }
167 
168 template<int Mode, bool LhsIsTriangular,
169  typename Lhs, bool LhsIsVector,
170  typename Rhs, bool RhsIsVector>
171 struct TriangularProduct;
172 
173 template<typename _MatrixType, unsigned int _Mode> class TriangularView
174  : public TriangularBase<TriangularView<_MatrixType, _Mode> >
175 {
176  public:
177 
179  typedef typename internal::traits<TriangularView>::Scalar Scalar;
180 
181  typedef _MatrixType MatrixType;
182  typedef typename internal::traits<TriangularView>::DenseMatrixType DenseMatrixType;
184 
185  protected:
186  typedef typename internal::traits<TriangularView>::MatrixTypeNested MatrixTypeNested;
187  typedef typename internal::traits<TriangularView>::MatrixTypeNestedNonRef MatrixTypeNestedNonRef;
188  typedef typename internal::traits<TriangularView>::MatrixTypeNestedCleaned MatrixTypeNestedCleaned;
189 
190  typedef typename internal::remove_all<typename MatrixType::ConjugateReturnType>::type MatrixConjugateReturnType;
191 
192  public:
193  using Base::evalToLazy;
194 
195 
196  typedef typename internal::traits<TriangularView>::StorageKind StorageKind;
197  typedef typename internal::traits<TriangularView>::Index Index;
198 
199  enum {
200  Mode = _Mode,
202  | (Mode & Lower ? Upper : 0)
203  | (Mode & (UnitDiag))
204  | (Mode & (ZeroDiag))
205  };
206 
207  inline TriangularView(const MatrixType& matrix) : m_matrix(matrix)
208  {}
209 
210  inline Index rows() const { return m_matrix.rows(); }
211  inline Index cols() const { return m_matrix.cols(); }
212  inline Index outerStride() const { return m_matrix.outerStride(); }
213  inline Index innerStride() const { return m_matrix.innerStride(); }
214 
216  template<typename Other> TriangularView& operator+=(const DenseBase<Other>& other) { return *this = m_matrix + other.derived(); }
218  template<typename Other> TriangularView& operator-=(const DenseBase<Other>& other) { return *this = m_matrix - other.derived(); }
220  TriangularView& operator*=(const typename internal::traits<MatrixType>::Scalar& other) { return *this = m_matrix * other; }
222  TriangularView& operator/=(const typename internal::traits<MatrixType>::Scalar& other) { return *this = m_matrix / other; }
223 
225  void fill(const Scalar& value) { setConstant(value); }
228  { return *this = MatrixType::Constant(rows(), cols(), value); }
233 
237  inline Scalar coeff(Index row, Index col) const
238  {
240  return m_matrix.coeff(row, col);
241  }
242 
247  {
249  return m_matrix.const_cast_derived().coeffRef(row, col);
250  }
251 
254 
256  template<typename OtherDerived>
258 
259  template<typename OtherDerived>
261 
263  { return *this = other.nestedExpression(); }
264 
265  template<typename OtherDerived>
266  void lazyAssign(const TriangularBase<OtherDerived>& other);
267 
268  template<typename OtherDerived>
269  void lazyAssign(const MatrixBase<OtherDerived>& other);
270 
273  { return m_matrix.conjugate(); }
276  { return m_matrix.conjugate(); }
277 
280  { return m_matrix.adjoint(); }
281 
284  {
286  return m_matrix.const_cast_derived().transpose();
287  }
290  {
291  return m_matrix.transpose();
292  }
293 
295  template<typename OtherDerived>
296  TriangularProduct<Mode,true,MatrixType,false,OtherDerived, OtherDerived::IsVectorAtCompileTime>
298  {
299  return TriangularProduct
300  <Mode,true,MatrixType,false,OtherDerived,OtherDerived::IsVectorAtCompileTime>
301  (m_matrix, rhs.derived());
302  }
303 
305  template<typename OtherDerived> friend
306  TriangularProduct<Mode,false,OtherDerived,OtherDerived::IsVectorAtCompileTime,MatrixType,false>
308  {
309  return TriangularProduct
310  <Mode,false,OtherDerived,OtherDerived::IsVectorAtCompileTime,MatrixType,false>
311  (lhs.derived(),rhs.m_matrix);
312  }
313 
314  #ifdef EIGEN2_SUPPORT
315  template<typename OtherDerived>
316  struct eigen2_product_return_type
317  {
319  typedef typename OtherDerived::PlainObject::DenseType OtherPlainObject;
321  typedef typename ProdRetType::PlainObject type;
322  };
323  template<typename OtherDerived>
324  const typename eigen2_product_return_type<OtherDerived>::type
325  operator*(const EigenBase<OtherDerived>& rhs) const
326  {
327  typename OtherDerived::PlainObject::DenseType rhsPlainObject;
328  rhs.evalTo(rhsPlainObject);
329  return this->toDenseMatrix() * rhsPlainObject;
330  }
331  template<typename OtherMatrixType>
332  bool isApprox(const TriangularView<OtherMatrixType, Mode>& other, typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision()) const
333  {
334  return this->toDenseMatrix().isApprox(other.toDenseMatrix(), precision);
335  }
336  template<typename OtherDerived>
337  bool isApprox(const MatrixBase<OtherDerived>& other, typename NumTraits<Scalar>::Real precision = NumTraits<Scalar>::dummy_precision()) const
338  {
339  return this->toDenseMatrix().isApprox(other, precision);
340  }
341  #endif // EIGEN2_SUPPORT
342 
343  template<int Side, typename Other>
344  inline const internal::triangular_solve_retval<Side,TriangularView, Other>
345  solve(const MatrixBase<Other>& other) const;
346 
347  template<int Side, typename OtherDerived>
348  void solveInPlace(const MatrixBase<OtherDerived>& other) const;
349 
350  template<typename Other>
351  inline const internal::triangular_solve_retval<OnTheLeft,TriangularView, Other>
352  solve(const MatrixBase<Other>& other) const
353  { return solve<OnTheLeft>(other); }
354 
355  template<typename OtherDerived>
356  void solveInPlace(const MatrixBase<OtherDerived>& other) const
357  { return solveInPlace<OnTheLeft>(other); }
358 
360  {
361  EIGEN_STATIC_ASSERT((Mode&UnitDiag)==0,PROGRAMMING_ERROR);
363  }
365  {
366  EIGEN_STATIC_ASSERT((Mode&UnitDiag)==0,PROGRAMMING_ERROR);
368  }
369 
370  template<typename OtherDerived>
371  void swap(TriangularBase<OtherDerived> const & other)
372  {
373  TriangularView<SwapWrapper<MatrixType>,Mode>(const_cast<MatrixType&>(m_matrix)).lazyAssign(other.derived());
374  }
375 
376  template<typename OtherDerived>
377  void swap(MatrixBase<OtherDerived> const & other)
378  {
379  SwapWrapper<MatrixType> swaper(const_cast<MatrixType&>(m_matrix));
380  TriangularView<SwapWrapper<MatrixType>,Mode>(swaper).lazyAssign(other.derived());
381  }
382 
384  {
385  if (Mode & UnitDiag)
386  return 1;
387  else if (Mode & ZeroDiag)
388  return 0;
389  else
390  return m_matrix.diagonal().prod();
391  }
392 
393  // TODO simplify the following:
394  template<typename ProductDerived, typename Lhs, typename Rhs>
396  {
397  setZero();
398  return assignProduct(other,1);
399  }
400 
401  template<typename ProductDerived, typename Lhs, typename Rhs>
403  {
404  return assignProduct(other,1);
405  }
406 
407  template<typename ProductDerived, typename Lhs, typename Rhs>
409  {
410  return assignProduct(other,-1);
411  }
412 
413 
414  template<typename ProductDerived>
416  {
417  setZero();
418  return assignProduct(other,other.alpha());
419  }
420 
421  template<typename ProductDerived>
423  {
424  return assignProduct(other,other.alpha());
425  }
426 
427  template<typename ProductDerived>
429  {
430  return assignProduct(other,-other.alpha());
431  }
432 
433  protected:
434 
435  template<typename ProductDerived, typename Lhs, typename Rhs>
437 
439 };
440 
441 /***************************************************************************
442 * Implementation of triangular evaluation/assignment
443 ***************************************************************************/
444 
445 namespace internal {
446 
447 template<typename Derived1, typename Derived2, unsigned int Mode, int UnrollCount, bool ClearOpposite>
448 struct triangular_assignment_selector
449 {
450  enum {
451  col = (UnrollCount-1) / Derived1::RowsAtCompileTime,
452  row = (UnrollCount-1) % Derived1::RowsAtCompileTime
453  };
454 
455  typedef typename Derived1::Scalar Scalar;
456 
457  static inline void run(Derived1 &dst, const Derived2 &src)
458  {
459  triangular_assignment_selector<Derived1, Derived2, Mode, UnrollCount-1, ClearOpposite>::run(dst, src);
460 
461  eigen_assert( Mode == Upper || Mode == Lower
462  || Mode == StrictlyUpper || Mode == StrictlyLower
463  || Mode == UnitUpper || Mode == UnitLower);
464  if((Mode == Upper && row <= col)
465  || (Mode == Lower && row >= col)
466  || (Mode == StrictlyUpper && row < col)
467  || (Mode == StrictlyLower && row > col)
468  || (Mode == UnitUpper && row < col)
469  || (Mode == UnitLower && row > col))
470  dst.copyCoeff(row, col, src);
471  else if(ClearOpposite)
472  {
473  if (Mode&UnitDiag && row==col)
474  dst.coeffRef(row, col) = Scalar(1);
475  else
476  dst.coeffRef(row, col) = Scalar(0);
477  }
478  }
479 };
480 
481 // prevent buggy user code from causing an infinite recursion
482 template<typename Derived1, typename Derived2, unsigned int Mode, bool ClearOpposite>
483 struct triangular_assignment_selector<Derived1, Derived2, Mode, 0, ClearOpposite>
484 {
485  static inline void run(Derived1 &, const Derived2 &) {}
486 };
487 
488 template<typename Derived1, typename Derived2, bool ClearOpposite>
489 struct triangular_assignment_selector<Derived1, Derived2, Upper, Dynamic, ClearOpposite>
490 {
491  typedef typename Derived1::Index Index;
492  typedef typename Derived1::Scalar Scalar;
493  static inline void run(Derived1 &dst, const Derived2 &src)
494  {
495  for(Index j = 0; j < dst.cols(); ++j)
496  {
497  Index maxi = (std::min)(j, dst.rows()-1);
498  for(Index i = 0; i <= maxi; ++i)
499  dst.copyCoeff(i, j, src);
500  if (ClearOpposite)
501  for(Index i = maxi+1; i < dst.rows(); ++i)
502  dst.coeffRef(i, j) = Scalar(0);
503  }
504  }
505 };
506 
507 template<typename Derived1, typename Derived2, bool ClearOpposite>
508 struct triangular_assignment_selector<Derived1, Derived2, Lower, Dynamic, ClearOpposite>
509 {
510  typedef typename Derived1::Index Index;
511  static inline void run(Derived1 &dst, const Derived2 &src)
512  {
513  for(Index j = 0; j < dst.cols(); ++j)
514  {
515  for(Index i = j; i < dst.rows(); ++i)
516  dst.copyCoeff(i, j, src);
517  Index maxi = (std::min)(j, dst.rows());
518  if (ClearOpposite)
519  for(Index i = 0; i < maxi; ++i)
520  dst.coeffRef(i, j) = static_cast<typename Derived1::Scalar>(0);
521  }
522  }
523 };
524 
525 template<typename Derived1, typename Derived2, bool ClearOpposite>
526 struct triangular_assignment_selector<Derived1, Derived2, StrictlyUpper, Dynamic, ClearOpposite>
527 {
528  typedef typename Derived1::Index Index;
529  static inline void run(Derived1 &dst, const Derived2 &src)
530  {
531  for(Index j = 0; j < dst.cols(); ++j)
532  {
533  Index maxi = (std::min)(j, dst.rows());
534  for(Index i = 0; i < maxi; ++i)
535  dst.copyCoeff(i, j, src);
536  if (ClearOpposite)
537  for(Index i = maxi; i < dst.rows(); ++i)
538  dst.coeffRef(i, j) = 0;
539  }
540  }
541 };
542 
543 template<typename Derived1, typename Derived2, bool ClearOpposite>
544 struct triangular_assignment_selector<Derived1, Derived2, StrictlyLower, Dynamic, ClearOpposite>
545 {
546  typedef typename Derived1::Index Index;
547  static inline void run(Derived1 &dst, const Derived2 &src)
548  {
549  for(Index j = 0; j < dst.cols(); ++j)
550  {
551  for(Index i = j+1; i < dst.rows(); ++i)
552  dst.copyCoeff(i, j, src);
553  Index maxi = (std::min)(j, dst.rows()-1);
554  if (ClearOpposite)
555  for(Index i = 0; i <= maxi; ++i)
556  dst.coeffRef(i, j) = static_cast<typename Derived1::Scalar>(0);
557  }
558  }
559 };
560 
561 template<typename Derived1, typename Derived2, bool ClearOpposite>
562 struct triangular_assignment_selector<Derived1, Derived2, UnitUpper, Dynamic, ClearOpposite>
563 {
564  typedef typename Derived1::Index Index;
565  static inline void run(Derived1 &dst, const Derived2 &src)
566  {
567  for(Index j = 0; j < dst.cols(); ++j)
568  {
569  Index maxi = (std::min)(j, dst.rows());
570  for(Index i = 0; i < maxi; ++i)
571  dst.copyCoeff(i, j, src);
572  if (ClearOpposite)
573  {
574  for(Index i = maxi+1; i < dst.rows(); ++i)
575  dst.coeffRef(i, j) = 0;
576  }
577  }
578  dst.diagonal().setOnes();
579  }
580 };
581 template<typename Derived1, typename Derived2, bool ClearOpposite>
582 struct triangular_assignment_selector<Derived1, Derived2, UnitLower, Dynamic, ClearOpposite>
583 {
584  typedef typename Derived1::Index Index;
585  static inline void run(Derived1 &dst, const Derived2 &src)
586  {
587  for(Index j = 0; j < dst.cols(); ++j)
588  {
589  Index maxi = (std::min)(j, dst.rows());
590  for(Index i = maxi+1; i < dst.rows(); ++i)
591  dst.copyCoeff(i, j, src);
592  if (ClearOpposite)
593  {
594  for(Index i = 0; i < maxi; ++i)
595  dst.coeffRef(i, j) = 0;
596  }
597  }
598  dst.diagonal().setOnes();
599  }
600 };
601 
602 } // end namespace internal
603 
604 // FIXME should we keep that possibility
605 template<typename MatrixType, unsigned int Mode>
606 template<typename OtherDerived>
607 inline TriangularView<MatrixType, Mode>&
609 {
610  if(OtherDerived::Flags & EvalBeforeAssigningBit)
611  {
612  typename internal::plain_matrix_type<OtherDerived>::type other_evaluated(other.rows(), other.cols());
613  other_evaluated.template triangularView<Mode>().lazyAssign(other.derived());
614  lazyAssign(other_evaluated);
615  }
616  else
617  lazyAssign(other.derived());
618  return *this;
619 }
620 
621 // FIXME should we keep that possibility
622 template<typename MatrixType, unsigned int Mode>
623 template<typename OtherDerived>
625 {
626  enum {
627  unroll = MatrixType::SizeAtCompileTime != Dynamic
628  && internal::traits<OtherDerived>::CoeffReadCost != Dynamic
629  && MatrixType::SizeAtCompileTime*internal::traits<OtherDerived>::CoeffReadCost/2 <= EIGEN_UNROLLING_LIMIT
630  };
631  eigen_assert(m_matrix.rows() == other.rows() && m_matrix.cols() == other.cols());
632 
633  internal::triangular_assignment_selector
634  <MatrixType, OtherDerived, int(Mode),
635  unroll ? int(MatrixType::SizeAtCompileTime) : Dynamic,
636  false // do not change the opposite triangular part
637  >::run(m_matrix.const_cast_derived(), other.derived());
638 }
639 
640 
641 
642 template<typename MatrixType, unsigned int Mode>
643 template<typename OtherDerived>
646 {
647  eigen_assert(Mode == int(OtherDerived::Mode));
648  if(internal::traits<OtherDerived>::Flags & EvalBeforeAssigningBit)
649  {
650  typename OtherDerived::DenseMatrixType other_evaluated(other.rows(), other.cols());
651  other_evaluated.template triangularView<Mode>().lazyAssign(other.derived().nestedExpression());
652  lazyAssign(other_evaluated);
653  }
654  else
655  lazyAssign(other.derived().nestedExpression());
656  return *this;
657 }
658 
659 template<typename MatrixType, unsigned int Mode>
660 template<typename OtherDerived>
662 {
663  enum {
664  unroll = MatrixType::SizeAtCompileTime != Dynamic
665  && internal::traits<OtherDerived>::CoeffReadCost != Dynamic
666  && MatrixType::SizeAtCompileTime * internal::traits<OtherDerived>::CoeffReadCost / 2
668  };
669  eigen_assert(m_matrix.rows() == other.rows() && m_matrix.cols() == other.cols());
670 
671  internal::triangular_assignment_selector
672  <MatrixType, OtherDerived, int(Mode),
673  unroll ? int(MatrixType::SizeAtCompileTime) : Dynamic,
674  false // preserve the opposite triangular part
675  >::run(m_matrix.const_cast_derived(), other.derived().nestedExpression());
676 }
677 
678 /***************************************************************************
679 * Implementation of TriangularBase methods
680 ***************************************************************************/
681 
684 template<typename Derived>
685 template<typename DenseDerived>
687 {
688  if(internal::traits<Derived>::Flags & EvalBeforeAssigningBit)
689  {
690  typename internal::plain_matrix_type<Derived>::type other_evaluated(rows(), cols());
691  evalToLazy(other_evaluated);
692  other.derived().swap(other_evaluated);
693  }
694  else
695  evalToLazy(other.derived());
696 }
697 
700 template<typename Derived>
701 template<typename DenseDerived>
703 {
704  enum {
705  unroll = DenseDerived::SizeAtCompileTime != Dynamic
706  && internal::traits<Derived>::CoeffReadCost != Dynamic
707  && DenseDerived::SizeAtCompileTime * internal::traits<Derived>::CoeffReadCost / 2
709  };
710  other.derived().resize(this->rows(), this->cols());
711 
712  internal::triangular_assignment_selector
713  <DenseDerived, typename internal::traits<Derived>::MatrixTypeNestedCleaned, Derived::Mode,
714  unroll ? int(DenseDerived::SizeAtCompileTime) : Dynamic,
715  true // clear the opposite triangular part
716  >::run(other.derived(), derived().nestedExpression());
717 }
718 
719 /***************************************************************************
720 * Implementation of TriangularView methods
721 ***************************************************************************/
722 
723 /***************************************************************************
724 * Implementation of MatrixBase methods
725 ***************************************************************************/
726 
727 #ifdef EIGEN2_SUPPORT
728 
729 // implementation of part<>(), including the SelfAdjoint case.
730 
731 namespace internal {
732 template<typename MatrixType, unsigned int Mode>
733 struct eigen2_part_return_type
734 {
736 };
737 
738 template<typename MatrixType>
739 struct eigen2_part_return_type<MatrixType, SelfAdjoint>
740 {
741  typedef SelfAdjointView<MatrixType, Upper> type;
742 };
743 }
744 
746 template<typename Derived>
747 template<unsigned int Mode>
748 const typename internal::eigen2_part_return_type<Derived, Mode>::type MatrixBase<Derived>::part() const
749 {
750  return derived();
751 }
752 
754 template<typename Derived>
755 template<unsigned int Mode>
756 typename internal::eigen2_part_return_type<Derived, Mode>::type MatrixBase<Derived>::part()
757 {
758  return derived();
759 }
760 #endif
761 
773 template<typename Derived>
774 template<unsigned int Mode>
775 typename MatrixBase<Derived>::template TriangularViewReturnType<Mode>::Type
777 {
778  return derived();
779 }
780 
782 template<typename Derived>
783 template<unsigned int Mode>
786 {
787  return derived();
788 }
789 
795 template<typename Derived>
797 {
798  RealScalar maxAbsOnUpperPart = static_cast<RealScalar>(-1);
799  for(Index j = 0; j < cols(); ++j)
800  {
801  Index maxi = (std::min)(j, rows()-1);
802  for(Index i = 0; i <= maxi; ++i)
803  {
804  RealScalar absValue = internal::abs(coeff(i,j));
805  if(absValue > maxAbsOnUpperPart) maxAbsOnUpperPart = absValue;
806  }
807  }
808  RealScalar threshold = maxAbsOnUpperPart * prec;
809  for(Index j = 0; j < cols(); ++j)
810  for(Index i = j+1; i < rows(); ++i)
811  if(internal::abs(coeff(i, j)) > threshold) return false;
812  return true;
813 }
814 
820 template<typename Derived>
822 {
823  RealScalar maxAbsOnLowerPart = static_cast<RealScalar>(-1);
824  for(Index j = 0; j < cols(); ++j)
825  for(Index i = j; i < rows(); ++i)
826  {
827  RealScalar absValue = internal::abs(coeff(i,j));
828  if(absValue > maxAbsOnLowerPart) maxAbsOnLowerPart = absValue;
829  }
830  RealScalar threshold = maxAbsOnLowerPart * prec;
831  for(Index j = 1; j < cols(); ++j)
832  {
833  Index maxi = (std::min)(j, rows()-1);
834  for(Index i = 0; i < maxi; ++i)
835  if(internal::abs(coeff(i, j)) > threshold) return false;
836  }
837  return true;
838 }
839 
840 } // end namespace Eigen
841 
842 #endif // EIGEN_TRIANGULARMATRIX_H