11 #ifndef EIGEN_COLPIVOTINGHOUSEHOLDERQR_H 12 #define EIGEN_COLPIVOTINGHOUSEHOLDERQR_H 17 template<
typename _MatrixType>
struct traits<ColPivHouseholderQR<_MatrixType> >
48 template<
typename _MatrixType>
class ColPivHouseholderQR
52 typedef _MatrixType MatrixType;
54 RowsAtCompileTime = MatrixType::RowsAtCompileTime,
55 ColsAtCompileTime = MatrixType::ColsAtCompileTime,
56 MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
57 MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
59 typedef typename MatrixType::Scalar Scalar;
60 typedef typename MatrixType::RealScalar RealScalar;
62 typedef typename MatrixType::StorageIndex StorageIndex;
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(),
91 m_isInitialized(false),
92 m_usePrescribedThreshold(false) {}
102 m_hCoeffs((
std::min)(rows,cols)),
103 m_colsPermutation(PermIndexType(cols)),
104 m_colsTranspositions(cols),
106 m_colNormsUpdated(cols),
107 m_colNormsDirect(cols),
108 m_isInitialized(false),
109 m_usePrescribedThreshold(false) {}
123 template<
typename InputType>
125 : m_qr(matrix.rows(), matrix.cols()),
126 m_hCoeffs((
std::min)(matrix.rows(),matrix.cols())),
127 m_colsPermutation(PermIndexType(matrix.cols())),
128 m_colsTranspositions(matrix.cols()),
129 m_temp(matrix.cols()),
130 m_colNormsUpdated(matrix.cols()),
131 m_colNormsDirect(matrix.cols()),
132 m_isInitialized(false),
133 m_usePrescribedThreshold(false)
144 template<
typename InputType>
146 : m_qr(matrix.derived()),
147 m_hCoeffs((
std::min)(matrix.rows(),matrix.cols())),
148 m_colsPermutation(PermIndexType(matrix.cols())),
149 m_colsTranspositions(matrix.cols()),
150 m_temp(matrix.cols()),
151 m_colNormsUpdated(matrix.cols()),
152 m_colNormsDirect(matrix.cols()),
153 m_isInitialized(false),
154 m_usePrescribedThreshold(false)
173 template<
typename Rhs>
177 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
181 HouseholderSequenceType householderQ()
const;
182 HouseholderSequenceType matrixQ()
const 184 return householderQ();
191 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
206 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
210 template<
typename InputType>
216 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
217 return m_colsPermutation;
233 typename MatrixType::RealScalar absDeterminant()
const;
247 typename MatrixType::RealScalar logAbsDeterminant()
const;
258 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
259 RealScalar premultiplied_threshold =
abs(m_maxpivot) * threshold();
261 for(
Index i = 0; i < m_nonzero_pivots; ++i)
262 result += (
abs(m_qr.coeff(i,i)) > premultiplied_threshold);
274 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
275 return cols() - rank();
287 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
288 return rank() == cols();
300 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
301 return rank() == rows();
312 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
313 return isInjective() && isSurjective();
323 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
327 inline Index rows()
const {
return m_qr.rows(); }
328 inline Index cols()
const {
return m_qr.cols(); }
334 const HCoeffsType&
hCoeffs()
const {
return m_hCoeffs; }
355 m_usePrescribedThreshold =
true;
356 m_prescribedThreshold = threshold;
370 m_usePrescribedThreshold =
false;
380 eigen_assert(m_isInitialized || m_usePrescribedThreshold);
381 return m_usePrescribedThreshold ? m_prescribedThreshold
396 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
397 return m_nonzero_pivots;
413 eigen_assert(m_isInitialized &&
"Decomposition is not initialized.");
417 #ifndef EIGEN_PARSED_BY_DOXYGEN 418 template<
typename RhsType,
typename DstType>
420 void _solve_impl(
const RhsType &rhs, DstType &dst)
const;
427 static void check_template_parameters()
429 EIGEN_STATIC_ASSERT_NON_INTEGER(Scalar);
432 void computeInPlace();
435 HCoeffsType m_hCoeffs;
436 PermutationType m_colsPermutation;
437 IntRowVectorType m_colsTranspositions;
438 RowVectorType m_temp;
439 RealRowVectorType m_colNormsUpdated;
440 RealRowVectorType m_colNormsDirect;
441 bool m_isInitialized, m_usePrescribedThreshold;
442 RealScalar m_prescribedThreshold, m_maxpivot;
443 Index m_nonzero_pivots;
447 template<
typename MatrixType>
451 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
452 eigen_assert(m_qr.rows() == m_qr.cols() &&
"You can't take the determinant of a non-square matrix!");
453 return abs(m_qr.diagonal().prod());
456 template<
typename MatrixType>
459 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
460 eigen_assert(m_qr.rows() == m_qr.cols() &&
"You can't take the determinant of a non-square matrix!");
461 return m_qr.diagonal().cwiseAbs().array().log().sum();
470 template<
typename MatrixType>
471 template<
typename InputType>
479 template<
typename MatrixType>
482 check_template_parameters();
489 Index rows = m_qr.rows();
490 Index cols = m_qr.cols();
491 Index size = m_qr.diagonalSize();
493 m_hCoeffs.resize(size);
497 m_colsTranspositions.resize(m_qr.cols());
498 Index number_of_transpositions = 0;
500 m_colNormsUpdated.resize(cols);
501 m_colNormsDirect.resize(cols);
502 for (
Index k = 0; k < cols; ++k) {
505 m_colNormsDirect.coeffRef(k) = m_qr.col(k).norm();
506 m_colNormsUpdated.coeffRef(k) = m_colNormsDirect.coeffRef(k);
512 m_nonzero_pivots = size;
513 m_maxpivot = RealScalar(0);
515 for(
Index k = 0; k < size; ++k)
518 Index biggest_col_index;
519 RealScalar biggest_col_sq_norm = numext::abs2(m_colNormsUpdated.tail(cols-k).maxCoeff(&biggest_col_index));
520 biggest_col_index += k;
524 if(m_nonzero_pivots==size && biggest_col_sq_norm < threshold_helper * RealScalar(rows-k))
525 m_nonzero_pivots = k;
528 m_colsTranspositions.coeffRef(k) = biggest_col_index;
529 if(k != biggest_col_index) {
530 m_qr.col(k).swap(m_qr.col(biggest_col_index));
531 std::swap(m_colNormsUpdated.coeffRef(k), m_colNormsUpdated.coeffRef(biggest_col_index));
532 std::swap(m_colNormsDirect.coeffRef(k), m_colNormsDirect.coeffRef(biggest_col_index));
533 ++number_of_transpositions;
538 m_qr.col(k).tail(rows-k).makeHouseholderInPlace(m_hCoeffs.coeffRef(k), beta);
541 m_qr.coeffRef(k,k) = beta;
544 if(
abs(beta) > m_maxpivot) m_maxpivot =
abs(beta);
547 m_qr.bottomRightCorner(rows-k, cols-k-1)
548 .applyHouseholderOnTheLeft(m_qr.col(k).tail(rows-k-1), m_hCoeffs.coeffRef(k), &m_temp.coeffRef(k+1));
551 for (
Index j = k + 1; j < cols; ++j) {
556 if (m_colNormsUpdated.coeffRef(j) != RealScalar(0)) {
557 RealScalar temp =
abs(m_qr.coeffRef(k, j)) / m_colNormsUpdated.coeffRef(j);
558 temp = (RealScalar(1) + temp) * (RealScalar(1) - temp);
559 temp = temp < RealScalar(0) ? RealScalar(0) : temp;
560 RealScalar temp2 = temp * numext::abs2<RealScalar>(m_colNormsUpdated.coeffRef(j) /
561 m_colNormsDirect.coeffRef(j));
562 if (temp2 <= norm_downdate_threshold) {
565 m_colNormsDirect.coeffRef(j) = m_qr.col(j).tail(rows - k - 1).norm();
566 m_colNormsUpdated.coeffRef(j) = m_colNormsDirect.coeffRef(j);
568 m_colNormsUpdated.coeffRef(j) *= numext::sqrt(temp);
574 m_colsPermutation.setIdentity(PermIndexType(cols));
575 for(PermIndexType k = 0; k < size; ++k)
576 m_colsPermutation.applyTranspositionOnTheRight(k, PermIndexType(m_colsTranspositions.coeff(k)));
578 m_det_pq = (number_of_transpositions%2) ? -1 : 1;
579 m_isInitialized =
true;
582 #ifndef EIGEN_PARSED_BY_DOXYGEN 583 template<
typename _MatrixType>
584 template<
typename RhsType,
typename DstType>
587 eigen_assert(rhs.rows() == rows());
589 const Index nonzero_pivots = nonzeroPivots();
591 if(nonzero_pivots == 0)
597 typename RhsType::PlainObject c(rhs);
601 .setLength(nonzero_pivots)
605 m_qr.topLeftCorner(nonzero_pivots, nonzero_pivots)
606 .template triangularView<Upper>()
607 .solveInPlace(c.topRows(nonzero_pivots));
609 for(
Index i = 0; i < nonzero_pivots; ++i) dst.row(m_colsPermutation.indices().coeff(i)) = c.row(i);
610 for(
Index i = nonzero_pivots; i < cols(); ++i) dst.row(m_colsPermutation.indices().coeff(i)).setZero();
616 template<
typename DstXprType,
typename MatrixType>
617 struct Assignment<DstXprType, Inverse<ColPivHouseholderQR<MatrixType> >, internal::assign_op<typename DstXprType::Scalar,typename ColPivHouseholderQR<MatrixType>::Scalar>, Dense2Dense>
621 static void run(DstXprType &dst,
const SrcXprType &src,
const internal::assign_op<typename DstXprType::Scalar,typename QrType::Scalar> &)
623 dst = src.nestedExpression().solve(MatrixType::Identity(src.rows(), src.cols()));
632 template<
typename MatrixType>
636 eigen_assert(m_isInitialized &&
"ColPivHouseholderQR is not initialized.");
644 template<
typename Derived>
653 #endif // EIGEN_COLPIVOTINGHOUSEHOLDERQR_H Index rank() const
Definition: ColPivHouseholderQR.h:255
bool isSurjective() const
Definition: ColPivHouseholderQR.h:298
ColPivHouseholderQR(EigenBase< InputType > &matrix)
Constructs a QR factorization from a given matrix.
Definition: ColPivHouseholderQR.h:145
ColPivHouseholderQR()
Default Constructor.
Definition: ColPivHouseholderQR.h:83
ColPivHouseholderQR(Index rows, Index cols)
Default Constructor with memory preallocation.
Definition: ColPivHouseholderQR.h:100
const HCoeffsType & hCoeffs() const
Definition: ColPivHouseholderQR.h:334
const ColPivHouseholderQR< PlainObject > colPivHouseholderQr() const
Definition: ColPivHouseholderQR.h:646
Namespace containing all symbols from the Eigen library.
Definition: Core:287
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:150
ComputationInfo info() const
Reports whether the QR factorization was succesful.
Definition: ColPivHouseholderQR.h:411
Derived & derived()
Definition: EigenBase.h:45
Complete orthogonal decomposition (COD) of a matrix.
Definition: ForwardDeclarations.h:257
HouseholderSequence< VectorsType, CoeffsType > householderSequence(const VectorsType &v, const CoeffsType &h)
Convenience function for constructing a Householder sequence.
Definition: HouseholderSequence.h:451
Definition: EigenBase.h:29
Sequence of Householder reflections acting on subspaces with decreasing size.
Definition: ForwardDeclarations.h:262
ColPivHouseholderQR(const EigenBase< InputType > &matrix)
Constructs a QR factorization from a given matrix.
Definition: ColPivHouseholderQR.h:124
Expression of the inverse of another expression.
Definition: Inverse.h:43
HouseholderSequenceType householderQ() const
Definition: ColPivHouseholderQR.h:634
RealScalar maxPivot() const
Definition: ColPivHouseholderQR.h:403
ColPivHouseholderQR & setThreshold(Default_t)
Definition: ColPivHouseholderQR.h:368
Householder rank-revealing QR decomposition of a matrix with column-pivoting.
Definition: ForwardDeclarations.h:255
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:33
const MatrixType & matrixR() const
Definition: ColPivHouseholderQR.h:204
bool isInjective() const
Definition: ColPivHouseholderQR.h:285
const Solve< ColPivHouseholderQR, Rhs > solve(const MatrixBase< Rhs > &b) const
Definition: ColPivHouseholderQR.h:175
bool isInvertible() const
Definition: ColPivHouseholderQR.h:310
Definition: Constants.h:432
MatrixType::RealScalar absDeterminant() const
Definition: ColPivHouseholderQR.h:448
ColPivHouseholderQR & setThreshold(const RealScalar &threshold)
Definition: ColPivHouseholderQR.h:353
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_abs_op< typename Derived::Scalar >, const Derived > abs(const Eigen::ArrayBase< Derived > &x)
Definition: Eigen_Colamd.h:50
RealScalar threshold() const
Definition: ColPivHouseholderQR.h:378
MatrixType::RealScalar logAbsDeterminant() const
Definition: ColPivHouseholderQR.h:457
Index dimensionOfKernel() const
Definition: ColPivHouseholderQR.h:272
Index nonzeroPivots() const
Definition: ColPivHouseholderQR.h:394
Pseudo expression representing a solving operation.
Definition: Solve.h:62
const MatrixType & matrixQR() const
Definition: ColPivHouseholderQR.h:189
ComputationInfo
Definition: Constants.h:430
const PermutationType & colsPermutation() const
Definition: ColPivHouseholderQR.h:214
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
const Inverse< ColPivHouseholderQR > inverse() const
Definition: ColPivHouseholderQR.h:321