* update sparse module wrt new diagonal matrix impl

* fix a bug is SparseMatrix

4 files changed, 58 insertions, 30 deletions

diff --git a/Eigen/src/Sparse/SparseDiagonalProduct.h b/Eigen/src/Sparse/SparseDiagonalProduct.h
index 3a148b09e..5fb149a2c 100644
--- a/Eigen/src/Sparse/SparseDiagonalProduct.h
+++ b/Eigen/src/Sparse/SparseDiagonalProduct.h
@@ -25,8 +25,9 @@
 #ifndef EIGEN_SPARSE_DIAGONAL_PRODUCT_H
 #define EIGEN_SPARSE_DIAGONAL_PRODUCT_H
 
-// the product a diagonal matrix with a sparse matrix can be easily
-// implemented using expression template. We have two very different cases:
+// The product of a diagonal matrix with a sparse matrix can be easily
+// implemented using expression template.
+// We have two consider very different cases:
 // 1 - diag * row-major sparse
 //     => each inner vector <=> scalar * sparse vector product
 //     => so we can reuse CwiseUnaryOp::InnerIterator
@@ -37,13 +38,21 @@
 // The two other cases are symmetric.
 
 template<typename Lhs, typename Rhs>
-struct ei_traits<SparseDiagonalProduct<Lhs, Rhs> > : ei_traits<SparseProduct<Lhs, Rhs, DiagonalProduct> >
+struct ei_traits<SparseDiagonalProduct<Lhs, Rhs> >
 {
   typedef typename ei_cleantype<Lhs>::type _Lhs;
   typedef typename ei_cleantype<Rhs>::type _Rhs;
+  typedef typename _Lhs::Scalar Scalar;
   enum {
+    RowsAtCompileTime = _Lhs::RowsAtCompileTime,
+    ColsAtCompileTime = _Rhs::ColsAtCompileTime,
+
+    MaxRowsAtCompileTime = _Lhs::MaxRowsAtCompileTime,
+    MaxColsAtCompileTime = _Rhs::MaxColsAtCompileTime,
+    
     SparseFlags = ei_is_diagonal<_Lhs>::ret ? int(_Rhs::Flags) : int(_Lhs::Flags),
-    Flags = SparseBit | (SparseFlags&RowMajorBit)
+    Flags = SparseBit | (SparseFlags&RowMajorBit),
+    CoeffReadCost = Dynamic
   };
 };
 
@@ -51,11 +60,16 @@ enum {SDP_IsDiagonal, SDP_IsSparseRowMajor, SDP_IsSparseColMajor};
 template<typename Lhs, typename Rhs, typename SparseDiagonalProductType, int RhsMode, int LhsMode> 
 class ei_sparse_diagonal_product_inner_iterator_selector;
 
-template<typename LhsNested, typename RhsNested>
-class SparseDiagonalProduct : public SparseMatrixBase<SparseDiagonalProduct<LhsNested,RhsNested> >, ei_no_assignment_operator
+template<typename Lhs, typename Rhs>
+class SparseDiagonalProduct
+  : public SparseMatrixBase<SparseDiagonalProduct<Lhs,Rhs> >,
+    ei_no_assignment_operator
 {
-    typedef typename ei_traits<SparseDiagonalProduct>::_LhsNested _LhsNested;
-    typedef typename ei_traits<SparseDiagonalProduct>::_RhsNested _RhsNested;
+    typedef typename Lhs::Nested LhsNested;
+    typedef typename Rhs::Nested RhsNested;
+
+    typedef typename ei_cleantype<LhsNested>::type _LhsNested;
+    typedef typename ei_cleantype<RhsNested>::type _RhsNested;
     
     enum {
       LhsMode = ei_is_diagonal<_LhsNested>::ret ? SDP_IsDiagonal
@@ -71,8 +85,8 @@ class SparseDiagonalProduct : public SparseMatrixBase<SparseDiagonalProduct<LhsN
     typedef ei_sparse_diagonal_product_inner_iterator_selector
                 <_LhsNested,_RhsNested,SparseDiagonalProduct,LhsMode,RhsMode> InnerIterator;
 
-    template<typename Lhs, typename Rhs>
-    EIGEN_STRONG_INLINE SparseDiagonalProduct(const Lhs& lhs, const Rhs& rhs)
+    template<typename _Lhs, typename _Rhs>
+    EIGEN_STRONG_INLINE SparseDiagonalProduct(const _Lhs& lhs, const _Rhs& rhs)
       : m_lhs(lhs), m_rhs(rhs)
     {
       ei_assert(lhs.cols() == rhs.rows() && "invalid sparse matrix * diagonal matrix product");
@@ -109,12 +123,12 @@ class ei_sparse_diagonal_product_inner_iterator_selector
   : public SparseCwiseBinaryOp<
       ei_scalar_product_op<typename Lhs::Scalar>,
       SparseInnerVectorSet<Rhs,1>,
-      typename Lhs::_CoeffsVectorType>::InnerIterator
+      typename Lhs::DiagonalVectorType>::InnerIterator
 {
     typedef typename SparseCwiseBinaryOp<
       ei_scalar_product_op<typename Lhs::Scalar>,
       SparseInnerVectorSet<Rhs,1>,
-      typename Lhs::_CoeffsVectorType>::InnerIterator Base;
+      typename Lhs::DiagonalVectorType>::InnerIterator Base;
   public:
     inline ei_sparse_diagonal_product_inner_iterator_selector(
               const SparseDiagonalProductType& expr, int outer)
@@ -141,12 +155,12 @@ class ei_sparse_diagonal_product_inner_iterator_selector
   : public SparseCwiseBinaryOp<
       ei_scalar_product_op<typename Rhs::Scalar>,
       SparseInnerVectorSet<Lhs,1>,
-      NestByValue<Transpose<typename Rhs::_CoeffsVectorType> > >::InnerIterator
+      NestByValue<Transpose<typename Rhs::DiagonalVectorType> > >::InnerIterator
 {
     typedef typename SparseCwiseBinaryOp<
       ei_scalar_product_op<typename Rhs::Scalar>,
       SparseInnerVectorSet<Lhs,1>,
-      NestByValue<Transpose<typename Rhs::_CoeffsVectorType> > >::InnerIterator Base;
+      NestByValue<Transpose<typename Rhs::DiagonalVectorType> > >::InnerIterator Base;
   public:
     inline ei_sparse_diagonal_product_inner_iterator_selector(
               const SparseDiagonalProductType& expr, int outer)
@@ -154,4 +168,14 @@ class ei_sparse_diagonal_product_inner_iterator_selector
     {}
 };
 
+// SparseMatrixBase functions
+
+template<typename Derived>
+template<typename OtherDerived>
+const SparseDiagonalProduct<Derived,OtherDerived>
+SparseMatrixBase<Derived>::operator*(const DiagonalBase<OtherDerived> &other) const
+{
+  return SparseDiagonalProduct<Derived,OtherDerived>(this->derived(), other.derived());
+}
+
 #endif // EIGEN_SPARSE_DIAGONAL_PRODUCT_H
diff --git a/Eigen/src/Sparse/SparseMatrix.h b/Eigen/src/Sparse/SparseMatrix.h
index 03bc60bb5..b751b5cb9 100644
--- a/Eigen/src/Sparse/SparseMatrix.h
+++ b/Eigen/src/Sparse/SparseMatrix.h
@@ -371,7 +371,7 @@ class SparseMatrix
       const int outerSize = IsRowMajor ? rows : cols;
       m_innerSize = IsRowMajor ? cols : rows;
       m_data.clear();
-      if (m_outerSize != outerSize)
+      if (m_outerSize != outerSize || m_outerSize==0)
       {
         delete[] m_outerIndex;
         m_outerIndex = new int [outerSize+1];
diff --git a/Eigen/src/Sparse/SparseMatrixBase.h b/Eigen/src/Sparse/SparseMatrixBase.h
index 8d9406cfb..8c53757b0 100644
--- a/Eigen/src/Sparse/SparseMatrixBase.h
+++ b/Eigen/src/Sparse/SparseMatrixBase.h
@@ -1,7 +1,7 @@
 // This file is part of Eigen, a lightweight C++ template library
 // for linear algebra.
 //
-// Copyright (C) 2008 Gael Guennebaud <g.gael@free.fr>
+// Copyright (C) 2008-2009 Gael Guennebaud <g.gael@free.fr>
 //
 // Eigen is free software; you can redistribute it and/or
 // modify it under the terms of the GNU Lesser General Public
@@ -37,6 +37,9 @@
   *
   */
 template<typename Derived> class SparseMatrixBase
+#ifndef EIGEN_PARSED_BY_DOXYGEN
+  : public MultiplierBase<Derived>
+#endif // not EIGEN_PARSED_BY_DOXYGEN
 {
   public:
 
@@ -301,10 +304,22 @@ template<typename Derived> class SparseMatrixBase
     { return matrix*scalar; }
 
 
+    // sparse * sparse
     template<typename OtherDerived>
     const typename SparseProductReturnType<Derived,OtherDerived>::Type
     operator*(const SparseMatrixBase<OtherDerived> &other) const;
 
+    // sparse * diagonal
+    template<typename OtherDerived>
+    const SparseDiagonalProduct<Derived,OtherDerived>
+    operator*(const DiagonalBase<OtherDerived> &other) const;
+
+    // diagonal * sparse
+    template<typename OtherDerived> friend 
+    const SparseDiagonalProduct<OtherDerived,Derived>
+    operator*(const DiagonalBase<OtherDerived> &lhs, const SparseMatrixBase& rhs)
+    { return SparseDiagonalProduct<OtherDerived,Derived>(lhs.derived(), rhs.derived()); }
+
     // dense * sparse (return a dense object)
     template<typename OtherDerived> friend
     const typename SparseProductReturnType<OtherDerived,Derived>::Type
diff --git a/Eigen/src/Sparse/SparseProduct.h b/Eigen/src/Sparse/SparseProduct.h
index 5222439d8..a9ad8aad4 100644
--- a/Eigen/src/Sparse/SparseProduct.h
+++ b/Eigen/src/Sparse/SparseProduct.h
@@ -1,7 +1,7 @@
 // This file is part of Eigen, a lightweight C++ template library
 // for linear algebra.
 //
-// Copyright (C) 2008 Gael Guennebaud <g.gael@free.fr>
+// Copyright (C) 2008-2009 Gael Guennebaud <g.gael@free.fr>
 //
 // Eigen is free software; you can redistribute it and/or
 // modify it under the terms of the GNU Lesser General Public
@@ -29,9 +29,7 @@ template<typename Lhs, typename Rhs> struct ei_sparse_product_mode
 {
   enum {
 
-    value = ei_is_diagonal<Lhs>::ret || ei_is_diagonal<Rhs>::ret
-          ? DiagonalProduct
-          :  (Rhs::Flags&Lhs::Flags&SparseBit)==SparseBit
+    value = (Rhs::Flags&Lhs::Flags&SparseBit)==SparseBit
           ? SparseTimeSparseProduct
           : (Lhs::Flags&SparseBit)==SparseBit
           ? SparseTimeDenseProduct
@@ -47,15 +45,6 @@ struct SparseProductReturnType
   typedef SparseProduct<LhsNested, RhsNested, ProductMode> Type;
 };
 
-template<typename Lhs, typename Rhs>
-struct SparseProductReturnType<Lhs,Rhs,DiagonalProduct>
-{
-  typedef const typename ei_nested<Lhs,Rhs::RowsAtCompileTime>::type LhsNested;
-  typedef const typename ei_nested<Rhs,Lhs::RowsAtCompileTime>::type RhsNested;
-
-  typedef SparseDiagonalProduct<LhsNested, RhsNested> Type;
-};
-
 // sparse product return type specialization
 template<typename Lhs, typename Rhs>
 struct SparseProductReturnType<Lhs,Rhs,SparseTimeSparseProduct>
@@ -106,7 +95,7 @@ struct ei_traits<SparseProduct<LhsNested, RhsNested, ProductMode> >
 //     RhsIsRowMajor = (RhsFlags & RowMajorBit)==RowMajorBit,
 
     EvalToRowMajor = (RhsFlags & LhsFlags & RowMajorBit),
-    ResultIsSparse = ProductMode==SparseTimeSparseProduct || ProductMode==DiagonalProduct,
+    ResultIsSparse = ProductMode==SparseTimeSparseProduct,
 
     RemovedBits = ~( (EvalToRowMajor ? 0 : RowMajorBit) | (ResultIsSparse ? 0 : SparseBit) ),