* update test to expose bug #57

* update createRandomMatrixOfRank to support fixed size
This commit is contained in:
Benoit Jacob 2009-09-28 09:40:18 -04:00
parent de942a44c2
commit 67bf7c90c5
3 changed files with 46 additions and 6 deletions

View File

@ -351,7 +351,6 @@ bool ColPivotingHouseholderQR<MatrixType>::solve(
}
const int rows = m_qr.rows();
const int cols = b.cols();
ei_assert(b.rows() == rows);
typename OtherDerived::PlainMatrixType c(b);

View File

@ -248,18 +248,22 @@ template<typename MatrixType>
void createRandomMatrixOfRank(int desired_rank, int rows, int cols, MatrixType& m)
{
typedef typename ei_traits<MatrixType>::Scalar Scalar;
typedef Matrix<Scalar, Dynamic, 1> VectorType;
enum { Rows = MatrixType::RowsAtCompileTime, Cols = MatrixType::ColsAtCompileTime };
MatrixType a = MatrixType::Random(rows,rows);
typedef Matrix<Scalar, Dynamic, 1> VectorType;
typedef Matrix<Scalar, Rows, Rows> MatrixAType;
typedef Matrix<Scalar, Cols, Cols> MatrixBType;
MatrixAType a = MatrixAType::Random(rows,rows);
MatrixType d = MatrixType::Identity(rows,cols);
MatrixType b = MatrixType::Random(cols,cols);
MatrixBType b = MatrixBType::Random(cols,cols);
// set the diagonal such that only desired_rank non-zero entries reamain
const int diag_size = std::min(d.rows(),d.cols());
d.diagonal().segment(desired_rank, diag_size-desired_rank) = VectorType::Zero(diag_size-desired_rank);
HouseholderQR<MatrixType> qra(a);
HouseholderQR<MatrixType> qrb(b);
HouseholderQR<MatrixAType> qra(a);
HouseholderQR<MatrixBType> qrb(b);
m = qra.matrixQ() * d * qrb.matrixQ();
}

View File

@ -63,6 +63,40 @@ template<typename MatrixType> void qr()
VERIFY(!qr.solve(m3, &m2));
}
template<typename MatrixType, int Cols2> void qr_fixedsize()
{
enum { Rows = MatrixType::RowsAtCompileTime, Cols = MatrixType::ColsAtCompileTime };
typedef typename MatrixType::Scalar Scalar;
int rank = ei_random<int>(1, std::min(int(Rows), int(Cols))-1);
Matrix<Scalar,Rows,Cols> m1;
createRandomMatrixOfRank(rank,Rows,Cols,m1);
ColPivotingHouseholderQR<Matrix<Scalar,Rows,Cols> > qr(m1);
VERIFY_IS_APPROX(rank, qr.rank());
VERIFY(Cols - qr.rank() == qr.dimensionOfKernel());
VERIFY(!qr.isInjective());
VERIFY(!qr.isInvertible());
VERIFY(!qr.isSurjective());
Matrix<Scalar,Rows,Cols> r = qr.matrixQR();
// FIXME need better way to construct trapezoid
for(int i = 0; i < Rows; i++) for(int j = 0; j < Cols; j++) if(i>j) r(i,j) = Scalar(0);
Matrix<Scalar,Rows,Cols> b = qr.matrixQ() * r;
Matrix<Scalar,Rows,Cols> c = MatrixType::Zero(Rows,Cols);
for(int i = 0; i < Cols; ++i) c.col(qr.colsPermutation().coeff(i)) = b.col(i);
VERIFY_IS_APPROX(m1, c);
Matrix<Scalar,Cols,Cols2> m2 = Matrix<Scalar,Cols,Cols2>::Random(Cols,Cols2);
Matrix<Scalar,Rows,Cols2> m3 = m1*m2;
m2 = Matrix<Scalar,Cols,Cols2>::Random(Cols,Cols2);
VERIFY(qr.solve(m3, &m2));
VERIFY_IS_APPROX(m3, m1*m2);
m3 = Matrix<Scalar,Rows,Cols2>::Random(Rows,Cols2);
VERIFY(!qr.solve(m3, &m2));
}
template<typename MatrixType> void qr_invertible()
{
typedef typename NumTraits<typename MatrixType::Scalar>::Real RealScalar;
@ -121,6 +155,9 @@ void test_qr_colpivoting()
CALL_SUBTEST( qr<MatrixXd>() );
CALL_SUBTEST( qr<MatrixXcd>() );
}
CALL_SUBTEST(( qr_fixedsize<Matrix<float,3,5>, 4 >() ));
CALL_SUBTEST(( qr_fixedsize<Matrix<double,6,2>, 3 >() ));
for(int i = 0; i < g_repeat; i++) {
CALL_SUBTEST( qr_invertible<MatrixXf>() );