From 34ca81b1bfe77b11effdd463099b30d12b6466c9 Mon Sep 17 00:00:00 2001
From: Gael Guennebaud <g.gael@free.fr>
Date: Mon, 2 Dec 2013 16:37:58 +0100
Subject: [PATCH] Add direct assignment of products

---
 Eigen/src/Core/AssignEvaluator.h   |  5 ++++
 Eigen/src/Core/ProductEvaluators.h | 40 ++++++++++++++++++++++++++++--
 test/evaluators.cpp                | 15 +++++++++--
 3 files changed, 56 insertions(+), 4 deletions(-)

diff --git a/Eigen/src/Core/AssignEvaluator.h b/Eigen/src/Core/AssignEvaluator.h
index 3babc16ba..6d3739407 100644
--- a/Eigen/src/Core/AssignEvaluator.h
+++ b/Eigen/src/Core/AssignEvaluator.h
@@ -679,6 +679,11 @@ void call_assignment(const NoAlias<Dst,StorageBase>& dst, const Src& src, const
 {
   Assignment<Dst,Src,Func>::run(dst.expression(), src, func);
 }
+template<typename Dst, template <typename> class StorageBase, typename Src, typename Func>
+void call_assignment(NoAlias<Dst,StorageBase>& dst, const Src& src, const Func& func)
+{
+  Assignment<Dst,Src,Func>::run(dst.expression(), src, func);
+}
 
 // Generic Dense to Dense assignment
 template< typename DstXprType, typename SrcXprType, typename Functor, typename Scalar>
diff --git a/Eigen/src/Core/ProductEvaluators.h b/Eigen/src/Core/ProductEvaluators.h
index 26145ceff..d0dba1644 100644
--- a/Eigen/src/Core/ProductEvaluators.h
+++ b/Eigen/src/Core/ProductEvaluators.h
@@ -45,6 +45,12 @@ struct evaluator<Product<Lhs, Rhs, Options> >
 template<typename Lhs, typename Rhs, int ProductType = internal::product_type<Lhs,Rhs>::value>
 struct dense_product_impl;
 
+template<typename Lhs, typename Rhs>
+struct evaluator_traits<Product<Lhs, Rhs, DefaultProduct> > 
+ : evaluator_traits_base<Product<Lhs, Rhs, DefaultProduct> >
+{
+  enum { AssumeAliasing = 1 };
+};
 
 // The evaluator for default dense products creates a temporary and call dense_product_impl
 template<typename Lhs, typename Rhs, int ProductTag>
@@ -58,9 +64,7 @@ struct product_evaluator<Product<Lhs, Rhs, DefaultProduct>, ProductTag, DenseSha
   product_evaluator(const XprType& xpr)
     : m_result(xpr.rows(), xpr.cols())
   {
-    
     ::new (static_cast<Base*>(this)) Base(m_result);
-    
     dense_product_impl<Lhs, Rhs>::evalTo(m_result, xpr.lhs(), xpr.rhs());
   }
   
@@ -68,6 +72,38 @@ protected:
   PlainObject m_result;
 };
 
+// Dense = Product
+template< typename DstXprType, typename Lhs, typename Rhs, typename Scalar>
+struct Assignment<DstXprType, Product<Lhs,Rhs,DefaultProduct>, internal::assign_op<Scalar>, Dense2Dense, Scalar>
+{
+  typedef Product<Lhs,Rhs,DefaultProduct> SrcXprType;
+  static void run(DstXprType &dst, const SrcXprType &src, const internal::assign_op<Scalar> &)
+  {
+    dense_product_impl<Lhs, Rhs>::evalTo(dst, src.lhs(), src.rhs());
+  }
+};
+
+// Dense += Product
+template< typename DstXprType, typename Lhs, typename Rhs, typename Scalar>
+struct Assignment<DstXprType, Product<Lhs,Rhs,DefaultProduct>, internal::add_assign_op<Scalar>, Dense2Dense, Scalar>
+{
+  typedef Product<Lhs,Rhs,DefaultProduct> SrcXprType;
+  static void run(DstXprType &dst, const SrcXprType &src, const internal::add_assign_op<Scalar> &)
+  {
+    dense_product_impl<Lhs, Rhs>::addTo(dst, src.lhs(), src.rhs());
+  }
+};
+
+// Dense -= Product
+template< typename DstXprType, typename Lhs, typename Rhs, typename Scalar>
+struct Assignment<DstXprType, Product<Lhs,Rhs,DefaultProduct>, internal::sub_assign_op<Scalar>, Dense2Dense, Scalar>
+{
+  typedef Product<Lhs,Rhs,DefaultProduct> SrcXprType;
+  static void run(DstXprType &dst, const SrcXprType &src, const internal::sub_assign_op<Scalar> &)
+  {
+    dense_product_impl<Lhs, Rhs>::subTo(dst, src.lhs(), src.rhs());
+  }
+};
 
 template<typename Lhs, typename Rhs>
 struct dense_product_impl<Lhs,Rhs,InnerProduct>
diff --git a/test/evaluators.cpp b/test/evaluators.cpp
index 7a20014dd..c07146260 100644
--- a/test/evaluators.cpp
+++ b/test/evaluators.cpp
@@ -47,14 +47,14 @@ namespace Eigen {
   void add_assign_using_evaluator(const DstXprType& dst, const SrcXprType& src)
   {
     typedef typename DstXprType::Scalar Scalar;
-    call_assignment(dst.const_cast_derived(), src.derived(), internal::add_assign_op<Scalar>());
+    call_assignment(const_cast<DstXprType&>(dst), src.derived(), internal::add_assign_op<Scalar>());
   }
 
   template<typename DstXprType, typename SrcXprType>
   void subtract_assign_using_evaluator(const DstXprType& dst, const SrcXprType& src)
   {
     typedef typename DstXprType::Scalar Scalar;
-    call_assignment(dst.const_cast_derived(), src.derived(), internal::sub_assign_op<Scalar>());
+    call_assignment(const_cast<DstXprType&>(dst), src.derived(), internal::sub_assign_op<Scalar>());
   }
 
   template<typename DstXprType, typename SrcXprType>
@@ -151,6 +151,17 @@ void test_evaluators()
     c = a*a;
     copy_using_evaluator(a, prod(a,a));
     VERIFY_IS_APPROX(a,c);
+    
+    // check compound assignment of products
+    d = c;
+    add_assign_using_evaluator(c.noalias(), prod(a,b));
+    d.noalias() += a*b;
+    VERIFY_IS_APPROX(c, d);
+    
+    d = c;
+    subtract_assign_using_evaluator(c.noalias(), prod(a,b));
+    d.noalias() -= a*b;
+    VERIFY_IS_APPROX(c, d);
   }
   
   {