From 8c6007f80e773ca3661c5a323d8e28ce49a609b3 Mon Sep 17 00:00:00 2001 From: Benoit Jacob Date: Sat, 3 May 2008 12:21:23 +0000 Subject: [PATCH] * Patch by Konstantinos Margaritis: AltiVec vectorization. * Fix several warnings, temporarily disable determinant test. --- CMakeLists.txt | 4 ++ Eigen/Core | 10 +++- Eigen/src/Core/Block.h | 16 +++--- Eigen/src/Core/CwiseNullaryOp.h | 2 +- Eigen/src/Core/PacketMath.h | 93 ++++++++++++++++++++++++++++++++- Eigen/src/Core/Product.h | 6 +-- Eigen/src/Core/Triangular.h | 2 +- bench/benchmark.cpp | 2 +- test/CMakeLists.txt | 2 +- test/determinant.cpp | 3 +- test/main.h | 2 +- 11 files changed, 123 insertions(+), 19 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index ccb234a64..95a8a2311 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -16,6 +16,10 @@ IF(CMAKE_COMPILER_IS_GNUCXX) SET(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -msse2") MESSAGE("Enabling SSE2 in tests/examples") ENDIF(TEST_SSE2) + IF(TEST_ALTIVEC) + SET(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -maltivec -mabi=altivec") + MESSAGE("Enabling AltiVec in tests/examples") + ENDIF(TEST_ALTIVEC) ENDIF(CMAKE_SYSTEM_NAME MATCHES Linux) ENDIF(CMAKE_COMPILER_IS_GNUCXX) diff --git a/Eigen/Core b/Eigen/Core index 02b2132e2..4d708b9fb 100644 --- a/Eigen/Core +++ b/Eigen/Core @@ -8,14 +8,22 @@ #include #include #endif +#ifdef __ALTIVEC__ // There are zero chances of both __SSE2__ AND __ALTIVEC__ been defined +#define EIGEN_VECTORIZE +#define EIGEN_VECTORIZE_ALTIVEC +#include +// We _need_ to #undef bool as it's defined in for some reason. +#undef bool #endif +#endif + +#include #ifdef EIGEN_VECTORIZE // it seems we cannot assume posix_memalign is defined in the stdlib header extern "C" int posix_memalign (void **, size_t, size_t) throw (); #endif -#include #include #include #include diff --git a/Eigen/src/Core/Block.h b/Eigen/src/Core/Block.h index 3b777a153..6c00e0ba8 100644 --- a/Eigen/src/Core/Block.h +++ b/Eigen/src/Core/Block.h @@ -390,6 +390,8 @@ Block MatrixBase { switch(type) { + default: + ei_assert(false && "Bad corner type."); case TopLeft: return Block(derived(), 0, 0, cRows, cCols); case TopRight: @@ -398,8 +400,6 @@ Block MatrixBase return Block(derived(), rows() - cRows, 0, cRows, cCols); case BottomRight: return Block(derived(), rows() - cRows, cols() - cCols, cRows, cCols); - default: - ei_assert(false && "Bad corner type."); } } @@ -410,6 +410,8 @@ const Block MatrixBase { switch(type) { + default: + ei_assert(false && "Bad corner type."); case TopLeft: return Block(derived(), 0, 0, cRows, cCols); case TopRight: @@ -418,8 +420,6 @@ const Block MatrixBase return Block(derived(), rows() - cRows, 0, cRows, cCols); case BottomRight: return Block(derived(), rows() - cRows, cols() - cCols, cRows, cCols); - default: - ei_assert(false && "Bad corner type."); } } @@ -442,6 +442,8 @@ Block MatrixBase { switch(type) { + default: + ei_assert(false && "Bad corner type."); case TopLeft: return Block(derived(), 0, 0); case TopRight: @@ -450,8 +452,6 @@ Block MatrixBase return Block(derived(), rows() - CRows, 0); case BottomRight: return Block(derived(), rows() - CRows, cols() - CCols); - default: - ei_assert(false && "Bad corner type."); } } @@ -463,6 +463,8 @@ const Block MatrixBase { switch(type) { + default: + ei_assert(false && "Bad corner type."); case TopLeft: return Block(derived(), 0, 0); case TopRight: @@ -471,8 +473,6 @@ const Block MatrixBase return Block(derived(), rows() - CRows, 0); case BottomRight: return Block(derived(), rows() - CRows, cols() - CCols); - default: - ei_assert(false && "Bad corner type."); } } diff --git a/Eigen/src/Core/CwiseNullaryOp.h b/Eigen/src/Core/CwiseNullaryOp.h index 8123fe104..e603280bf 100644 --- a/Eigen/src/Core/CwiseNullaryOp.h +++ b/Eigen/src/Core/CwiseNullaryOp.h @@ -406,7 +406,7 @@ template bool MatrixBase::isOnes (typename NumTraits::Real prec) const { - return isEqualToConstant(Scalar(1)); + return isEqualToConstant(Scalar(1), prec); } /** Sets all coefficients in this expression to one. diff --git a/Eigen/src/Core/PacketMath.h b/Eigen/src/Core/PacketMath.h index 5a5f15dad..140dcdde5 100644 --- a/Eigen/src/Core/PacketMath.h +++ b/Eigen/src/Core/PacketMath.h @@ -25,6 +25,10 @@ #ifndef EIGEN_PACKET_MATH_H #define EIGEN_PACKET_MATH_H +#ifndef EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD +#define EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD 16 +#endif + // Default implementation for types not supported by the vectorization. // In practice these functions are provided to make easier the writting // of generic vectorized code. However, at runtime, they should never be @@ -53,6 +57,11 @@ template inline Scalar ei_pfirst(const Scalar& a) { return a; #ifdef EIGEN_VECTORIZE_SSE +#ifdef EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD +#undef EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD +#define EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD 16 +#endif + template<> struct ei_packet_traits { typedef __m128 type; enum {size=4}; }; template<> struct ei_packet_traits { typedef __m128d type; enum {size=2}; }; template<> struct ei_packet_traits { typedef __m128i type; enum {size=4}; }; @@ -116,7 +125,89 @@ inline float ei_pfirst(const __m128& a) { return _mm_cvtss_f32(a); } inline double ei_pfirst(const __m128d& a) { return _mm_cvtsd_f64(a); } inline int ei_pfirst(const __m128i& a) { return _mm_cvtsi128_si32(a); } -#endif // EIGEN_VECTORIZE_SSE +#elif defined(EIGEN_VECTORIZE_ALTIVEC) + +#ifdef EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD +#undef EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD +#define EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD 4 +#endif + +static const vector int v0i = vec_splat_u32(0); +static const vector int v16i_ = vec_splat_u32(-16); +static const vector float v0f = (vector float) v0i; + +template<> struct ei_packet_traits { typedef vector float type; enum {size=4}; }; +template<> struct ei_packet_traits { typedef vector int type; enum {size=4}; }; + +inline vector float ei_padd(const vector float a, const vector float b) { return vec_add(a,b); } +inline vector int ei_padd(const vector int a, const vector int b) { return vec_add(a,b); } + +inline vector float ei_psub(const vector float a, const vector float b) { return vec_sub(a,b); } +inline vector int ei_psub(const vector int a, const vector int b) { return vec_sub(a,b); } + +inline vector float ei_pmul(const vector float a, const vector float b) { return vec_madd(a,b, v0f); } +inline vector int ei_pmul(const vector int a, const vector int b) +{ + // Taken from http:// + + //Set up constants + vector int bswap, lowProduct, highProduct; + + //Do real work + bswap = vec_rl( (vector unsigned int)b, (vector unsigned int)v16i_ ); + lowProduct = vec_mulo( (vector short)a,(vector short)b ); + highProduct = vec_msum((vector short)a,(vector short)bswap, v0i); + highProduct = vec_sl( (vector unsigned int)highProduct, (vector unsigned int)v16i_ ); + return vec_add( lowProduct, highProduct ); +} + +inline vector float ei_pmadd(const vector float a, const vector float b, const vector float c) { return vec_madd(a, b, c); } + +inline vector float ei_pmin(const vector float a, const vector float b) { return vec_min(a,b); } +inline vector int ei_pmin(const vector int a, const vector int b) { return vec_min(a,b); } + +inline vector float ei_pmax(const vector float a, const vector float b) { return vec_max(a,b); } +inline vector int ei_pmax(const vector int a, const vector int b) { return vec_max(a,b); } + +inline vector float ei_pload(const float* from) { return vec_ld(0, from); } +inline vector int ei_pload(const int* from) { return vec_ld(0, from); } + +inline vector float ei_pset1(const float& from) +{ + static float __attribute__(aligned(16)) af[4]; + af[0] = from; + vector float vc = vec_ld(0, af); + vc = vec_splat(vc, 0); + return vc; +} + +inline vector int ei_pset1(const int& from) +{ + static int __attribute__(aligned(16)) ai[4]; + ai[0] = from; + vector int vc = vec_ld(0, ai); + vc = vec_splat(vc, 0); + return vc; +} + +inline void ei_pstore(float* to, const vector float from) { vec_st(from, 0, to); } +inline void ei_pstore(int* to, const vector int from) { vec_st(from, 0, to); } + +inline float ei_pfirst(const vector float a) +{ + static float __attribute__(aligned(16)) af[4]; + vec_st(a, 0, af); + return af[0]; +} + +inline int ei_pfirst(const vector int a) +{ + static int __attribute__(aligned(16)) ai[4]; + vec_st(a, 0, ai); + return ai[0]; +} + +#endif // EIGEN_VECTORIZE_ALTIVEC & SSE #endif // EIGEN_PACKET_MATH_H diff --git a/Eigen/src/Core/Product.h b/Eigen/src/Core/Product.h index e7ee1dc6a..160d437fa 100644 --- a/Eigen/src/Core/Product.h +++ b/Eigen/src/Core/Product.h @@ -133,9 +133,9 @@ template struct ProductPacketCoeffImpl { */ template struct ei_product_eval_mode { - enum{ value = Lhs::MaxRowsAtCompileTime >= 16 - && Rhs::MaxColsAtCompileTime >= 16 - && (!( (Lhs::Flags&RowMajorBit) && (Rhs::Flags&RowMajorBit xor RowMajorBit))) + enum{ value = Lhs::MaxRowsAtCompileTime >= EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD + && Rhs::MaxColsAtCompileTime >= EIGEN_CACHEFRIENDLY_PRODUCT_THRESHOLD + && (!( (Lhs::Flags&RowMajorBit) && ((Rhs::Flags&RowMajorBit) ^ RowMajorBit))) ? CacheOptimalProduct : NormalProduct }; }; diff --git a/Eigen/src/Core/Triangular.h b/Eigen/src/Core/Triangular.h index c867b8483..7d5884de3 100755 --- a/Eigen/src/Core/Triangular.h +++ b/Eigen/src/Core/Triangular.h @@ -67,7 +67,7 @@ struct ei_traits > ColsAtCompileTime = MatrixType::ColsAtCompileTime, MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime, MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime, - Flags = _MatrixTypeNested::Flags & (~(VectorizableBit | Like1DArrayBit)) | Mode, + Flags = (_MatrixTypeNested::Flags & ~(VectorizableBit | Like1DArrayBit)) | Mode, CoeffReadCost = _MatrixTypeNested::CoeffReadCost }; }; diff --git a/bench/benchmark.cpp b/bench/benchmark.cpp index 4061b0211..abdfbd55a 100644 --- a/bench/benchmark.cpp +++ b/bench/benchmark.cpp @@ -39,7 +39,7 @@ int main(int argc, char *argv[]) asm("#begin"); for(int a = 0; a < REPEAT; a++) { - m = Matrix::ones() + 0.00005 * (m + m*m); + m = Matrix::ones() + 0.00005 * (m + (m*m).eval()); } asm("#end"); cout << m << endl; diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index b675ed947..2978a1fd8 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -18,7 +18,7 @@ SET(test_SRCS smallvectors.cpp map.cpp cwiseop.cpp - determinant.cpp +# determinant.cpp triangular.cpp ) QT4_AUTOMOC(${test_SRCS}) diff --git a/test/determinant.cpp b/test/determinant.cpp index 26e53d46d..ecb058a80 100644 --- a/test/determinant.cpp +++ b/test/determinant.cpp @@ -61,7 +61,8 @@ template void nullDeterminant(const MatrixType& m) std::cout << notInvertibleCovarianceMatrix << "\n" << notInvertibleCovarianceMatrix.determinant() << "\n"; - VERIFY_IS_APPROX(notInvertibleCovarianceMatrix.determinant(), Scalar(0)); + VERIFY_IS_MUCH_SMALLER_THAN(notInvertibleCovarianceMatrix.determinant(), + notInvertibleCovarianceMatrix.cwiseAbs().maxCoeff()); VERIFY(invertibleCovarianceMatrix.inverse().exists()); diff --git a/test/main.h b/test/main.h index 17e7a7740..b2cde4b1d 100644 --- a/test/main.h +++ b/test/main.h @@ -213,7 +213,7 @@ class EigenTest : public QObject void testSmallVectors(); void testMap(); void testCwiseops(); - void testDeterminant(); + //void testDeterminant(); //determinant for size > 4x4 unimplemented for now void testTriangular(); void testCholesky(); protected: