42 template<
typename _MatrixType>
class HouseholderQR
46 typedef _MatrixType MatrixType;
48 RowsAtCompileTime = MatrixType::RowsAtCompileTime,
49 ColsAtCompileTime = MatrixType::ColsAtCompileTime,
50 Options = MatrixType::Options,
51 MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
52 MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime
54 typedef typename MatrixType::Scalar Scalar;
55 typedef typename MatrixType::RealScalar RealScalar;
56 typedef typename MatrixType::Index Index;
57 typedef Matrix<Scalar, RowsAtCompileTime, RowsAtCompileTime, (MatrixType::Flags&RowMajorBit) ? RowMajor : ColMajor, MaxRowsAtCompileTime, MaxRowsAtCompileTime> MatrixQType;
58 typedef typename internal::plain_diag_type<MatrixType>::type HCoeffsType;
59 typedef typename internal::plain_row_type<MatrixType>::type RowVectorType;
60 typedef HouseholderSequence<MatrixType,typename internal::remove_all<typename HCoeffsType::ConjugateReturnType>::type> HouseholderSequenceType;
68 HouseholderQR() : m_qr(), m_hCoeffs(), m_temp(), m_isInitialized(false) {}
78 m_hCoeffs((
std::min)(rows,cols)),
80 m_isInitialized(false) {}
95 : m_qr(matrix.rows(), matrix.cols()),
96 m_hCoeffs((
std::min)(matrix.rows(),matrix.cols())),
97 m_temp(matrix.cols()),
98 m_isInitialized(false)
120 template<
typename Rhs>
121 inline const internal::solve_retval<HouseholderQR, Rhs>
124 eigen_assert(m_isInitialized &&
"HouseholderQR is not initialized.");
125 return internal::solve_retval<HouseholderQR, Rhs>(*
this, b.derived());
138 eigen_assert(m_isInitialized &&
"HouseholderQR is not initialized.");
139 return HouseholderSequenceType(m_qr, m_hCoeffs.conjugate());
147 eigen_assert(m_isInitialized &&
"HouseholderQR is not initialized.");
182 inline Index rows()
const {
return m_qr.rows(); }
183 inline Index cols()
const {
return m_qr.cols(); }
189 const HCoeffsType&
hCoeffs()
const {
return m_hCoeffs; }
193 HCoeffsType m_hCoeffs;
194 RowVectorType m_temp;
195 bool m_isInitialized;
198 template<
typename MatrixType>
202 eigen_assert(m_isInitialized &&
"HouseholderQR is not initialized.");
203 eigen_assert(m_qr.rows() == m_qr.cols() &&
"You can't take the determinant of a non-square matrix!");
204 return abs(m_qr.diagonal().prod());
207 template<
typename MatrixType>
210 eigen_assert(m_isInitialized &&
"HouseholderQR is not initialized.");
211 eigen_assert(m_qr.rows() == m_qr.cols() &&
"You can't take the determinant of a non-square matrix!");
212 return m_qr.diagonal().cwiseAbs().array().log().sum();
218 template<
typename MatrixQR,
typename HCoeffs>
219 void householder_qr_inplace_unblocked(MatrixQR& mat, HCoeffs&
hCoeffs,
typename MatrixQR::Scalar* tempData = 0)
221 typedef typename MatrixQR::Index Index;
222 typedef typename MatrixQR::Scalar Scalar;
223 typedef typename MatrixQR::RealScalar RealScalar;
224 Index rows = mat.rows();
225 Index cols = mat.cols();
226 Index size = (std::min)(rows,cols);
228 eigen_assert(hCoeffs.size() == size);
235 tempData = tempVector.data();
238 for(Index k = 0; k < size; ++k)
240 Index remainingRows = rows - k;
241 Index remainingCols = cols - k - 1;
244 mat.col(k).tail(remainingRows).makeHouseholderInPlace(hCoeffs.coeffRef(k), beta);
245 mat.coeffRef(k,k) = beta;
248 mat.bottomRightCorner(remainingRows, remainingCols)
249 .applyHouseholderOnTheLeft(mat.col(k).tail(remainingRows-1), hCoeffs.coeffRef(k), tempData+k+1);
254 template<
typename MatrixQR,
typename HCoeffs>
255 void householder_qr_inplace_blocked(MatrixQR& mat, HCoeffs& hCoeffs,
256 typename MatrixQR::Index maxBlockSize=32,
257 typename MatrixQR::Scalar* tempData = 0)
259 typedef typename MatrixQR::Index Index;
260 typedef typename MatrixQR::Scalar Scalar;
263 Index rows = mat.rows();
264 Index cols = mat.cols();
265 Index size = (std::min)(rows, cols);
272 tempData = tempVector.data();
275 Index blockSize = (std::min)(maxBlockSize,size);
278 for (k = 0; k < size; k += blockSize)
280 Index bs = (std::min)(size-k,blockSize);
281 Index tcols = cols - k - bs;
282 Index brows = rows-k;
292 BlockType A11_21 = mat.block(k,k,brows,bs);
295 householder_qr_inplace_unblocked(A11_21, hCoeffsSegment, tempData);
299 BlockType A21_22 = mat.block(k,k+bs,brows,tcols);
300 apply_block_householder_on_the_left(A21_22,A11_21,hCoeffsSegment.adjoint());
305 template<
typename _MatrixType,
typename Rhs>
306 struct solve_retval<HouseholderQR<_MatrixType>, Rhs>
307 : solve_retval_base<HouseholderQR<_MatrixType>, Rhs>
311 template<
typename Dest>
void evalTo(Dest& dst)
const 313 const Index rows = dec().rows(), cols = dec().cols();
314 const Index rank = (std::min)(rows, cols);
315 eigen_assert(rhs().rows() == rows);
317 typename Rhs::PlainObject c(rhs());
322 dec().
hCoeffs().head(rank)).transpose()
326 .topLeftCorner(rank, rank)
327 .template triangularView<Upper>()
328 .solveInPlace(c.topRows(rank));
330 dst.topRows(rank) = c.topRows(rank);
331 dst.bottomRows(cols-rank).setZero();
343 template<
typename MatrixType>
346 Index rows = matrix.rows();
347 Index cols = matrix.cols();
348 Index size = (std::min)(rows,cols);
351 m_hCoeffs.resize(size);
355 internal::householder_qr_inplace_blocked(m_qr, m_hCoeffs, 48, m_temp.data());
357 m_isInitialized =
true;
365 template<
typename Derived>
HouseholderQR(const MatrixType &matrix)
Constructs a QR factorization from a given matrix.
Definition: HouseholderQR.h:94
const MatrixType & matrixQR() const
Definition: HouseholderQR.h:145
Definition: StdDeque.h:58
const internal::solve_retval< HouseholderQR, Rhs > solve(const MatrixBase< Rhs > &b) const
Definition: HouseholderQR.h:122
HouseholderSequence< VectorsType, CoeffsType > householderSequence(const VectorsType &v, const CoeffsType &h)
Convenience function for constructing a Householder sequence.
Definition: HouseholderSequence.h:422
HouseholderSequenceType householderQ() const
Definition: HouseholderQR.h:136
MatrixType::RealScalar absDeterminant() const
Definition: HouseholderQR.h:199
HouseholderQR()
Default Constructor.
Definition: HouseholderQR.h:68
const HCoeffsType & hCoeffs() const
Definition: HouseholderQR.h:189
HouseholderQR & compute(const MatrixType &matrix)
Definition: HouseholderQR.h:344
Definition: Eigen_Colamd.h:54
MatrixType::RealScalar logAbsDeterminant() const
Definition: HouseholderQR.h:208
Expression of a fixed-size or dynamic-size block.
Definition: Block.h:102
const HouseholderQR< PlainObject > householderQr() const
Definition: HouseholderQR.h:367
Householder QR decomposition of a matrix.
Definition: ForwardDeclarations.h:221
HouseholderQR(Index rows, Index cols)
Default Constructor with memory preallocation.
Definition: HouseholderQR.h:76
void resize(Index nbRows, Index nbCols)
Definition: PlainObjectBase.h:235
The matrix class, also used for vectors and row-vectors.
Definition: Matrix.h:127
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48