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;
58 typedef typename internal::plain_diag_type<MatrixType>::type HCoeffsType;
59 typedef typename internal::plain_row_type<MatrixType>::type RowVectorType;
68 HouseholderQR() : m_qr(), m_hCoeffs(), m_temp(), m_isInitialized(false) {}
78 m_hCoeffs((std::min)(rows,cols)),
80 m_isInitialized(false) {}
83 : m_qr(matrix.rows(), matrix.cols()),
84 m_hCoeffs((std::min)(matrix.rows(),matrix.cols())),
85 m_temp(matrix.cols()),
86 m_isInitialized(false)
108 template<
typename Rhs>
109 inline const internal::solve_retval<HouseholderQR, Rhs>
112 eigen_assert(m_isInitialized &&
"HouseholderQR is not initialized.");
113 return internal::solve_retval<HouseholderQR, Rhs>(*
this, b.derived());
116 HouseholderSequenceType householderQ()
const
118 eigen_assert(m_isInitialized &&
"HouseholderQR is not initialized.");
119 return HouseholderSequenceType(m_qr, m_hCoeffs.conjugate());
127 eigen_assert(m_isInitialized &&
"HouseholderQR is not initialized.");
162 inline Index rows()
const {
return m_qr.rows(); }
163 inline Index cols()
const {
return m_qr.cols(); }
164 const HCoeffsType& hCoeffs()
const {
return m_hCoeffs; }
168 HCoeffsType m_hCoeffs;
169 RowVectorType m_temp;
170 bool m_isInitialized;
173 template<
typename MatrixType>
176 eigen_assert(m_isInitialized &&
"HouseholderQR is not initialized.");
177 eigen_assert(m_qr.rows() == m_qr.cols() &&
"You can't take the determinant of a non-square matrix!");
181 template<
typename MatrixType>
184 eigen_assert(m_isInitialized &&
"HouseholderQR is not initialized.");
185 eigen_assert(m_qr.rows() == m_qr.cols() &&
"You can't take the determinant of a non-square matrix!");
186 return m_qr.diagonal().cwiseAbs().array().log().sum();
192 template<
typename MatrixQR,
typename HCoeffs>
193 void householder_qr_inplace_unblocked(MatrixQR& mat, HCoeffs& hCoeffs,
typename MatrixQR::Scalar* tempData = 0)
195 typedef typename MatrixQR::Index Index;
196 typedef typename MatrixQR::Scalar Scalar;
197 typedef typename MatrixQR::RealScalar RealScalar;
198 Index rows = mat.rows();
199 Index cols = mat.cols();
200 Index size = (std::min)(rows,cols);
202 eigen_assert(hCoeffs.size() == size);
209 tempData = tempVector.data();
212 for(Index k = 0; k < size; ++k)
214 Index remainingRows = rows - k;
215 Index remainingCols = cols - k - 1;
218 mat.col(k).tail(remainingRows).makeHouseholderInPlace(hCoeffs.coeffRef(k), beta);
219 mat.coeffRef(k,k) = beta;
222 mat.bottomRightCorner(remainingRows, remainingCols)
223 .applyHouseholderOnTheLeft(mat.col(k).tail(remainingRows-1), hCoeffs.coeffRef(k), tempData+k+1);
228 template<
typename MatrixQR,
typename HCoeffs>
229 void householder_qr_inplace_blocked(MatrixQR& mat, HCoeffs& hCoeffs,
230 typename MatrixQR::Index maxBlockSize=32,
231 typename MatrixQR::Scalar* tempData = 0)
233 typedef typename MatrixQR::Index Index;
234 typedef typename MatrixQR::Scalar Scalar;
235 typedef typename MatrixQR::RealScalar RealScalar;
236 typedef Block<MatrixQR,Dynamic,Dynamic> BlockType;
238 Index rows = mat.rows();
239 Index cols = mat.cols();
240 Index size = (std::min)(rows, cols);
242 typedef Matrix<Scalar,Dynamic,1,ColMajor,MatrixQR::MaxColsAtCompileTime,1> TempType;
246 tempVector.resize(cols);
247 tempData = tempVector.data();
250 Index blockSize = (std::min)(maxBlockSize,size);
253 for (k = 0; k < size; k += blockSize)
255 Index bs = (std::min)(size-k,blockSize);
256 Index tcols = cols - k - bs;
257 Index brows = rows-k;
267 BlockType A11_21 = mat.block(k,k,brows,bs);
268 Block<HCoeffs,Dynamic,1> hCoeffsSegment = hCoeffs.segment(k,bs);
270 householder_qr_inplace_unblocked(A11_21, hCoeffsSegment, tempData);
274 BlockType A21_22 = mat.block(k,k+bs,brows,tcols);
275 apply_block_householder_on_the_left(A21_22,A11_21,hCoeffsSegment.adjoint());
280 template<
typename _MatrixType,
typename Rhs>
281 struct solve_retval<HouseholderQR<_MatrixType>, Rhs>
282 : solve_retval_base<HouseholderQR<_MatrixType>, Rhs>
284 EIGEN_MAKE_SOLVE_HELPERS(HouseholderQR<_MatrixType>,Rhs)
286 template<typename Dest>
void evalTo(Dest& dst)
const
288 const Index rows = dec().rows(), cols = dec().cols();
289 const Index rank = (std::min)(rows, cols);
290 eigen_assert(rhs().rows() == rows);
292 typename Rhs::PlainObject c(rhs());
296 dec().matrixQR().leftCols(rank),
297 dec().hCoeffs().head(rank)).transpose()
301 .topLeftCorner(rank, rank)
302 .template triangularView<Upper>()
303 .solveInPlace(c.topRows(rank));
305 dst.topRows(rank) = c.topRows(rank);
306 dst.bottomRows(cols-rank).setZero();
312 template<
typename MatrixType>
313 HouseholderQR<MatrixType>& HouseholderQR<MatrixType>::compute(
const MatrixType& matrix)
315 Index rows = matrix.rows();
316 Index cols = matrix.cols();
317 Index size = (std::min)(rows,cols);
320 m_hCoeffs.resize(size);
324 internal::householder_qr_inplace_blocked(m_qr, m_hCoeffs, 48, m_temp.data());
326 m_isInitialized =
true;
334 template<
typename Derived>
335 const HouseholderQR<typename MatrixBase<Derived>::PlainObject>