11 #ifndef EIGEN_COLPIVOTINGHOUSEHOLDERQR_H 12 #define EIGEN_COLPIVOTINGHOUSEHOLDERQR_H 17 template<
typename _MatrixType>
struct traits<ColPivHouseholderQR<_MatrixType> >
46 template<
typename _MatrixType>
class ColPivHouseholderQR
50 typedef _MatrixType MatrixType;
52 RowsAtCompileTime = MatrixType::RowsAtCompileTime,
53 ColsAtCompileTime = MatrixType::ColsAtCompileTime,
54 Options = MatrixType::Options,
55 MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
56 MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
58 typedef typename MatrixType::Scalar Scalar;
59 typedef typename MatrixType::RealScalar RealScalar;
61 typedef typename MatrixType::StorageIndex StorageIndex;
62 typedef Matrix<Scalar, RowsAtCompileTime, RowsAtCompileTime, Options, MaxRowsAtCompileTime, MaxRowsAtCompileTime> MatrixQType;
63 typedef typename internal::plain_diag_type<MatrixType>::type HCoeffsType;
64 typedef PermutationMatrix<ColsAtCompileTime, MaxColsAtCompileTime> PermutationType;
65 typedef typename internal::plain_row_type<MatrixType, Index>::type IntRowVectorType;
66 typedef typename internal::plain_row_type<MatrixType>::type RowVectorType;
67 typedef typename internal::plain_row_type<MatrixType, RealScalar>::type RealRowVectorType;
68 typedef HouseholderSequence<MatrixType,typename internal::remove_all<typename HCoeffsType::ConjugateReturnType>::type> HouseholderSequenceType;
69 typedef typename MatrixType::PlainObject PlainObject;
73 typedef typename PermutationType::StorageIndex PermIndexType;
87 m_colsTranspositions(),
90 m_isInitialized(false),
91 m_usePrescribedThreshold(false) {}
101 m_hCoeffs((
std::min)(rows,cols)),
102 m_colsPermutation(PermIndexType(cols)),
103 m_colsTranspositions(cols),
106 m_isInitialized(false),
107 m_usePrescribedThreshold(false) {}
121 template<
typename InputType>
123 : m_qr(matrix.rows(), matrix.cols()),
124 m_hCoeffs((
std::min)(matrix.rows(),matrix.cols())),
125 m_colsPermutation(PermIndexType(matrix.cols())),
126 m_colsTranspositions(matrix.cols()),
127 m_temp(matrix.cols()),
128 m_colSqNorms(matrix.cols()),
129 m_isInitialized(false),
130 m_usePrescribedThreshold(false)
152 template<
typename Rhs>
156 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
160 HouseholderSequenceType householderQ()
const;
161 HouseholderSequenceType matrixQ()
const 163 return householderQ();
170 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
185 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
189 template<
typename InputType>
195 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
196 return m_colsPermutation;
212 typename MatrixType::RealScalar absDeterminant()
const;
226 typename MatrixType::RealScalar logAbsDeterminant()
const;
237 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
238 RealScalar premultiplied_threshold = abs(m_maxpivot) * threshold();
240 for(Index i = 0; i < m_nonzero_pivots; ++i)
241 result += (abs(m_qr.coeff(i,i)) > premultiplied_threshold);
253 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
254 return cols() - rank();
266 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
267 return rank() == cols();
279 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
280 return rank() == rows();
291 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
292 return isInjective() && isSurjective();
302 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
306 inline Index rows()
const {
return m_qr.rows(); }
307 inline Index cols()
const {
return m_qr.cols(); }
313 const HCoeffsType&
hCoeffs()
const {
return m_hCoeffs; }
334 m_usePrescribedThreshold =
true;
335 m_prescribedThreshold = threshold;
349 m_usePrescribedThreshold =
false;
359 eigen_assert(m_isInitialized || m_usePrescribedThreshold);
360 return m_usePrescribedThreshold ? m_prescribedThreshold
375 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
376 return m_nonzero_pivots;
392 eigen_assert(m_isInitialized &&
"Decomposition is not initialized.");
396 #ifndef EIGEN_PARSED_BY_DOXYGEN 397 template<
typename RhsType,
typename DstType>
399 void _solve_impl(
const RhsType &rhs, DstType &dst)
const;
404 static void check_template_parameters()
406 EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar);
409 void computeInPlace();
412 HCoeffsType m_hCoeffs;
413 PermutationType m_colsPermutation;
414 IntRowVectorType m_colsTranspositions;
415 RowVectorType m_temp;
416 RealRowVectorType m_colSqNorms;
417 bool m_isInitialized, m_usePrescribedThreshold;
418 RealScalar m_prescribedThreshold, m_maxpivot;
419 Index m_nonzero_pivots;
423 template<
typename MatrixType>
427 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
428 eigen_assert(m_qr.rows() == m_qr.cols() &&
"You can't take the determinant of a non-square matrix!");
429 return abs(m_qr.diagonal().prod());
432 template<
typename MatrixType>
435 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
436 eigen_assert(m_qr.rows() == m_qr.cols() &&
"You can't take the determinant of a non-square matrix!");
437 return m_qr.diagonal().cwiseAbs().array().log().sum();
446 template<
typename MatrixType>
447 template<
typename InputType>
450 check_template_parameters();
462 template<
typename MatrixType>
466 Index rows = m_qr.rows();
467 Index cols = m_qr.cols();
468 Index size = m_qr.diagonalSize();
470 m_hCoeffs.resize(size);
474 m_colsTranspositions.resize(m_qr.cols());
475 Index number_of_transpositions = 0;
477 m_colSqNorms.resize(cols);
478 for(Index k = 0; k < cols; ++k)
479 m_colSqNorms.coeffRef(k) = m_qr.col(k).squaredNorm();
483 m_nonzero_pivots = size;
484 m_maxpivot = RealScalar(0);
486 for(Index k = 0; k < size; ++k)
489 Index biggest_col_index;
490 RealScalar biggest_col_sq_norm = m_colSqNorms.tail(cols-k).maxCoeff(&biggest_col_index);
491 biggest_col_index += k;
497 biggest_col_sq_norm = m_qr.col(biggest_col_index).tail(rows-k).squaredNorm();
500 m_colSqNorms.coeffRef(biggest_col_index) = biggest_col_sq_norm;
504 if(m_nonzero_pivots==size && biggest_col_sq_norm < threshold_helper * RealScalar(rows-k))
505 m_nonzero_pivots = k;
508 m_colsTranspositions.coeffRef(k) = biggest_col_index;
509 if(k != biggest_col_index) {
510 m_qr.col(k).swap(m_qr.col(biggest_col_index));
511 std::swap(m_colSqNorms.coeffRef(k), m_colSqNorms.coeffRef(biggest_col_index));
512 ++number_of_transpositions;
517 m_qr.col(k).tail(rows-k).makeHouseholderInPlace(m_hCoeffs.coeffRef(k), beta);
520 m_qr.coeffRef(k,k) = beta;
523 if(abs(beta) > m_maxpivot) m_maxpivot = abs(beta);
526 m_qr.bottomRightCorner(rows-k, cols-k-1)
527 .applyHouseholderOnTheLeft(m_qr.col(k).tail(rows-k-1), m_hCoeffs.coeffRef(k), &m_temp.coeffRef(k+1));
530 m_colSqNorms.tail(cols-k-1) -= m_qr.row(k).tail(cols-k-1).cwiseAbs2();
533 m_colsPermutation.setIdentity(PermIndexType(cols));
534 for(PermIndexType k = 0; k < size; ++k)
535 m_colsPermutation.applyTranspositionOnTheRight(k, PermIndexType(m_colsTranspositions.coeff(k)));
537 m_det_pq = (number_of_transpositions%2) ? -1 : 1;
538 m_isInitialized =
true;
541 #ifndef EIGEN_PARSED_BY_DOXYGEN 542 template<
typename _MatrixType>
543 template<
typename RhsType,
typename DstType>
546 eigen_assert(rhs.rows() == rows());
548 const Index nonzero_pivots = nonzeroPivots();
550 if(nonzero_pivots == 0)
556 typename RhsType::PlainObject c(rhs);
560 .setLength(nonzero_pivots)
564 m_qr.topLeftCorner(nonzero_pivots, nonzero_pivots)
565 .template triangularView<Upper>()
566 .solveInPlace(c.topRows(nonzero_pivots));
568 for(Index i = 0; i < nonzero_pivots; ++i) dst.row(m_colsPermutation.indices().coeff(i)) = c.row(i);
569 for(Index i = nonzero_pivots; i < cols(); ++i) dst.row(m_colsPermutation.indices().coeff(i)).setZero();
575 template<
typename DstXprType,
typename MatrixType,
typename Scalar>
576 struct Assignment<DstXprType, Inverse<ColPivHouseholderQR<MatrixType> >, internal::assign_op<Scalar>, Dense2Dense, Scalar>
580 static void run(DstXprType &dst,
const SrcXprType &src,
const internal::assign_op<Scalar> &)
582 dst = src.nestedExpression().solve(MatrixType::Identity(src.rows(), src.cols()));
591 template<
typename MatrixType>
595 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
604 template<
typename Derived>
614 #endif // EIGEN_COLPIVOTINGHOUSEHOLDERQR_H const PermutationType & colsPermutation() const
Definition: ColPivHouseholderQR.h:193
Index nonzeroPivots() const
Definition: ColPivHouseholderQR.h:373
const Inverse< ColPivHouseholderQR > inverse() const
Definition: ColPivHouseholderQR.h:300
bool isInjective() const
Definition: ColPivHouseholderQR.h:264
ColPivHouseholderQR()
Default Constructor.
Definition: ColPivHouseholderQR.h:83
ColPivHouseholderQR(Index rows, Index cols)
Default Constructor with memory preallocation.
Definition: ColPivHouseholderQR.h:99
const Solve< ColPivHouseholderQR, Rhs > solve(const MatrixBase< Rhs > &b) const
Definition: ColPivHouseholderQR.h:154
MatrixType::RealScalar logAbsDeterminant() const
Definition: ColPivHouseholderQR.h:433
Index rank() const
Definition: ColPivHouseholderQR.h:234
bool isInvertible() const
Definition: ColPivHouseholderQR.h:289
Definition: StdDeque.h:58
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:107
Derived & derived()
Definition: EigenBase.h:44
HouseholderSequence< VectorsType, CoeffsType > householderSequence(const VectorsType &v, const CoeffsType &h)
Convenience function for constructing a Householder sequence.
Definition: HouseholderSequence.h:452
ComputationInfo info() const
Reports whether the QR factorization was succesful.
Definition: ColPivHouseholderQR.h:390
Definition: EigenBase.h:28
Sequence of Householder reflections acting on subspaces with decreasing size.
Definition: ForwardDeclarations.h:259
ColPivHouseholderQR(const EigenBase< InputType > &matrix)
Constructs a QR factorization from a given matrix.
Definition: ColPivHouseholderQR.h:122
Expression of the inverse of another expression.
Definition: Inverse.h:43
ColPivHouseholderQR & setThreshold(Default_t)
Definition: ColPivHouseholderQR.h:347
RealScalar threshold() const
Definition: ColPivHouseholderQR.h:357
const MatrixType & matrixR() const
Definition: ColPivHouseholderQR.h:183
Householder rank-revealing QR decomposition of a matrix with column-pivoting.
Definition: ForwardDeclarations.h:253
const MatrixType & matrixQR() const
Definition: ColPivHouseholderQR.h:168
HouseholderSequenceType householderQ() const
Definition: ColPivHouseholderQR.h:593
const HCoeffsType & hCoeffs() const
Definition: ColPivHouseholderQR.h:313
const ColPivHouseholderQR< PlainObject > colPivHouseholderQr() const
Definition: ColPivHouseholderQR.h:606
Definition: Constants.h:432
ColPivHouseholderQR & setThreshold(const RealScalar &threshold)
Definition: ColPivHouseholderQR.h:332
RealScalar maxPivot() const
Definition: ColPivHouseholderQR.h:382
Definition: Eigen_Colamd.h:54
MatrixType::RealScalar absDeterminant() const
Definition: ColPivHouseholderQR.h:424
bool isSurjective() const
Definition: ColPivHouseholderQR.h:277
Index cols() const
Definition: EigenBase.h:61
Pseudo expression representing a solving operation.
Definition: Solve.h:62
ComputationInfo
Definition: Constants.h:430
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
Index dimensionOfKernel() const
Definition: ColPivHouseholderQR.h:251