diff --git a/Eigen/Core b/Eigen/Core index 5260ce7e8..56005f836 100644 --- a/Eigen/Core +++ b/Eigen/Core @@ -93,7 +93,14 @@ #endif // include files - + #if (defined __GNUC__) && (defined __MINGW32__) + #include + //including intrin.h works around a MINGW bug http://sourceforge.net/tracker/?func=detail&atid=102435&aid=2962480&group_id=2435 + //in essence, intrin.h is included by windows.h and also declares intrinsics (just as emmintrin.h etc. below do). However, + //intrin.h uses an extern "C" declaration, and g++ thus complains of duplicate declarations with conflicting linkage. The linkage for intrinsics + //doesn't matter, but at that stage the compiler doesn't know; so, to avoid compile errors when windows.h is included after Eigen/Core, + //include intrin here. + #endif #include #include #ifdef EIGEN_VECTORIZE_SSE3 diff --git a/Eigen/src/CMakeLists.txt b/Eigen/src/CMakeLists.txt index 761c83334..52696a803 100644 --- a/Eigen/src/CMakeLists.txt +++ b/Eigen/src/CMakeLists.txt @@ -1,12 +1,6 @@ -ADD_SUBDIRECTORY(Core) -ADD_SUBDIRECTORY(LU) -ADD_SUBDIRECTORY(QR) -ADD_SUBDIRECTORY(SVD) -ADD_SUBDIRECTORY(Cholesky) -ADD_SUBDIRECTORY(Geometry) -ADD_SUBDIRECTORY(Sparse) -ADD_SUBDIRECTORY(Jacobi) -ADD_SUBDIRECTORY(Householder) -ADD_SUBDIRECTORY(Eigenvalues) -ADD_SUBDIRECTORY(misc) -ADD_SUBDIRECTORY(plugins) +file(GLOB Eigen_src_subdirectories "*") +foreach(f ${Eigen_src_subdirectories}) + if(NOT f MATCHES ".txt") + add_subdirectory(${f}) + endif() +endforeach() diff --git a/Eigen/src/Core/CwiseBinaryOp.h b/Eigen/src/Core/CwiseBinaryOp.h index 171140c27..5def0db2a 100644 --- a/Eigen/src/Core/CwiseBinaryOp.h +++ b/Eigen/src/Core/CwiseBinaryOp.h @@ -216,7 +216,7 @@ EIGEN_STRONG_INLINE Derived & MatrixBase::operator-=(const MatrixBase &other) { SelfCwiseBinaryOp, Derived, OtherDerived> tmp(derived()); - tmp = other; + tmp = other.derived(); return derived(); } diff --git a/Eigen/src/Core/NoAlias.h b/Eigen/src/Core/NoAlias.h index d34f83b7b..53ad3bfee 100644 --- a/Eigen/src/Core/NoAlias.h +++ b/Eigen/src/Core/NoAlias.h @@ -43,6 +43,7 @@ template class StorageBase> class NoAlias { + typedef typename ExpressionType::Scalar Scalar; public: NoAlias(ExpressionType& expression) : m_expression(expression) {} @@ -50,17 +51,31 @@ class NoAlias * \sa MatrixBase::lazyAssign() */ template EIGEN_STRONG_INLINE ExpressionType& operator=(const StorageBase& other) - { return m_expression.lazyAssign(other.derived()); } + { return ei_assign_selector::run(m_expression,other.derived()); } /** \sa MatrixBase::operator+= */ template EIGEN_STRONG_INLINE ExpressionType& operator+=(const StorageBase& other) - { return m_expression.lazyAssign(m_expression + other.derived()); } + { + typedef SelfCwiseBinaryOp, ExpressionType, OtherDerived> SelfAdder; + SelfAdder tmp(m_expression); + typedef typename ei_nested::type OtherDerivedNested; + typedef typename ei_cleantype::type _OtherDerivedNested; + ei_assign_selector::run(tmp,OtherDerivedNested(other.derived())); + return m_expression; + } /** \sa MatrixBase::operator-= */ template EIGEN_STRONG_INLINE ExpressionType& operator-=(const StorageBase& other) - { return m_expression.lazyAssign(m_expression - other.derived()); } + { + typedef SelfCwiseBinaryOp, ExpressionType, OtherDerived> SelfAdder; + SelfAdder tmp(m_expression); + typedef typename ei_nested::type OtherDerivedNested; + typedef typename ei_cleantype::type _OtherDerivedNested; + ei_assign_selector::run(tmp,OtherDerivedNested(other.derived())); + return m_expression; + } #ifndef EIGEN_PARSED_BY_DOXYGEN template diff --git a/Eigen/src/Core/SelfCwiseBinaryOp.h b/Eigen/src/Core/SelfCwiseBinaryOp.h index 5100f6b25..f77589747 100644 --- a/Eigen/src/Core/SelfCwiseBinaryOp.h +++ b/Eigen/src/Core/SelfCwiseBinaryOp.h @@ -62,8 +62,6 @@ template class SelfCwiseBinaryOp typedef typename ei_packet_traits::type Packet; - using Base::operator=; - inline SelfCwiseBinaryOp(Lhs& xpr, const BinaryOp& func = BinaryOp()) : m_matrix(xpr), m_functor(func) {} inline Index rows() const { return m_matrix.rows(); } @@ -142,6 +140,15 @@ template class SelfCwiseBinaryOp #endif return *this; } + + // overloaded to honor evaluation of special matrices + // maybe another solution would be to not use SelfCwiseBinaryOp + // at first... + SelfCwiseBinaryOp& operator=(const Rhs& _rhs) + { + typename ei_nested::type rhs(_rhs); + return Base::operator=(rhs); + } protected: Lhs& m_matrix; diff --git a/Eigen/src/Eigen2Support/CMakeLists.txt b/Eigen/src/Eigen2Support/CMakeLists.txt new file mode 100644 index 000000000..2d635042e --- /dev/null +++ b/Eigen/src/Eigen2Support/CMakeLists.txt @@ -0,0 +1,6 @@ +FILE(GLOB Eigen_Eigen2Support_SRCS "*.h") + +INSTALL(FILES + ${Eigen_Eigen2Support_SRCS} + DESTINATION ${INCLUDE_INSTALL_DIR}/Eigen/src/Eigen2Support COMPONENT Devel + ) diff --git a/Eigen/src/Geometry/OrthoMethods.h b/Eigen/src/Geometry/OrthoMethods.h index cb8d3458a..d03d85beb 100644 --- a/Eigen/src/Geometry/OrthoMethods.h +++ b/Eigen/src/Geometry/OrthoMethods.h @@ -54,7 +54,7 @@ MatrixBase::cross(const MatrixBase& other) const template< int Arch,typename VectorLhs,typename VectorRhs, typename Scalar = typename VectorLhs::Scalar, - int Vectorizable = (VectorLhs::Flags&VectorRhs::Flags)&PacketAccessBit> + bool Vectorizable = (VectorLhs::Flags&VectorRhs::Flags)&PacketAccessBit> struct ei_cross3_impl { inline static typename ei_plain_matrix_type::type run(const VectorLhs& lhs, const VectorRhs& rhs) diff --git a/Eigen/src/Geometry/Quaternion.h b/Eigen/src/Geometry/Quaternion.h index 78365cd38..5e52d5b5a 100644 --- a/Eigen/src/Geometry/Quaternion.h +++ b/Eigen/src/Geometry/Quaternion.h @@ -1,7 +1,7 @@ // This file is part of Eigen, a lightweight C++ template library // for linear algebra. // -// Copyright (C) 2008 Gael Guennebaud +// Copyright (C) 2008-2010 Gael Guennebaud // Copyright (C) 2009 Mathieu Gautier // // Eigen is free software; you can redistribute it and/or @@ -180,6 +180,10 @@ public: return typename ei_cast_return_type >::type( coeffs().template cast()); } + +#ifdef EIGEN_QUATERNIONBASE_PLUGIN +# include EIGEN_QUATERNIONBASE_PLUGIN +#endif }; /*************************************************************************** @@ -395,7 +399,8 @@ QuaternionBase::_transformVector(Vector3 v) const // It appears to be much faster than the common algorithm found // in the litterature (30 versus 39 flops). It also requires two // Vector3 as temporaries. - Vector3 uv = Scalar(2) * this->vec().cross(v); + Vector3 uv = this->vec().cross(v); + uv += uv; return v + this->w() * uv + this->vec().cross(uv); } diff --git a/Eigen/src/Geometry/arch/Geometry_SSE.h b/Eigen/src/Geometry/arch/Geometry_SSE.h index 798d81c91..7d82be694 100644 --- a/Eigen/src/Geometry/arch/Geometry_SSE.h +++ b/Eigen/src/Geometry/arch/Geometry_SSE.h @@ -54,8 +54,8 @@ struct ei_cross3_impl inline static typename ei_plain_matrix_type::type run(const VectorLhs& lhs, const VectorRhs& rhs) { - __m128 a = lhs.coeffs().packet(0); - __m128 b = rhs.coeffs().packet(0); + __m128 a = lhs.template packet(0); + __m128 b = rhs.template packet(0); __m128 mul1=_mm_mul_ps(ei_vec4f_swizzle1(a,1,2,0,3),ei_vec4f_swizzle1(b,2,0,1,3)); __m128 mul2=_mm_mul_ps(ei_vec4f_swizzle1(a,2,0,1,3),ei_vec4f_swizzle1(b,1,2,0,3)); typename ei_plain_matrix_type::type res; diff --git a/Eigen/src/StlSupport/CMakeLists.txt b/Eigen/src/StlSupport/CMakeLists.txt new file mode 100644 index 000000000..0f094f637 --- /dev/null +++ b/Eigen/src/StlSupport/CMakeLists.txt @@ -0,0 +1,6 @@ +FILE(GLOB Eigen_StlSupport_SRCS "*.h") + +INSTALL(FILES + ${Eigen_StlSupport_SRCS} + DESTINATION ${INCLUDE_INSTALL_DIR}/Eigen/src/StlSupport COMPONENT Devel + ) diff --git a/bench/geometry.cpp b/bench/geometry.cpp new file mode 100644 index 000000000..b187a515f --- /dev/null +++ b/bench/geometry.cpp @@ -0,0 +1,126 @@ + +#include +#include +#include + +using namespace std; +using namespace Eigen; + +#ifndef SCALAR +#define SCALAR float +#endif + +#ifndef SIZE +#define SIZE 8 +#endif + +typedef SCALAR Scalar; +typedef NumTraits::Real RealScalar; +typedef Matrix A; +typedef Matrix B; +typedef Matrix C; +typedef Matrix M; + +template +EIGEN_DONT_INLINE void transform(const Transformation& t, Data& data) +{ + EIGEN_ASM_COMMENT("begin"); + data = t * data; + EIGEN_ASM_COMMENT("end"); +} + +template +EIGEN_DONT_INLINE void transform(const Quaternion& t, Data& data) +{ + EIGEN_ASM_COMMENT("begin quat"); + for(int i=0;i struct ToRotationMatrixWrapper +{ + enum {Dim = T::Dim}; + typedef typename T::Scalar Scalar; + ToRotationMatrixWrapper(const T& o) : object(o) {} + T object; +}; + +template +EIGEN_DONT_INLINE void transform(const ToRotationMatrixWrapper& t, Data& data) +{ + EIGEN_ASM_COMMENT("begin quat via mat"); + data = t.object.toRotationMatrix() * data; + EIGEN_ASM_COMMENT("end quat via mat"); +} + +template +EIGEN_DONT_INLINE void transform(const Transform& t, Data& data) +{ + data = (t * data.colwise().homogeneous()).template block(0,0); +} + +template struct get_dim { enum { Dim = T::Dim }; }; +template +struct get_dim > { enum { Dim = R }; }; + +template +struct bench_impl +{ + static EIGEN_DONT_INLINE void run(const Transformation& t) + { + Matrix::Dim,N> data; + data.setRandom(); + bench_impl::run(t); + BenchTimer timer; + BENCH(timer,10,100000,transform(t,data)); + cout.width(9); + cout << timer.best() << " "; + } +}; + + +template +struct bench_impl +{ + static EIGEN_DONT_INLINE void run(const Transformation&) {} +}; + +template +EIGEN_DONT_INLINE void bench(const std::string& msg, const Transformation& t) +{ + cout << msg << " "; + bench_impl::run(t); + std::cout << "\n"; +} + +int main(int argc, char ** argv) +{ + Matrix mat34; mat34.setRandom(); + Transform iso3(mat34); + Transform aff3(mat34); + Transform caff3(mat34); + Transform proj3(mat34); + Quaternion quat;quat.setIdentity(); + ToRotationMatrixWrapper > quatmat(quat); + Matrix mat33; mat33.setRandom(); + + cout.precision(4); + std::cout + << "N "; + for(int i=0;i void basicStuff(const MatrixType& m) typedef typename MatrixType::Index Index; typedef typename MatrixType::Scalar Scalar; typedef Matrix VectorType; + typedef Matrix SquareMatrixType; Index rows = m.rows(); Index cols = m.cols(); @@ -47,6 +48,7 @@ template void basicStuff(const MatrixType& m) VectorType v1 = VectorType::Random(rows), v2 = VectorType::Random(rows), vzero = VectorType::Zero(rows); + SquareMatrixType sm1 = SquareMatrixType::Random(rows,rows), sm2(rows,rows); Scalar x = ei_random(); @@ -121,6 +123,27 @@ template void basicStuff(const MatrixType& m) m1 = m2; VERIFY(m1==m2); VERIFY(!(m1!=m2)); + + // check automatic transposition + sm2.setZero(); + for(typename MatrixType::Index i=0;i void basicStuffComplex(const MatrixType& m) @@ -177,14 +200,14 @@ void test_basicstuff() for(int i = 0; i < g_repeat; i++) { CALL_SUBTEST_1( basicStuff(Matrix()) ); CALL_SUBTEST_2( basicStuff(Matrix4d()) ); - CALL_SUBTEST_3( basicStuff(MatrixXcf(3, 3)) ); - CALL_SUBTEST_4( basicStuff(MatrixXi(8, 12)) ); - CALL_SUBTEST_5( basicStuff(MatrixXcd(20, 20)) ); + CALL_SUBTEST_3( basicStuff(MatrixXcf(ei_random(1,100), ei_random(1,100))) ); + CALL_SUBTEST_4( basicStuff(MatrixXi(ei_random(1,100), ei_random(1,100))) ); + CALL_SUBTEST_5( basicStuff(MatrixXcd(ei_random(1,100), ei_random(1,100))) ); CALL_SUBTEST_6( basicStuff(Matrix()) ); - CALL_SUBTEST_7( basicStuff(Matrix(10,10)) ); + CALL_SUBTEST_7( basicStuff(Matrix(ei_random(1,100),ei_random(1,100))) ); - CALL_SUBTEST_3( basicStuffComplex(MatrixXcf(21, 17)) ); - CALL_SUBTEST_5( basicStuffComplex(MatrixXcd(2, 3)) ); + CALL_SUBTEST_3( basicStuffComplex(MatrixXcf(ei_random(1,100), ei_random(1,100))) ); + CALL_SUBTEST_5( basicStuffComplex(MatrixXcd(ei_random(1,100), ei_random(1,100))) ); } CALL_SUBTEST_2(casting()); diff --git a/test/cholesky.cpp b/test/cholesky.cpp index 0edf9a793..46140bb11 100644 --- a/test/cholesky.cpp +++ b/test/cholesky.cpp @@ -168,6 +168,21 @@ template void cholesky(const MatrixType& m) } } + // test some special use cases of SelfCwiseBinaryOp: + MatrixType m1 = MatrixType::Random(rows,cols), m2(rows,cols); + m2 = m1; + m2 += symmLo.template selfadjointView().llt().solve(matB); + VERIFY_IS_APPROX(m2, m1 + symmLo.template selfadjointView().llt().solve(matB)); + m2 = m1; + m2 -= symmLo.template selfadjointView().llt().solve(matB); + VERIFY_IS_APPROX(m2, m1 - symmLo.template selfadjointView().llt().solve(matB)); + m2 = m1; + m2.noalias() += symmLo.template selfadjointView().llt().solve(matB); + VERIFY_IS_APPROX(m2, m1 + symmLo.template selfadjointView().llt().solve(matB)); + m2 = m1; + m2.noalias() -= symmLo.template selfadjointView().llt().solve(matB); + VERIFY_IS_APPROX(m2, m1 - symmLo.template selfadjointView().llt().solve(matB)); + } template void cholesky_cplx(const MatrixType& m) diff --git a/test/linearstructure.cpp b/test/linearstructure.cpp index a0b8982dd..b5c58bdaa 100644 --- a/test/linearstructure.cpp +++ b/test/linearstructure.cpp @@ -27,7 +27,7 @@ template void linearStructure(const MatrixType& m) { /* this test covers the following files: - Sum.h Difference.h Opposite.h ScalarMultiple.h + CwiseUnaryOp.h, CwiseBinaryOp.h, SelfCwiseBinaryOp.h */ typedef typename MatrixType::Index Index; typedef typename MatrixType::Scalar Scalar; diff --git a/unsupported/Eigen/CholmodSupport b/unsupported/Eigen/CholmodSupport index a5700a799..10fa169f2 100644 --- a/unsupported/Eigen/CholmodSupport +++ b/unsupported/Eigen/CholmodSupport @@ -5,11 +5,9 @@ #include "src/Core/util/DisableMSVCWarnings.h" -#ifdef EIGEN_CHOLMOD_SUPPORT - extern "C" { - #include - } -#endif +extern "C" { + #include +} namespace Eigen { diff --git a/unsupported/Eigen/TaucsSupport b/unsupported/Eigen/TaucsSupport index 585a81810..dca2c095b 100644 --- a/unsupported/Eigen/TaucsSupport +++ b/unsupported/Eigen/TaucsSupport @@ -39,9 +39,7 @@ namespace Eigen { * \endcode */ -#ifdef EIGEN_TAUCS_SUPPORT # include "src/SparseExtra/TaucsSupport.h" -#endif } // namespace Eigen