FullPivLU.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) 2006-2009 Benoit Jacob <jacob.benoit.1@gmail.com>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9 
10 #ifndef EIGEN_LU_H
11 #define EIGEN_LU_H
12 
13 namespace Eigen {
14 
15 namespace internal {
16 template<typename _MatrixType> struct traits<FullPivLU<_MatrixType> >
17  : traits<_MatrixType>
18 {
19  typedef MatrixXpr XprKind;
21  typedef int StorageIndex;
22  enum { Flags = 0 };
23 };
24 
25 } // end namespace internal
26 
60 template<typename _MatrixType> class FullPivLU
61  : public SolverBase<FullPivLU<_MatrixType> >
62 {
63  public:
64  typedef _MatrixType MatrixType;
66  friend class SolverBase<FullPivLU>;
67 
69  enum {
70  MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
71  MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
72  };
77  typedef typename MatrixType::PlainObject PlainObject;
78 
85  FullPivLU();
86 
94 
100  template<typename InputType>
101  explicit FullPivLU(const EigenBase<InputType>& matrix);
102 
109  template<typename InputType>
111 
119  template<typename InputType>
121  m_lu = matrix.derived();
122  computeInPlace();
123  return *this;
124  }
125 
132  inline const MatrixType& matrixLU() const
133  {
134  eigen_assert(m_isInitialized && "LU is not initialized.");
135  return m_lu;
136  }
137 
145  inline Index nonzeroPivots() const
146  {
147  eigen_assert(m_isInitialized && "LU is not initialized.");
148  return m_nonzero_pivots;
149  }
150 
154  RealScalar maxPivot() const { return m_maxpivot; }
155 
160  EIGEN_DEVICE_FUNC inline const PermutationPType& permutationP() const
161  {
162  eigen_assert(m_isInitialized && "LU is not initialized.");
163  return m_p;
164  }
165 
170  inline const PermutationQType& permutationQ() const
171  {
172  eigen_assert(m_isInitialized && "LU is not initialized.");
173  return m_q;
174  }
175 
191  {
192  eigen_assert(m_isInitialized && "LU is not initialized.");
194  }
195 
216  image(const MatrixType& originalMatrix) const
217  {
218  eigen_assert(m_isInitialized && "LU is not initialized.");
219  return internal::image_retval<FullPivLU>(*this, originalMatrix);
220  }
221 
222  #ifdef EIGEN_PARSED_BY_DOXYGEN
223 
242  template<typename Rhs>
243  inline const Solve<FullPivLU, Rhs>
244  solve(const MatrixBase<Rhs>& b) const;
245  #endif
246 
250  inline RealScalar rcond() const
251  {
252  eigen_assert(m_isInitialized && "PartialPivLU is not initialized.");
253  return internal::rcond_estimate_helper(m_l1_norm, *this);
254  }
255 
272 
290  FullPivLU& setThreshold(const RealScalar& threshold)
291  {
292  m_usePrescribedThreshold = true;
293  m_prescribedThreshold = threshold;
294  return *this;
295  }
296 
306  {
307  m_usePrescribedThreshold = false;
308  return *this;
309  }
310 
316  {
317  eigen_assert(m_isInitialized || m_usePrescribedThreshold);
318  return m_usePrescribedThreshold ? m_prescribedThreshold
319  // this formula comes from experimenting (see "LU precision tuning" thread on the list)
320  // and turns out to be identical to Higham's formula used already in LDLt.
321  : NumTraits<Scalar>::epsilon() * RealScalar(m_lu.diagonalSize());
322  }
323 
330  inline Index rank() const
331  {
332  using std::abs;
333  eigen_assert(m_isInitialized && "LU is not initialized.");
334  RealScalar premultiplied_threshold = abs(m_maxpivot) * threshold();
335  Index result = 0;
336  for(Index i = 0; i < m_nonzero_pivots; ++i)
337  result += (abs(m_lu.coeff(i,i)) > premultiplied_threshold);
338  return result;
339  }
340 
347  inline Index dimensionOfKernel() const
348  {
349  eigen_assert(m_isInitialized && "LU is not initialized.");
350  return cols() - rank();
351  }
352 
360  inline bool isInjective() const
361  {
362  eigen_assert(m_isInitialized && "LU is not initialized.");
363  return rank() == cols();
364  }
365 
373  inline bool isSurjective() const
374  {
375  eigen_assert(m_isInitialized && "LU is not initialized.");
376  return rank() == rows();
377  }
378 
385  inline bool isInvertible() const
386  {
387  eigen_assert(m_isInitialized && "LU is not initialized.");
388  return isInjective() && (m_lu.rows() == m_lu.cols());
389  }
390 
398  inline const Inverse<FullPivLU> inverse() const
399  {
400  eigen_assert(m_isInitialized && "LU is not initialized.");
401  eigen_assert(m_lu.rows() == m_lu.cols() && "You can't take the inverse of a non-square matrix!");
402  return Inverse<FullPivLU>(*this);
403  }
404 
405  MatrixType reconstructedMatrix() const;
406 
408  inline Index rows() const EIGEN_NOEXCEPT { return m_lu.rows(); }
410  inline Index cols() const EIGEN_NOEXCEPT { return m_lu.cols(); }
411 
412  #ifndef EIGEN_PARSED_BY_DOXYGEN
413  template<typename RhsType, typename DstType>
414  void _solve_impl(const RhsType &rhs, DstType &dst) const;
415 
416  template<bool Conjugate, typename RhsType, typename DstType>
417  void _solve_impl_transposed(const RhsType &rhs, DstType &dst) const;
418  #endif
419 
420  protected:
421 
423  {
425  }
426 
427  void computeInPlace();
428 
429  MatrixType m_lu;
430  PermutationPType m_p;
431  PermutationQType m_q;
432  IntColVectorType m_rowsTranspositions;
433  IntRowVectorType m_colsTranspositions;
437  signed char m_det_pq;
438  bool m_isInitialized, m_usePrescribedThreshold;
439 };
440 
441 template<typename MatrixType>
443  : m_isInitialized(false), m_usePrescribedThreshold(false)
444 {
445 }
446 
447 template<typename MatrixType>
449  : m_lu(rows, cols),
450  m_p(rows),
451  m_q(cols),
452  m_rowsTranspositions(rows),
453  m_colsTranspositions(cols),
454  m_isInitialized(false),
456 {
457 }
458 
459 template<typename MatrixType>
460 template<typename InputType>
462  : m_lu(matrix.rows(), matrix.cols()),
463  m_p(matrix.rows()),
464  m_q(matrix.cols()),
465  m_rowsTranspositions(matrix.rows()),
466  m_colsTranspositions(matrix.cols()),
467  m_isInitialized(false),
469 {
470  compute(matrix.derived());
471 }
472 
473 template<typename MatrixType>
474 template<typename InputType>
476  : m_lu(matrix.derived()),
477  m_p(matrix.rows()),
478  m_q(matrix.cols()),
479  m_rowsTranspositions(matrix.rows()),
480  m_colsTranspositions(matrix.cols()),
481  m_isInitialized(false),
483 {
484  computeInPlace();
485 }
486 
487 template<typename MatrixType>
489 {
491 
492  // the permutations are stored as int indices, so just to be sure:
494 
495  m_l1_norm = m_lu.cwiseAbs().colwise().sum().maxCoeff();
496 
497  const Index size = m_lu.diagonalSize();
498  const Index rows = m_lu.rows();
499  const Index cols = m_lu.cols();
500 
501  // will store the transpositions, before we accumulate them at the end.
502  // can't accumulate on-the-fly because that will be done in reverse order for the rows.
503  m_rowsTranspositions.resize(m_lu.rows());
504  m_colsTranspositions.resize(m_lu.cols());
505  Index number_of_transpositions = 0; // number of NONTRIVIAL transpositions, i.e. m_rowsTranspositions[i]!=i
506 
507  m_nonzero_pivots = size; // the generic case is that in which all pivots are nonzero (invertible case)
508  m_maxpivot = RealScalar(0);
509 
510  for(Index k = 0; k < size; ++k)
511  {
512  // First, we need to find the pivot.
513 
514  // biggest coefficient in the remaining bottom-right corner (starting at row k, col k)
515  Index row_of_biggest_in_corner, col_of_biggest_in_corner;
517  typedef typename Scoring::result_type Score;
518  Score biggest_in_corner;
519  biggest_in_corner = m_lu.bottomRightCorner(rows-k, cols-k)
520  .unaryExpr(Scoring())
521  .maxCoeff(&row_of_biggest_in_corner, &col_of_biggest_in_corner);
522  row_of_biggest_in_corner += k; // correct the values! since they were computed in the corner,
523  col_of_biggest_in_corner += k; // need to add k to them.
524 
525  if(biggest_in_corner==Score(0))
526  {
527  // before exiting, make sure to initialize the still uninitialized transpositions
528  // in a sane state without destroying what we already have.
529  m_nonzero_pivots = k;
530  for(Index i = k; i < size; ++i)
531  {
532  m_rowsTranspositions.coeffRef(i) = internal::convert_index<StorageIndex>(i);
533  m_colsTranspositions.coeffRef(i) = internal::convert_index<StorageIndex>(i);
534  }
535  break;
536  }
537 
538  RealScalar abs_pivot = internal::abs_knowing_score<Scalar>()(m_lu(row_of_biggest_in_corner, col_of_biggest_in_corner), biggest_in_corner);
539  if(abs_pivot > m_maxpivot) m_maxpivot = abs_pivot;
540 
541  // Now that we've found the pivot, we need to apply the row/col swaps to
542  // bring it to the location (k,k).
543 
544  m_rowsTranspositions.coeffRef(k) = internal::convert_index<StorageIndex>(row_of_biggest_in_corner);
545  m_colsTranspositions.coeffRef(k) = internal::convert_index<StorageIndex>(col_of_biggest_in_corner);
546  if(k != row_of_biggest_in_corner) {
547  m_lu.row(k).swap(m_lu.row(row_of_biggest_in_corner));
548  ++number_of_transpositions;
549  }
550  if(k != col_of_biggest_in_corner) {
551  m_lu.col(k).swap(m_lu.col(col_of_biggest_in_corner));
552  ++number_of_transpositions;
553  }
554 
555  // Now that the pivot is at the right location, we update the remaining
556  // bottom-right corner by Gaussian elimination.
557 
558  if(k<rows-1)
559  m_lu.col(k).tail(rows-k-1) /= m_lu.coeff(k,k);
560  if(k<size-1)
561  m_lu.block(k+1,k+1,rows-k-1,cols-k-1).noalias() -= m_lu.col(k).tail(rows-k-1) * m_lu.row(k).tail(cols-k-1);
562  }
563 
564  // the main loop is over, we still have to accumulate the transpositions to find the
565  // permutations P and Q
566 
567  m_p.setIdentity(rows);
568  for(Index k = size-1; k >= 0; --k)
570 
571  m_q.setIdentity(cols);
572  for(Index k = 0; k < size; ++k)
574 
575  m_det_pq = (number_of_transpositions%2) ? -1 : 1;
576 
577  m_isInitialized = true;
578 }
579 
580 template<typename MatrixType>
582 {
583  eigen_assert(m_isInitialized && "LU is not initialized.");
584  eigen_assert(m_lu.rows() == m_lu.cols() && "You can't take the determinant of a non-square matrix!");
585  return Scalar(m_det_pq) * Scalar(m_lu.diagonal().prod());
586 }
587 
591 template<typename MatrixType>
593 {
594  eigen_assert(m_isInitialized && "LU is not initialized.");
595  const Index smalldim = (std::min)(m_lu.rows(), m_lu.cols());
596  // LU
597  MatrixType res(m_lu.rows(),m_lu.cols());
598  // FIXME the .toDenseMatrix() should not be needed...
599  res = m_lu.leftCols(smalldim)
600  .template triangularView<UnitLower>().toDenseMatrix()
601  * m_lu.topRows(smalldim)
602  .template triangularView<Upper>().toDenseMatrix();
603 
604  // P^{-1}(LU)
605  res = m_p.inverse() * res;
606 
607  // (P^{-1}LU)Q^{-1}
608  res = res * m_q.inverse();
609 
610  return res;
611 }
612 
613 /********* Implementation of kernel() **************************************************/
614 
615 namespace internal {
616 template<typename _MatrixType>
617 struct kernel_retval<FullPivLU<_MatrixType> >
618  : kernel_retval_base<FullPivLU<_MatrixType> >
619 {
621 
622  enum { MaxSmallDimAtCompileTime = EIGEN_SIZE_MIN_PREFER_FIXED(
623  MatrixType::MaxColsAtCompileTime,
624  MatrixType::MaxRowsAtCompileTime)
625  };
626 
627  template<typename Dest> void evalTo(Dest& dst) const
628  {
629  using std::abs;
630  const Index cols = dec().matrixLU().cols(), dimker = cols - rank();
631  if(dimker == 0)
632  {
633  // The Kernel is just {0}, so it doesn't have a basis properly speaking, but let's
634  // avoid crashing/asserting as that depends on floating point calculations. Let's
635  // just return a single column vector filled with zeros.
636  dst.setZero();
637  return;
638  }
639 
640  /* Let us use the following lemma:
641  *
642  * Lemma: If the matrix A has the LU decomposition PAQ = LU,
643  * then Ker A = Q(Ker U).
644  *
645  * Proof: trivial: just keep in mind that P, Q, L are invertible.
646  */
647 
648  /* Thus, all we need to do is to compute Ker U, and then apply Q.
649  *
650  * U is upper triangular, with eigenvalues sorted so that any zeros appear at the end.
651  * Thus, the diagonal of U ends with exactly
652  * dimKer zero's. Let us use that to construct dimKer linearly
653  * independent vectors in Ker U.
654  */
655 
657  RealScalar premultiplied_threshold = dec().maxPivot() * dec().threshold();
658  Index p = 0;
659  for(Index i = 0; i < dec().nonzeroPivots(); ++i)
660  if(abs(dec().matrixLU().coeff(i,i)) > premultiplied_threshold)
661  pivots.coeffRef(p++) = i;
662  eigen_internal_assert(p == rank());
663 
664  // we construct a temporaty trapezoid matrix m, by taking the U matrix and
665  // permuting the rows and cols to bring the nonnegligible pivots to the top of
666  // the main diagonal. We need that to be able to apply our triangular solvers.
667  // FIXME when we get triangularView-for-rectangular-matrices, this can be simplified
668  Matrix<typename MatrixType::Scalar, Dynamic, Dynamic, MatrixType::Options,
669  MaxSmallDimAtCompileTime, MatrixType::MaxColsAtCompileTime>
670  m(dec().matrixLU().block(0, 0, rank(), cols));
671  for(Index i = 0; i < rank(); ++i)
672  {
673  if(i) m.row(i).head(i).setZero();
674  m.row(i).tail(cols-i) = dec().matrixLU().row(pivots.coeff(i)).tail(cols-i);
675  }
676  m.block(0, 0, rank(), rank());
677  m.block(0, 0, rank(), rank()).template triangularView<StrictlyLower>().setZero();
678  for(Index i = 0; i < rank(); ++i)
679  m.col(i).swap(m.col(pivots.coeff(i)));
680 
681  // ok, we have our trapezoid matrix, we can apply the triangular solver.
682  // notice that the math behind this suggests that we should apply this to the
683  // negative of the RHS, but for performance we just put the negative sign elsewhere, see below.
684  m.topLeftCorner(rank(), rank())
685  .template triangularView<Upper>().solveInPlace(
686  m.topRightCorner(rank(), dimker)
687  );
688 
689  // now we must undo the column permutation that we had applied!
690  for(Index i = rank()-1; i >= 0; --i)
691  m.col(i).swap(m.col(pivots.coeff(i)));
692 
693  // see the negative sign in the next line, that's what we were talking about above.
694  for(Index i = 0; i < rank(); ++i) dst.row(dec().permutationQ().indices().coeff(i)) = -m.row(i).tail(dimker);
695  for(Index i = rank(); i < cols; ++i) dst.row(dec().permutationQ().indices().coeff(i)).setZero();
696  for(Index k = 0; k < dimker; ++k) dst.coeffRef(dec().permutationQ().indices().coeff(rank()+k), k) = Scalar(1);
697  }
698 };
699 
700 /***** Implementation of image() *****************************************************/
701 
702 template<typename _MatrixType>
703 struct image_retval<FullPivLU<_MatrixType> >
704  : image_retval_base<FullPivLU<_MatrixType> >
705 {
707 
708  enum { MaxSmallDimAtCompileTime = EIGEN_SIZE_MIN_PREFER_FIXED(
709  MatrixType::MaxColsAtCompileTime,
710  MatrixType::MaxRowsAtCompileTime)
711  };
712 
713  template<typename Dest> void evalTo(Dest& dst) const
714  {
715  using std::abs;
716  if(rank() == 0)
717  {
718  // The Image is just {0}, so it doesn't have a basis properly speaking, but let's
719  // avoid crashing/asserting as that depends on floating point calculations. Let's
720  // just return a single column vector filled with zeros.
721  dst.setZero();
722  return;
723  }
724 
726  RealScalar premultiplied_threshold = dec().maxPivot() * dec().threshold();
727  Index p = 0;
728  for(Index i = 0; i < dec().nonzeroPivots(); ++i)
729  if(abs(dec().matrixLU().coeff(i,i)) > premultiplied_threshold)
730  pivots.coeffRef(p++) = i;
731  eigen_internal_assert(p == rank());
732 
733  for(Index i = 0; i < rank(); ++i)
734  dst.col(i) = originalMatrix().col(dec().permutationQ().indices().coeff(pivots.coeff(i)));
735  }
736 };
737 
738 /***** Implementation of solve() *****************************************************/
739 
740 } // end namespace internal
741 
742 #ifndef EIGEN_PARSED_BY_DOXYGEN
743 template<typename _MatrixType>
744 template<typename RhsType, typename DstType>
745 void FullPivLU<_MatrixType>::_solve_impl(const RhsType &rhs, DstType &dst) const
746 {
747  /* The decomposition PAQ = LU can be rewritten as A = P^{-1} L U Q^{-1}.
748  * So we proceed as follows:
749  * Step 1: compute c = P * rhs.
750  * Step 2: replace c by the solution x to Lx = c. Exists because L is invertible.
751  * Step 3: replace c by the solution x to Ux = c. May or may not exist.
752  * Step 4: result = Q * c;
753  */
754 
755  const Index rows = this->rows(),
756  cols = this->cols(),
757  nonzero_pivots = this->rank();
758  const Index smalldim = (std::min)(rows, cols);
759 
760  if(nonzero_pivots == 0)
761  {
762  dst.setZero();
763  return;
764  }
765 
766  typename RhsType::PlainObject c(rhs.rows(), rhs.cols());
767 
768  // Step 1
769  c = permutationP() * rhs;
770 
771  // Step 2
772  m_lu.topLeftCorner(smalldim,smalldim)
773  .template triangularView<UnitLower>()
774  .solveInPlace(c.topRows(smalldim));
775  if(rows>cols)
776  c.bottomRows(rows-cols) -= m_lu.bottomRows(rows-cols) * c.topRows(cols);
777 
778  // Step 3
779  m_lu.topLeftCorner(nonzero_pivots, nonzero_pivots)
780  .template triangularView<Upper>()
781  .solveInPlace(c.topRows(nonzero_pivots));
782 
783  // Step 4
784  for(Index i = 0; i < nonzero_pivots; ++i)
785  dst.row(permutationQ().indices().coeff(i)) = c.row(i);
786  for(Index i = nonzero_pivots; i < m_lu.cols(); ++i)
787  dst.row(permutationQ().indices().coeff(i)).setZero();
788 }
789 
790 template<typename _MatrixType>
791 template<bool Conjugate, typename RhsType, typename DstType>
792 void FullPivLU<_MatrixType>::_solve_impl_transposed(const RhsType &rhs, DstType &dst) const
793 {
794  /* The decomposition PAQ = LU can be rewritten as A = P^{-1} L U Q^{-1},
795  * and since permutations are real and unitary, we can write this
796  * as A^T = Q U^T L^T P,
797  * So we proceed as follows:
798  * Step 1: compute c = Q^T rhs.
799  * Step 2: replace c by the solution x to U^T x = c. May or may not exist.
800  * Step 3: replace c by the solution x to L^T x = c.
801  * Step 4: result = P^T c.
802  * If Conjugate is true, replace "^T" by "^*" above.
803  */
804 
805  const Index rows = this->rows(), cols = this->cols(),
806  nonzero_pivots = this->rank();
807  const Index smalldim = (std::min)(rows, cols);
808 
809  if(nonzero_pivots == 0)
810  {
811  dst.setZero();
812  return;
813  }
814 
815  typename RhsType::PlainObject c(rhs.rows(), rhs.cols());
816 
817  // Step 1
818  c = permutationQ().inverse() * rhs;
819 
820  // Step 2
821  m_lu.topLeftCorner(nonzero_pivots, nonzero_pivots)
822  .template triangularView<Upper>()
823  .transpose()
824  .template conjugateIf<Conjugate>()
825  .solveInPlace(c.topRows(nonzero_pivots));
826 
827  // Step 3
828  m_lu.topLeftCorner(smalldim, smalldim)
829  .template triangularView<UnitLower>()
830  .transpose()
831  .template conjugateIf<Conjugate>()
832  .solveInPlace(c.topRows(smalldim));
833 
834  // Step 4
835  PermutationPType invp = permutationP().inverse().eval();
836  for(Index i = 0; i < smalldim; ++i)
837  dst.row(invp.indices().coeff(i)) = c.row(i);
838  for(Index i = smalldim; i < rows; ++i)
839  dst.row(invp.indices().coeff(i)).setZero();
840 }
841 
842 #endif
843 
844 namespace internal {
845 
846 
847 /***** Implementation of inverse() *****************************************************/
848 template<typename DstXprType, typename MatrixType>
849 struct Assignment<DstXprType, Inverse<FullPivLU<MatrixType> >, internal::assign_op<typename DstXprType::Scalar,typename FullPivLU<MatrixType>::Scalar>, Dense2Dense>
850 {
853  static void run(DstXprType &dst, const SrcXprType &src, const internal::assign_op<typename DstXprType::Scalar,typename MatrixType::Scalar> &)
854  {
855  dst = src.nestedExpression().solve(MatrixType::Identity(src.rows(), src.cols()));
856  }
857 };
858 } // end namespace internal
859 
860 /******* MatrixBase methods *****************************************************************/
861 
868 template<typename Derived>
871 {
872  return FullPivLU<PlainObject>(eval());
873 }
874 
875 } // end namespace Eigen
876 
877 #endif // EIGEN_LU_H
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
Definition: Inverse.h:57
internal::plain_row_type< MatrixType, StorageIndex >::type IntRowVectorType
Definition: FullPivLU.h:73
PermutationMatrix< ColsAtCompileTime, MaxColsAtCompileTime > PermutationQType
Definition: FullPivLU.h:75
Matrix3f m
RealScalar m_maxpivot
Definition: FullPivLU.h:436
SCALAR Scalar
Definition: bench_gemm.cpp:46
bool m_usePrescribedThreshold
Definition: FullPivLU.h:438
#define EIGEN_GENERIC_PUBLIC_INTERFACE(Derived)
Definition: Macros.h:1264
const internal::image_retval< FullPivLU > image(const MatrixType &originalMatrix) const
Definition: FullPivLU.h:216
SolverBase< FullPivLU > Base
Definition: FullPivLU.h:65
m m block(1, 0, 2, 2)<< 4
Scalar * b
Definition: benchVecAdd.cpp:17
internal::traits< FullPivLU< _MatrixType > >::Scalar Scalar
Definition: SolverBase.h:73
FullPivLU & setThreshold(const RealScalar &threshold)
Definition: FullPivLU.h:290
RealScalar m_prescribedThreshold
Definition: FullPivLU.h:436
InverseReturnType inverse() const
#define min(a, b)
Definition: datatypes.h:19
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
Definition: FullPivLU.h:408
void determinant(const MatrixType &m)
Definition: determinant.cpp:14
FullPivLU & compute(const EigenBase< InputType > &matrix)
Definition: FullPivLU.h:120
Scalar Scalar * c
Definition: benchVecAdd.cpp:17
const PermutationQType & permutationQ() const
Definition: FullPivLU.h:170
_MatrixType MatrixType
Definition: FullPivLU.h:64
Derived & applyTranspositionOnTheRight(Index i, Index j)
RealScalar rcond() const
Definition: FullPivLU.h:250
Namespace containing all symbols from the Eigen library.
Definition: jet.h:637
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE FixedSegmentReturnType< internal::get_fixed_value< NType >::value >::Type tail(NType n)
MatrixXf MatrixType
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:232
Default_t
Definition: Constants.h:362
void _solve_impl(const RhsType &rhs, DstType &dst) const
Definition: FullPivLU.h:745
Decomposition::RealScalar rcond_estimate_helper(typename Decomposition::RealScalar matrix_norm, const Decomposition &dec)
Reciprocal condition number estimator.
RealScalar m_l1_norm
Definition: FullPivLU.h:435
Eigen::Index Index
The interface type of indices.
Definition: EigenBase.h:39
const FullPivLU< PlainObject > fullPivLu() const
Definition: FullPivLU.h:870
PermutationQType m_q
Definition: FullPivLU.h:431
#define EIGEN_SIZE_MIN_PREFER_FIXED(a, b)
Definition: Macros.h:1302
signed char m_det_pq
Definition: FullPivLU.h:437
Index nonzeroPivots() const
Definition: FullPivLU.h:145
static double epsilon
Definition: testRot3.cpp:37
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE Scalar & coeffRef(Index rowId, Index colId)
Index rank() const
Definition: FullPivLU.h:330
Expression of the inverse of another expression.
Definition: Inverse.h:43
cout<< "Here is the matrix m:"<< endl<< m<< endl;Matrix< ptrdiff_t, 3, 1 > res
RealScalar threshold() const
Definition: FullPivLU.h:315
internal::traits< MatrixType >::Scalar determinant() const
Definition: FullPivLU.h:581
const MatrixType & matrixLU() const
Definition: FullPivLU.h:132
Values result
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index size() const EIGEN_NOEXCEPT
Definition: EigenBase.h:67
Index m_nonzero_pivots
Definition: FullPivLU.h:434
IntRowVectorType m_colsTranspositions
Definition: FullPivLU.h:433
ConstTransposeReturnType transpose() const
Definition: SolverBase.h:121
#define EIGEN_NOEXCEPT
Definition: Macros.h:1418
FullPivLU & setThreshold(Default_t)
Definition: FullPivLU.h:305
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const Scalar & coeff(Index rowId, Index colId) const
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:74
const Inverse< FullPivLU > inverse() const
Definition: FullPivLU.h:398
#define eigen_assert(x)
Definition: Macros.h:1037
static void run(DstXprType &dst, const SrcXprType &src, const internal::assign_op< typename DstXprType::Scalar, typename MatrixType::Scalar > &)
Definition: FullPivLU.h:853
static void check_template_parameters()
Definition: FullPivLU.h:422
bool isInjective() const
Definition: FullPivLU.h:360
void _solve_impl_transposed(const RhsType &rhs, DstType &dst) const
Definition: FullPivLU.h:792
#define EIGEN_STATIC_ASSERT_NON_INTEGER(TYPE)
Definition: StaticAssert.h:187
bool isSurjective() const
Definition: FullPivLU.h:373
#define EIGEN_CONSTEXPR
Definition: Macros.h:787
MatrixType reconstructedMatrix() const
Definition: FullPivLU.h:592
PermutationMatrix< RowsAtCompileTime, MaxRowsAtCompileTime > PermutationPType
Definition: FullPivLU.h:76
NumTraits< Scalar >::Real RealScalar
Definition: bench_gemm.cpp:47
bool isInvertible() const
Definition: FullPivLU.h:385
#define EIGEN_DEVICE_FUNC
Definition: Macros.h:976
internal::plain_col_type< MatrixType, StorageIndex >::type IntColVectorType
Definition: FullPivLU.h:74
#define EIGEN_MAKE_KERNEL_HELPERS(DecompositionType)
Definition: Kernel.h:66
const internal::kernel_retval< FullPivLU > kernel() const
Definition: FullPivLU.h:190
PermutationPType m_p
Definition: FullPivLU.h:430
LU decomposition of a matrix with complete pivoting, and related features.
FullPivLU()
Default Constructor.
Definition: FullPivLU.h:442
float * p
IntColVectorType m_rowsTranspositions
Definition: FullPivLU.h:432
#define EIGEN_MAKE_IMAGE_HELPERS(DecompositionType)
Definition: Image.h:67
void computeInPlace()
Definition: FullPivLU.h:488
EIGEN_DEVICE_FUNC const PermutationPType & permutationP() const
Definition: FullPivLU.h:160
internal::nested_eval< T, 1 >::type eval(const T &xpr)
const int Dynamic
Definition: Constants.h:22
Pseudo expression representing a solving operation.
Definition: Solve.h:62
Index dimensionOfKernel() const
Definition: FullPivLU.h:347
#define eigen_internal_assert(x)
Definition: Macros.h:1043
Generic expression where a coefficient-wise unary operator is applied to an expression.
Definition: CwiseUnaryOp.h:55
MatrixType::PlainObject PlainObject
Definition: FullPivLU.h:77
Map< Matrix< T, Dynamic, Dynamic, ColMajor >, 0, OuterStride<> > matrix(T *data, int rows, int cols, int stride)
The matrix class, also used for vectors and row-vectors.
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index cols() const EIGEN_NOEXCEPT
Definition: FullPivLU.h:410
#define abs(x)
Definition: datatypes.h:17
A base class for matrix decomposition and solvers.
Definition: SolverBase.h:68
EIGEN_DEVICE_FUNC Derived & derived()
Definition: EigenBase.h:46
bool m_isInitialized
Definition: FullPivLU.h:438
MatrixType m_lu
Definition: FullPivLU.h:429
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
RealScalar maxPivot() const
Definition: FullPivLU.h:154
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index cols() const EIGEN_NOEXCEPT
Definition: Inverse.h:58
EIGEN_DEVICE_FUNC const XprTypeNestedCleaned & nestedExpression() const
Definition: Inverse.h:60
v setZero(3)


gtsam
Author(s):
autogenerated on Tue Jul 4 2023 02:34:14