Fixed a couple of issues introduced in previous commits.

Added a test for Triangular.
This commit is contained in:
Gael Guennebaud 2008-04-26 20:28:27 +00:00
parent b4c974d059
commit 1ec2d21ca5
7 changed files with 53 additions and 29 deletions

View File

@ -104,7 +104,7 @@ bool Vectorize = (Derived::Flags & OtherDerived::Flags & VectorizableBit)
||((Derived::Flags&RowMajorBit) ||((Derived::Flags&RowMajorBit)
? Derived::ColsAtCompileTime!=Dynamic && (Derived::ColsAtCompileTime%ei_packet_traits<typename Derived::Scalar>::size==0) ? Derived::ColsAtCompileTime!=Dynamic && (Derived::ColsAtCompileTime%ei_packet_traits<typename Derived::Scalar>::size==0)
: Derived::RowsAtCompileTime!=Dynamic && (Derived::RowsAtCompileTime%ei_packet_traits<typename Derived::Scalar>::size==0)) ), : Derived::RowsAtCompileTime!=Dynamic && (Derived::RowsAtCompileTime%ei_packet_traits<typename Derived::Scalar>::size==0)) ),
bool TriangularAssign = false> bool TriangularAssign = Derived::Flags & (NullLowerBit | NullUpperBit)>
struct ei_assignment_impl; struct ei_assignment_impl;
template<typename Derived> template<typename Derived>
@ -113,9 +113,7 @@ Derived& MatrixBase<Derived>
::lazyAssign(const MatrixBase<OtherDerived>& other) ::lazyAssign(const MatrixBase<OtherDerived>& other)
{ {
// std::cout << "lazyAssign = " << Derived::Flags << " " << OtherDerived::Flags << "\n"; // std::cout << "lazyAssign = " << Derived::Flags << " " << OtherDerived::Flags << "\n";
ei_assignment_impl<Derived, OtherDerived, ei_assignment_impl<Derived, OtherDerived>::execute(derived(),other.derived());
Derived::Flags & (NullLowerBit | NullUpperBit)>
::execute(derived(),other.derived());
return derived(); return derived();
} }
@ -145,7 +143,7 @@ Derived& MatrixBase<Derived>
} }
template <typename Derived, typename OtherDerived> template <typename Derived, typename OtherDerived>
struct ei_assignment_impl<Derived, OtherDerived, false> struct ei_assignment_impl<Derived, OtherDerived, false, false>
{ {
static void execute(Derived & dst, const OtherDerived & src) static void execute(Derived & dst, const OtherDerived & src)
{ {
@ -179,7 +177,7 @@ struct ei_assignment_impl<Derived, OtherDerived, false>
}; };
template <typename Derived, typename OtherDerived> template <typename Derived, typename OtherDerived>
struct ei_assignment_impl<Derived, OtherDerived, true> struct ei_assignment_impl<Derived, OtherDerived, true, false>
{ {
static void execute(Derived & dst, const OtherDerived & src) static void execute(Derived & dst, const OtherDerived & src)
{ {

View File

@ -52,7 +52,7 @@ struct ei_traits<CwiseUnaryOp<UnaryOp, MatrixType> >
MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime, MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
Flags = (MatrixType::Flags & ( Flags = (MatrixType::Flags & (
DefaultLostFlagMask | Like1DArrayBit DefaultLostFlagMask | Like1DArrayBit
| ei_functor_traits<UnaryOp>::IsVectorizable ? VectorizableBit : 0)), | (ei_functor_traits<UnaryOp>::IsVectorizable ? VectorizableBit : 0))),
CoeffReadCost = MatrixType::CoeffReadCost + ei_functor_traits<UnaryOp>::Cost CoeffReadCost = MatrixType::CoeffReadCost + ei_functor_traits<UnaryOp>::Cost
}; };
}; };

View File

@ -61,10 +61,10 @@ struct ei_traits<Triangular<Mode, MatrixType> >
{ {
typedef typename MatrixType::Scalar Scalar; typedef typename MatrixType::Scalar Scalar;
enum { enum {
RowsAtCompileTime = MatrixType::SizeAtCompileTime, RowsAtCompileTime = MatrixType::RowsAtCompileTime,
ColsAtCompileTime = MatrixType::SizeAtCompileTime, ColsAtCompileTime = MatrixType::ColsAtCompileTime,
MaxRowsAtCompileTime = MatrixType::MaxSizeAtCompileTime, MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
MaxColsAtCompileTime = MatrixType::MaxSizeAtCompileTime, MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
Flags = MatrixType::Flags & (~(VectorizableBit | Like1DArrayBit)) | Mode, Flags = MatrixType::Flags & (~(VectorizableBit | Like1DArrayBit)) | Mode,
CoeffReadCost = MatrixType::CoeffReadCost CoeffReadCost = MatrixType::CoeffReadCost
}; };
@ -84,16 +84,18 @@ template<int Mode, typename MatrixType> class Triangular
assert(matrix.rows()==matrix.cols()); assert(matrix.rows()==matrix.cols());
} }
EIGEN_INHERIT_ASSIGNMENT_OPERATORS(Triangular)
/** Overloaded to keep a Triangular expression */ /** Overloaded to keep a Triangular expression */
Triangular<(Upper | Lower) xor Mode, Transpose<MatrixType> > transpose() Triangular<(Upper | Lower) xor Mode, Temporary<Transpose<MatrixType> > > transpose()
{ {
return Triangular<(Upper | Lower) xor Mode, Transpose<MatrixType> >((m_matrix.transpose())); return Triangular<(Upper | Lower) xor Mode, Temporary<Transpose<MatrixType> > >((m_matrix.transpose().temporary()));
} }
/** Overloaded to keep a Triangular expression */ /** Overloaded to keep a Triangular expression */
const Triangular<(Upper | Lower) xor Mode, Transpose<MatrixType> > transpose() const const Triangular<(Upper | Lower) xor Mode, Temporary<Transpose<MatrixType> > > transpose() const
{ {
return Triangular<(Upper | Lower) xor Mode, Transpose<MatrixType> >((m_matrix.transpose())); return Triangular<(Upper | Lower) xor Mode, Temporary<Transpose<MatrixType> > >((m_matrix.transpose().temporary()));
} }
/** \returns the product of the inverse of *this with \a other. /** \returns the product of the inverse of *this with \a other.
@ -121,7 +123,7 @@ template<int Mode, typename MatrixType> class Triangular
res(0,c) = other(0,c)/_coeff(0, 0); res(0,c) = other(0,c)/_coeff(0, 0);
for (int i=1 ; i<_rows() ; ++i) for (int i=1 ; i<_rows() ; ++i)
{ {
Scalar tmp = other(i,c) - ((this->row(i).start(i)).transpose() * res.col(c).start(i))(0,0); Scalar tmp = other(i,c) - ((this->row(i).start(i)) * res.col(c).start(i))(0,0);
if (Flags & UnitDiagBit) if (Flags & UnitDiagBit)
res(i,c) = tmp; res(i,c) = tmp;
else else
@ -137,7 +139,7 @@ template<int Mode, typename MatrixType> class Triangular
res(_cols()-1,c) = other(_cols()-1, c)/_coeff(_rows()-1, _cols()-1); res(_cols()-1,c) = other(_cols()-1, c)/_coeff(_rows()-1, _cols()-1);
for (int i=_rows()-2 ; i>=0 ; --i) for (int i=_rows()-2 ; i>=0 ; --i)
{ {
Scalar tmp = other(i,c) - ((this->row(i).end(_cols()-i-1)).transpose() * res.col(c).end(_cols()-i-1))(0,0); Scalar tmp = other(i,c) - ((this->row(i).end(_cols()-i-1)) * res.col(c).end(_cols()-i-1))(0,0);
if (Flags & UnitDiagBit) if (Flags & UnitDiagBit)
res(i,c) = tmp; res(i,c) = tmp;
else else
@ -155,8 +157,8 @@ template<int Mode, typename MatrixType> class Triangular
Scalar& _coeffRef(int row, int col) Scalar& _coeffRef(int row, int col)
{ {
assert( ((! Flags & Lower) && row<=col) || (Flags & Lower && col<=row)); ei_assert( ((! (Flags & Lower)) && row<=col) || (Flags & Lower && col<=row));
return m_matrix.coeffRef(row, col); return m_matrix.const_cast_derived().coeffRef(row, col);
} }
Scalar _coeff(int row, int col) const Scalar _coeff(int row, int col) const

View File

@ -18,6 +18,17 @@ USING_PART_OF_NAMESPACE_EIGEN
int main(int argc, char *argv[]) int main(int argc, char *argv[])
{ {
Matrix4i m1, m2, m3;
m1.setRandom();
m2.setConstant(2);
int s1 = 2;
m3 = m1;
std::cout << m1 << "\n\n";
std::cout << m2 << "\n\n";
m3 = m1.cwiseProduct(m2);
std::cout << m3 << "\n==\n" << m1*s1 << "\n\n";
// v(1,2,3,4);
// std::cout << v * 2 << "\n";
Matrix<SCALAR,MATSIZE,MATSIZE> I = Matrix<SCALAR,MATSIZE,MATSIZE>::ones(); Matrix<SCALAR,MATSIZE,MATSIZE> I = Matrix<SCALAR,MATSIZE,MATSIZE>::ones();
Matrix<SCALAR,MATSIZE,MATSIZE> m; Matrix<SCALAR,MATSIZE,MATSIZE> m;
for(int i = 0; i < MATSIZE; i++) for(int i = 0; i < MATSIZE; i++)

View File

@ -7,6 +7,7 @@ FIND_PACKAGE(Qt4 REQUIRED)
INCLUDE_DIRECTORIES( ${QT_INCLUDE_DIR} ) INCLUDE_DIRECTORIES( ${QT_INCLUDE_DIR} )
SET(test_SRCS SET(test_SRCS
triangular.cpp
main.cpp main.cpp
basicstuff.cpp basicstuff.cpp
linearstructure.cpp linearstructure.cpp

View File

@ -40,21 +40,32 @@ template<typename MatrixType> void nullDeterminant(const MatrixType& m)
typedef Matrix<Scalar, MatrixType::ColsAtCompileTime, MatrixType::ColsAtCompileTime> SquareMatrixType; typedef Matrix<Scalar, MatrixType::ColsAtCompileTime, MatrixType::ColsAtCompileTime> SquareMatrixType;
typedef Matrix<Scalar, MatrixType::ColsAtCompileTime, 1> VectorType; typedef Matrix<Scalar, MatrixType::ColsAtCompileTime, 1> VectorType;
MatrixType d(rows, cols); MatrixType dinv(rows, cols), dnotinv(rows, cols);
// build a ill-conditionned matrix with a nul determinant dinv.col(0).setOnes();
d.col(0).setOnes(); dinv.block(0,1, rows, cols-2).setRandom();
d.block(0,1, rows, cols-2).setRandom();
d.col(cols-1).setOnes(); dnotinv.col(0).setOnes();
dnotinv.block(0,1, rows, cols-2).setRandom();
dnotinv.col(cols-1).setOnes();
for (int i=0 ; i<rows ; ++i) for (int i=0 ; i<rows ; ++i)
d.row(i).block(0,1,1,cols-2) = d.row(i).block(0,1,1,cols-2).normalized(); {
dnotinv.row(i).block(0,1,1,cols-2) = ei_random<Scalar>(99.999999,100.00000001)*dnotinv.row(i).block(0,1,1,cols-2).normalized();
dnotinv(i,cols-1) = dnotinv.row(i).block(0,1,1,cols-2).norm2();
dinv(i,cols-1) = dinv.row(i).block(0,1,1,cols-2).norm2();
}
SquareMatrixType covarianceMatrix = d.transpose() * d; SquareMatrixType invertibleCovarianceMatrix = dinv.transpose() * dinv;
SquareMatrixType notInvertibleCovarianceMatrix = dnotinv.transpose() * dnotinv;
// std::cout << covarianceMatrix << "\n" << covarianceMatrix.determinant() << "\n"; std::cout << notInvertibleCovarianceMatrix << "\n" << notInvertibleCovarianceMatrix.determinant() << "\n";
VERIFY_IS_APPROX(covarianceMatrix.determinant(), Scalar(0)); VERIFY_IS_APPROX(notInvertibleCovarianceMatrix.determinant(), Scalar(0));
VERIFY(invertibleCovarianceMatrix.inverse().exists());
VERIFY(!notInvertibleCovarianceMatrix.inverse().exists());
} }
void EigenTest::testDeterminant() void EigenTest::testDeterminant()

View File

@ -32,7 +32,7 @@
#include <ctime> #include <ctime>
#include <iostream> #include <iostream>
#define DEFAULT_REPEAT 50 #define DEFAULT_REPEAT 10
#ifndef EIGEN_NO_ASSERTION_CHECKING #ifndef EIGEN_NO_ASSERTION_CHECKING
@ -214,6 +214,7 @@ class EigenTest : public QObject
void testMap(); void testMap();
void testCwiseops(); void testCwiseops();
void testDeterminant(); void testDeterminant();
void testTriangular();
protected: protected:
int m_repeat; int m_repeat;
}; };