mirror/eigen
2
0
Fork 0
mirror of https://gitlab.com/libeigen/eigen.git synced 2024-12-21 07:19:46 +08:00
Code Issues Packages Projects Releases Wiki Activity

bug #206 - part 2: For HouseholderSequence objects, added non-allocating versions of evalTo() and applyThisOnTheRight/Left that take additional working vector parameters.

Browse Source
This commit is contained in:
Adolfo Rodriguez Tsourouksdissian 2011-10-30 23:55:16 -04:00
parent bca18a13ea
commit 7bf0e8cd82
1 changed files with 36 additions and 13 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

49
Eigen/src/Householder/HouseholderSequence.h
View File

@ -237,13 +237,20 @@ template<typename VectorsType, typename CoeffsType, int Side> class HouseholderS
ConjugateReturnType inverse() const { return adjoint(); } ConjugateReturnType inverse() const { return adjoint(); }
/** \internal */ /** \internal */
template<typename DestType> void evalTo(DestType& dst) const template<typename DestType> inline void evalTo(DestType& dst) const
{ {
Index vecs = m_length;
// FIXME find a way to pass this temporary if the user wants to
Matrix<Scalar, DestType::RowsAtCompileTime, 1, Matrix<Scalar, DestType::RowsAtCompileTime, 1,
AutoAlign|ColMajor, DestType::MaxRowsAtCompileTime, 1> temp(rows()); AutoAlign|ColMajor, DestType::MaxRowsAtCompileTime, 1> workspace(rows());
if( internal::is_same<typename internal::remove_all<VectorsType>::type,DestType>::value evalTo(dst, workspace);
}
/** \internal */
template<typename Dest, typename Workspace>
void evalTo(Dest& dst, Workspace& workspace) const
{
workspace.resize(rows());
Index vecs = m_length;
if( internal::is_same<typename internal::remove_all<VectorsType>::type,Dest>::value
&& internal::extract_data(dst) == internal::extract_data(m_vectors)) && internal::extract_data(dst) == internal::extract_data(m_vectors))
{ {
// in-place // in-place
@ -254,10 +261,10 @@ template<typename VectorsType, typename CoeffsType, int Side> class HouseholderS
Index cornerSize = rows() - k - m_shift; Index cornerSize = rows() - k - m_shift;
if(m_trans) if(m_trans)
dst.bottomRightCorner(cornerSize, cornerSize) dst.bottomRightCorner(cornerSize, cornerSize)
.applyHouseholderOnTheRight(essentialVector(k), m_coeffs.coeff(k), &temp.coeffRef(0)); .applyHouseholderOnTheRight(essentialVector(k), m_coeffs.coeff(k), workspace.data());
else else
dst.bottomRightCorner(cornerSize, cornerSize) dst.bottomRightCorner(cornerSize, cornerSize)
.applyHouseholderOnTheLeft(essentialVector(k), m_coeffs.coeff(k), &temp.coeffRef(0)); .applyHouseholderOnTheLeft(essentialVector(k), m_coeffs.coeff(k), workspace.data());
// clear the off diagonal vector // clear the off diagonal vector
dst.col(k).tail(rows()-k-1).setZero(); dst.col(k).tail(rows()-k-1).setZero();
@ -274,10 +281,10 @@ template<typename VectorsType, typename CoeffsType, int Side> class HouseholderS
Index cornerSize = rows() - k - m_shift; Index cornerSize = rows() - k - m_shift;
if(m_trans) if(m_trans)
dst.bottomRightCorner(cornerSize, cornerSize) dst.bottomRightCorner(cornerSize, cornerSize)
.applyHouseholderOnTheRight(essentialVector(k), m_coeffs.coeff(k), &temp.coeffRef(0)); .applyHouseholderOnTheRight(essentialVector(k), m_coeffs.coeff(k), &workspace.coeffRef(0));
else else
dst.bottomRightCorner(cornerSize, cornerSize) dst.bottomRightCorner(cornerSize, cornerSize)
.applyHouseholderOnTheLeft(essentialVector(k), m_coeffs.coeff(k), &temp.coeffRef(0)); .applyHouseholderOnTheLeft(essentialVector(k), m_coeffs.coeff(k), &workspace.coeffRef(0));
} }
} }
} }
@ -285,24 +292,40 @@ template<typename VectorsType, typename CoeffsType, int Side> class HouseholderS
/** \internal */ /** \internal */
template<typename Dest> inline void applyThisOnTheRight(Dest& dst) const template<typename Dest> inline void applyThisOnTheRight(Dest& dst) const
{ {
Matrix<Scalar,1,Dest::RowsAtCompileTime> temp(dst.rows()); Matrix<Scalar,1,Dest::RowsAtCompileTime,RowMajor,1,Dest::MaxRowsAtCompileTime> workspace(dst.rows());
applyThisOnTheRight(dst, workspace);
}
/** \internal */
template<typename Dest, typename Workspace>
inline void applyThisOnTheRight(Dest& dst, Workspace& workspace) const
{
workspace.resize(dst.rows());
for(Index k = 0; k < m_length; ++k) for(Index k = 0; k < m_length; ++k)
{ {
Index actual_k = m_trans ? m_length-k-1 : k; Index actual_k = m_trans ? m_length-k-1 : k;
dst.rightCols(rows()-m_shift-actual_k) dst.rightCols(rows()-m_shift-actual_k)
.applyHouseholderOnTheRight(essentialVector(actual_k), m_coeffs.coeff(actual_k), &temp.coeffRef(0)); .applyHouseholderOnTheRight(essentialVector(actual_k), m_coeffs.coeff(actual_k), workspace.data());
} }
} }
/** \internal */ /** \internal */
template<typename Dest> inline void applyThisOnTheLeft(Dest& dst) const template<typename Dest> inline void applyThisOnTheLeft(Dest& dst) const
{ {
Matrix<Scalar,1,Dest::ColsAtCompileTime> temp(dst.cols()); Matrix<Scalar,1,Dest::ColsAtCompileTime,RowMajor,1,Dest::MaxColsAtCompileTime> workspace(dst.cols());
applyThisOnTheLeft(dst, workspace);
}
/** \internal */
template<typename Dest, typename Workspace>
inline void applyThisOnTheLeft(Dest& dst, Workspace& workspace) const
{
workspace.resize(dst.cols());
for(Index k = 0; k < m_length; ++k) for(Index k = 0; k < m_length; ++k)
{ {
Index actual_k = m_trans ? k : m_length-k-1; Index actual_k = m_trans ? k : m_length-k-1;
dst.bottomRows(rows()-m_shift-actual_k) dst.bottomRows(rows()-m_shift-actual_k)
.applyHouseholderOnTheLeft(essentialVector(actual_k), m_coeffs.coeff(actual_k), &temp.coeffRef(0)); .applyHouseholderOnTheLeft(essentialVector(actual_k), m_coeffs.coeff(actual_k), workspace.data());
} }
} }