mirror of
https://gitlab.com/libeigen/eigen.git
synced 2025-09-22 22:33:15 +08:00
Fix HouseholderSequence::conjugate() and ::adjoint() and add respective unit tests.
This commit is contained in:
parent
9f11f80db1
commit
55365566b2
@ -82,7 +82,7 @@ template<typename _MatrixType> class HessenbergDecomposition
|
|||||||
typedef Matrix<Scalar, SizeMinusOne, 1, Options & ~RowMajor, MaxSizeMinusOne, 1> CoeffVectorType;
|
typedef Matrix<Scalar, SizeMinusOne, 1, Options & ~RowMajor, MaxSizeMinusOne, 1> CoeffVectorType;
|
||||||
|
|
||||||
/** \brief Return type of matrixQ() */
|
/** \brief Return type of matrixQ() */
|
||||||
typedef typename HouseholderSequence<MatrixType,CoeffVectorType>::ConjugateReturnType HouseholderSequenceType;
|
typedef HouseholderSequence<MatrixType,typename internal::remove_all<typename CoeffVectorType::ConjugateReturnType>::type> HouseholderSequenceType;
|
||||||
|
|
||||||
typedef internal::HessenbergDecompositionMatrixHReturnType<MatrixType> MatrixHReturnType;
|
typedef internal::HessenbergDecompositionMatrixHReturnType<MatrixType> MatrixHReturnType;
|
||||||
|
|
||||||
|
@ -96,7 +96,7 @@ template<typename _MatrixType> class Tridiagonalization
|
|||||||
>::type SubDiagonalReturnType;
|
>::type SubDiagonalReturnType;
|
||||||
|
|
||||||
/** \brief Return type of matrixQ() */
|
/** \brief Return type of matrixQ() */
|
||||||
typedef typename HouseholderSequence<MatrixType,CoeffVectorType>::ConjugateReturnType HouseholderSequenceType;
|
typedef HouseholderSequence<MatrixType,typename internal::remove_all<typename CoeffVectorType::ConjugateReturnType>::type> HouseholderSequenceType;
|
||||||
|
|
||||||
/** \brief Default constructor.
|
/** \brief Default constructor.
|
||||||
*
|
*
|
||||||
|
@ -125,7 +125,9 @@ template<typename VectorsType, typename CoeffsType, int Side> class HouseholderS
|
|||||||
typedef typename VectorsType::Index Index;
|
typedef typename VectorsType::Index Index;
|
||||||
|
|
||||||
typedef HouseholderSequence<
|
typedef HouseholderSequence<
|
||||||
VectorsType,
|
typename internal::conditional<NumTraits<Scalar>::IsComplex,
|
||||||
|
typename internal::remove_all<typename VectorsType::ConjugateReturnType>::type,
|
||||||
|
VectorsType>::type,
|
||||||
typename internal::conditional<NumTraits<Scalar>::IsComplex,
|
typename internal::conditional<NumTraits<Scalar>::IsComplex,
|
||||||
typename internal::remove_all<typename CoeffsType::ConjugateReturnType>::type,
|
typename internal::remove_all<typename CoeffsType::ConjugateReturnType>::type,
|
||||||
CoeffsType>::type,
|
CoeffsType>::type,
|
||||||
@ -206,7 +208,7 @@ template<typename VectorsType, typename CoeffsType, int Side> class HouseholderS
|
|||||||
/** \brief Complex conjugate of the Householder sequence. */
|
/** \brief Complex conjugate of the Householder sequence. */
|
||||||
ConjugateReturnType conjugate() const
|
ConjugateReturnType conjugate() const
|
||||||
{
|
{
|
||||||
return ConjugateReturnType(m_vectors, m_coeffs.conjugate())
|
return ConjugateReturnType(m_vectors.conjugate(), m_coeffs.conjugate())
|
||||||
.setTrans(m_trans)
|
.setTrans(m_trans)
|
||||||
.setLength(m_length)
|
.setLength(m_length)
|
||||||
.setShift(m_shift);
|
.setShift(m_shift);
|
||||||
|
@ -55,7 +55,7 @@ template<typename _MatrixType> class ColPivHouseholderQR
|
|||||||
typedef typename internal::plain_row_type<MatrixType, Index>::type IntRowVectorType;
|
typedef typename internal::plain_row_type<MatrixType, Index>::type IntRowVectorType;
|
||||||
typedef typename internal::plain_row_type<MatrixType>::type RowVectorType;
|
typedef typename internal::plain_row_type<MatrixType>::type RowVectorType;
|
||||||
typedef typename internal::plain_row_type<MatrixType, RealScalar>::type RealRowVectorType;
|
typedef typename internal::plain_row_type<MatrixType, RealScalar>::type RealRowVectorType;
|
||||||
typedef typename HouseholderSequence<MatrixType,HCoeffsType>::ConjugateReturnType HouseholderSequenceType;
|
typedef HouseholderSequence<MatrixType,typename internal::remove_all<typename HCoeffsType::ConjugateReturnType>::type> HouseholderSequenceType;
|
||||||
|
|
||||||
private:
|
private:
|
||||||
|
|
||||||
|
@ -57,7 +57,7 @@ template<typename _MatrixType> class HouseholderQR
|
|||||||
typedef Matrix<Scalar, RowsAtCompileTime, RowsAtCompileTime, (MatrixType::Flags&RowMajorBit) ? RowMajor : ColMajor, MaxRowsAtCompileTime, MaxRowsAtCompileTime> MatrixQType;
|
typedef Matrix<Scalar, RowsAtCompileTime, RowsAtCompileTime, (MatrixType::Flags&RowMajorBit) ? RowMajor : ColMajor, MaxRowsAtCompileTime, MaxRowsAtCompileTime> MatrixQType;
|
||||||
typedef typename internal::plain_diag_type<MatrixType>::type HCoeffsType;
|
typedef typename internal::plain_diag_type<MatrixType>::type HCoeffsType;
|
||||||
typedef typename internal::plain_row_type<MatrixType>::type RowVectorType;
|
typedef typename internal::plain_row_type<MatrixType>::type RowVectorType;
|
||||||
typedef typename HouseholderSequence<MatrixType,HCoeffsType>::ConjugateReturnType HouseholderSequenceType;
|
typedef HouseholderSequence<MatrixType,typename internal::remove_all<typename HCoeffsType::ConjugateReturnType>::type> HouseholderSequenceType;
|
||||||
|
|
||||||
/**
|
/**
|
||||||
* \brief Default Constructor.
|
* \brief Default Constructor.
|
||||||
|
@ -39,7 +39,7 @@ template<typename _MatrixType> class UpperBidiagonalization
|
|||||||
CwiseUnaryOp<internal::scalar_conjugate_op<Scalar>, const Diagonal<const MatrixType,0> >
|
CwiseUnaryOp<internal::scalar_conjugate_op<Scalar>, const Diagonal<const MatrixType,0> >
|
||||||
> HouseholderUSequenceType;
|
> HouseholderUSequenceType;
|
||||||
typedef HouseholderSequence<
|
typedef HouseholderSequence<
|
||||||
const MatrixType,
|
const typename internal::remove_all<typename MatrixType::ConjugateReturnType>::type,
|
||||||
Diagonal<const MatrixType,1>,
|
Diagonal<const MatrixType,1>,
|
||||||
OnTheRight
|
OnTheRight
|
||||||
> HouseholderVSequenceType;
|
> HouseholderVSequenceType;
|
||||||
@ -74,7 +74,7 @@ template<typename _MatrixType> class UpperBidiagonalization
|
|||||||
const HouseholderVSequenceType householderV() // const here gives nasty errors and i'm lazy
|
const HouseholderVSequenceType householderV() // const here gives nasty errors and i'm lazy
|
||||||
{
|
{
|
||||||
eigen_assert(m_isInitialized && "UpperBidiagonalization is not initialized.");
|
eigen_assert(m_isInitialized && "UpperBidiagonalization is not initialized.");
|
||||||
return HouseholderVSequenceType(m_householder, m_householder.const_derived().template diagonal<1>())
|
return HouseholderVSequenceType(m_householder.conjugate(), m_householder.const_derived().template diagonal<1>())
|
||||||
.setLength(m_householder.cols()-1)
|
.setLength(m_householder.cols()-1)
|
||||||
.setShift(1);
|
.setShift(1);
|
||||||
}
|
}
|
||||||
|
@ -89,12 +89,29 @@ template<typename MatrixType> void householder(const MatrixType& m)
|
|||||||
hseq.setLength(hc.size()).setShift(shift);
|
hseq.setLength(hc.size()).setShift(shift);
|
||||||
VERIFY(hseq.length() == hc.size());
|
VERIFY(hseq.length() == hc.size());
|
||||||
VERIFY(hseq.shift() == shift);
|
VERIFY(hseq.shift() == shift);
|
||||||
|
|
||||||
MatrixType m5 = m2;
|
MatrixType m5 = m2;
|
||||||
m5.block(shift,0,brows,cols).template triangularView<StrictlyLower>().setZero();
|
m5.block(shift,0,brows,cols).template triangularView<StrictlyLower>().setZero();
|
||||||
VERIFY_IS_APPROX(hseq * m5, m1); // test applying hseq directly
|
VERIFY_IS_APPROX(hseq * m5, m1); // test applying hseq directly
|
||||||
m3 = hseq;
|
m3 = hseq;
|
||||||
VERIFY_IS_APPROX(m3 * m5, m1); // test evaluating hseq to a dense matrix, then applying
|
VERIFY_IS_APPROX(m3 * m5, m1); // test evaluating hseq to a dense matrix, then applying
|
||||||
|
|
||||||
|
SquareMatrixType hseq_mat = hseq;
|
||||||
|
SquareMatrixType hseq_mat_conj = hseq.conjugate();
|
||||||
|
SquareMatrixType hseq_mat_adj = hseq.adjoint();
|
||||||
|
SquareMatrixType hseq_mat_trans = hseq.transpose();
|
||||||
|
SquareMatrixType m6 = SquareMatrixType::Random(rows, rows);
|
||||||
|
VERIFY_IS_APPROX(hseq_mat.adjoint(), hseq_mat_adj);
|
||||||
|
VERIFY_IS_APPROX(hseq_mat.conjugate(), hseq_mat_conj);
|
||||||
|
VERIFY_IS_APPROX(hseq_mat.transpose(), hseq_mat_trans);
|
||||||
|
VERIFY_IS_APPROX(hseq_mat * m6, hseq_mat * m6);
|
||||||
|
VERIFY_IS_APPROX(hseq_mat.adjoint() * m6, hseq_mat_adj * m6);
|
||||||
|
VERIFY_IS_APPROX(hseq_mat.conjugate() * m6, hseq_mat_conj * m6);
|
||||||
|
VERIFY_IS_APPROX(hseq_mat.transpose() * m6, hseq_mat_trans * m6);
|
||||||
|
VERIFY_IS_APPROX(m6 * hseq_mat, m6 * hseq_mat);
|
||||||
|
VERIFY_IS_APPROX(m6 * hseq_mat.adjoint(), m6 * hseq_mat_adj);
|
||||||
|
VERIFY_IS_APPROX(m6 * hseq_mat.conjugate(), m6 * hseq_mat_conj);
|
||||||
|
VERIFY_IS_APPROX(m6 * hseq_mat.transpose(), m6 * hseq_mat_trans);
|
||||||
|
|
||||||
// test householder sequence on the right with a shift
|
// test householder sequence on the right with a shift
|
||||||
|
|
||||||
|
@ -16,6 +16,7 @@ template<typename MatrixType> void upperbidiag(const MatrixType& m)
|
|||||||
const typename MatrixType::Index cols = m.cols();
|
const typename MatrixType::Index cols = m.cols();
|
||||||
|
|
||||||
typedef Matrix<typename MatrixType::RealScalar, MatrixType::RowsAtCompileTime, MatrixType::ColsAtCompileTime> RealMatrixType;
|
typedef Matrix<typename MatrixType::RealScalar, MatrixType::RowsAtCompileTime, MatrixType::ColsAtCompileTime> RealMatrixType;
|
||||||
|
typedef Matrix<typename MatrixType::Scalar, MatrixType::ColsAtCompileTime, MatrixType::RowsAtCompileTime> TransposeMatrixType;
|
||||||
|
|
||||||
MatrixType a = MatrixType::Random(rows,cols);
|
MatrixType a = MatrixType::Random(rows,cols);
|
||||||
internal::UpperBidiagonalization<MatrixType> ubd(a);
|
internal::UpperBidiagonalization<MatrixType> ubd(a);
|
||||||
@ -24,6 +25,8 @@ template<typename MatrixType> void upperbidiag(const MatrixType& m)
|
|||||||
b.block(0,0,cols,cols) = ubd.bidiagonal();
|
b.block(0,0,cols,cols) = ubd.bidiagonal();
|
||||||
MatrixType c = ubd.householderU() * b * ubd.householderV().adjoint();
|
MatrixType c = ubd.householderU() * b * ubd.householderV().adjoint();
|
||||||
VERIFY_IS_APPROX(a,c);
|
VERIFY_IS_APPROX(a,c);
|
||||||
|
TransposeMatrixType d = ubd.householderV() * b.adjoint() * ubd.householderU().adjoint();
|
||||||
|
VERIFY_IS_APPROX(a.adjoint(),d);
|
||||||
}
|
}
|
||||||
|
|
||||||
void test_upperbidiagonalization()
|
void test_upperbidiagonalization()
|
||||||
|
Loading…
x
Reference in New Issue
Block a user