From 3411146984316c97f56543333a46f47aeb7f9d35 Mon Sep 17 00:00:00 2001 From: Sergey Sharybin Date: Fri, 21 Mar 2014 16:04:53 +0600 Subject: Update Eigen to version 3.2.1 Main purpose of this is to have SparseLU solver which we can use now as a replacement to opennl library. --- extern/Eigen3/Eigen/src/Core/GeneralProduct.h | 100 ++++++++++++++++---------- 1 file changed, 61 insertions(+), 39 deletions(-) (limited to 'extern/Eigen3/Eigen/src/Core/GeneralProduct.h') diff --git a/extern/Eigen3/Eigen/src/Core/GeneralProduct.h b/extern/Eigen3/Eigen/src/Core/GeneralProduct.h index bfc2a67b12f..2a59d94645e 100644 --- a/extern/Eigen3/Eigen/src/Core/GeneralProduct.h +++ b/extern/Eigen3/Eigen/src/Core/GeneralProduct.h @@ -222,7 +222,29 @@ class GeneralProduct ***********************************************************************/ namespace internal { -template struct outer_product_selector; + +// Column major +template +EIGEN_DONT_INLINE void outer_product_selector_run(const ProductType& prod, Dest& dest, const Func& func, const false_type&) +{ + typedef typename Dest::Index Index; + // FIXME make sure lhs is sequentially stored + // FIXME not very good if rhs is real and lhs complex while alpha is real too + const Index cols = dest.cols(); + for (Index j=0; j +EIGEN_DONT_INLINE void outer_product_selector_run(const ProductType& prod, Dest& dest, const Func& func, const true_type&) { + typedef typename Dest::Index Index; + // FIXME make sure rhs is sequentially stored + // FIXME not very good if lhs is real and rhs complex while alpha is real too + const Index rows = dest.rows(); + for (Index i=0; i struct traits > @@ -235,6 +257,8 @@ template class GeneralProduct : public ProductBase, Lhs, Rhs> { + template struct IsRowMajor : internal::conditional<(int(T::Flags)&RowMajorBit), internal::true_type, internal::false_type>::type {}; + public: EIGEN_PRODUCT_PUBLIC_INTERFACE(GeneralProduct) @@ -243,41 +267,39 @@ class GeneralProduct EIGEN_STATIC_ASSERT((internal::is_same::value), YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY) } - - template void scaleAndAddTo(Dest& dest, Scalar alpha) const - { - internal::outer_product_selector<(int(Dest::Flags)&RowMajorBit) ? RowMajor : ColMajor>::run(*this, dest, alpha); + + struct set { template void operator()(const Dst& dst, const Src& src) const { dst.const_cast_derived() = src; } }; + struct add { template void operator()(const Dst& dst, const Src& src) const { dst.const_cast_derived() += src; } }; + struct sub { template void operator()(const Dst& dst, const Src& src) const { dst.const_cast_derived() -= src; } }; + struct adds { + Scalar m_scale; + adds(const Scalar& s) : m_scale(s) {} + template void operator()(const Dst& dst, const Src& src) const { + dst.const_cast_derived() += m_scale * src; + } + }; + + template + inline void evalTo(Dest& dest) const { + internal::outer_product_selector_run(*this, dest, set(), IsRowMajor()); + } + + template + inline void addTo(Dest& dest) const { + internal::outer_product_selector_run(*this, dest, add(), IsRowMajor()); } -}; - -namespace internal { -template<> struct outer_product_selector { - template - static EIGEN_DONT_INLINE void run(const ProductType& prod, Dest& dest, typename ProductType::Scalar alpha) { - typedef typename Dest::Index Index; - // FIXME make sure lhs is sequentially stored - // FIXME not very good if rhs is real and lhs complex while alpha is real too - const Index cols = dest.cols(); - for (Index j=0; j + inline void subTo(Dest& dest) const { + internal::outer_product_selector_run(*this, dest, sub(), IsRowMajor()); + } -template<> struct outer_product_selector { - template - static EIGEN_DONT_INLINE void run(const ProductType& prod, Dest& dest, typename ProductType::Scalar alpha) { - typedef typename Dest::Index Index; - // FIXME make sure rhs is sequentially stored - // FIXME not very good if lhs is real and rhs complex while alpha is real too - const Index rows = dest.rows(); - for (Index i=0; i void scaleAndAddTo(Dest& dest, const Scalar& alpha) const + { + internal::outer_product_selector_run(*this, dest, adds(alpha), IsRowMajor()); + } }; -} // end namespace internal - /*********************************************************************** * Implementation of General Matrix Vector Product ***********************************************************************/ @@ -311,7 +333,7 @@ class GeneralProduct typedef typename Lhs::Scalar LhsScalar; typedef typename Rhs::Scalar RhsScalar; - GeneralProduct(const Lhs& lhs, const Rhs& rhs) : Base(lhs,rhs) + GeneralProduct(const Lhs& a_lhs, const Rhs& a_rhs) : Base(a_lhs,a_rhs) { // EIGEN_STATIC_ASSERT((internal::is_same::value), // YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY) @@ -320,7 +342,7 @@ class GeneralProduct enum { Side = Lhs::IsVectorAtCompileTime ? OnTheLeft : OnTheRight }; typedef typename internal::conditional::type MatrixType; - template void scaleAndAddTo(Dest& dst, Scalar alpha) const + template void scaleAndAddTo(Dest& dst, const Scalar& alpha) const { eigen_assert(m_lhs.rows() == dst.rows() && m_rhs.cols() == dst.cols()); internal::gemv_selector struct gemv_selector { template - static void run(const ProductType& prod, Dest& dest, typename ProductType::Scalar alpha) + static void run(const ProductType& prod, Dest& dest, const typename ProductType::Scalar& alpha) { Transpose destT(dest); enum { OtherStorageOrder = StorageOrder == RowMajor ? ColMajor : RowMajor }; @@ -384,7 +406,7 @@ struct gemv_static_vector_if template<> struct gemv_selector { template - static inline void run(const ProductType& prod, Dest& dest, typename ProductType::Scalar alpha) + static inline void run(const ProductType& prod, Dest& dest, const typename ProductType::Scalar& alpha) { typedef typename ProductType::Index Index; typedef typename ProductType::LhsScalar LhsScalar; @@ -413,7 +435,7 @@ template<> struct gemv_selector gemv_static_vector_if static_dest; - bool alphaIsCompatible = (!ComplexByReal) || (imag(actualAlpha)==RealScalar(0)); + bool alphaIsCompatible = (!ComplexByReal) || (numext::imag(actualAlpha)==RealScalar(0)); bool evalToDest = EvalToDestAtCompileTime && alphaIsCompatible; RhsScalar compatibleAlpha = get_factor::run(actualAlpha); @@ -457,7 +479,7 @@ template<> struct gemv_selector template<> struct gemv_selector { template - static void run(const ProductType& prod, Dest& dest, typename ProductType::Scalar alpha) + static void run(const ProductType& prod, Dest& dest, const typename ProductType::Scalar& alpha) { typedef typename ProductType::LhsScalar LhsScalar; typedef typename ProductType::RhsScalar RhsScalar; @@ -508,7 +530,7 @@ template<> struct gemv_selector template<> struct gemv_selector { template - static void run(const ProductType& prod, Dest& dest, typename ProductType::Scalar alpha) + static void run(const ProductType& prod, Dest& dest, const typename ProductType::Scalar& alpha) { typedef typename Dest::Index Index; // TODO makes sure dest is sequentially stored in memory, otherwise use a temp @@ -521,7 +543,7 @@ template<> struct gemv_selector template<> struct gemv_selector { template - static void run(const ProductType& prod, Dest& dest, typename ProductType::Scalar alpha) + static void run(const ProductType& prod, Dest& dest, const typename ProductType::Scalar& alpha) { typedef typename Dest::Index Index; // TODO makes sure rhs is sequentially stored in memory, otherwise use a temp -- cgit v1.2.3