* optimize sum() for sparse matrices and vectors

* fix the row()/col() functions of some InnerVector

Eigen/src/Sparse/DynamicSparseMatrix.h

@@ -67,10 +67,10 @@ class DynamicSparseMatrix
     // EIGEN_SPARSE_INHERIT_ASSIGNMENT_OPERATOR(DynamicSparseMatrix, +=)
     // EIGEN_SPARSE_INHERIT_ASSIGNMENT_OPERATOR(DynamicSparseMatrix, -=)
     typedef MappedSparseMatrix<Scalar,Flags> Map;
+    using Base::IsRowMajor;
 
   protected:
 
-    enum { IsRowMajor = Base::IsRowMajor };
     typedef DynamicSparseMatrix<Scalar,(Flags&~RowMajorBit)|(IsRowMajor?RowMajorBit:0)> TransposedSparseMatrix;
 
     int m_innerSize;

Eigen/src/Sparse/SparseBlock.h

@@ -43,16 +43,21 @@ template<typename MatrixType, int Size>
 class SparseInnerVectorSet : ei_no_assignment_operator,
   public SparseMatrixBase<SparseInnerVectorSet<MatrixType, Size> >
 {
-    enum { IsRowMajor = ei_traits<SparseInnerVectorSet>::IsRowMajor };
   public:
 
+    enum { IsRowMajor = ei_traits<SparseInnerVectorSet>::IsRowMajor };
+
     EIGEN_SPARSE_GENERIC_PUBLIC_INTERFACE(SparseInnerVectorSet)
     class InnerIterator: public MatrixType::InnerIterator
     {
       public:
         inline InnerIterator(const SparseInnerVectorSet& xpr, int outer)
-          : MatrixType::InnerIterator(xpr.m_matrix, xpr.m_outerStart + outer)
+          : MatrixType::InnerIterator(xpr.m_matrix, xpr.m_outerStart + outer), m_outer(outer)
         {}
+        inline int row() const { return IsRowMajor ? m_outer : this->index(); }
+        inline int col() const { return IsRowMajor ? this->index() : m_outer; }
+      protected:
+        int m_outer;
     };
 
     inline SparseInnerVectorSet(const MatrixType& matrix, int outerStart, int outerSize)
@@ -100,16 +105,21 @@ class SparseInnerVectorSet<DynamicSparseMatrix<_Scalar, _Options>, Size>
   : public SparseMatrixBase<SparseInnerVectorSet<DynamicSparseMatrix<_Scalar, _Options>, Size> >
 {
     typedef DynamicSparseMatrix<_Scalar, _Options> MatrixType;
-    enum { IsRowMajor = ei_traits<SparseInnerVectorSet>::IsRowMajor };
   public:
 
+    enum { IsRowMajor = ei_traits<SparseInnerVectorSet>::IsRowMajor };
+
     EIGEN_SPARSE_GENERIC_PUBLIC_INTERFACE(SparseInnerVectorSet)
     class InnerIterator: public MatrixType::InnerIterator
     {
       public:
         inline InnerIterator(const SparseInnerVectorSet& xpr, int outer)
-          : MatrixType::InnerIterator(xpr.m_matrix, xpr.m_outerStart + outer)
+          : MatrixType::InnerIterator(xpr.m_matrix, xpr.m_outerStart + outer), m_outer(outer)
         {}
+        inline int row() const { return IsRowMajor ? m_outer : this->index(); }
+        inline int col() const { return IsRowMajor ? this->index() : m_outer; }
+      protected:
+        int m_outer;
     };
 
     inline SparseInnerVectorSet(const MatrixType& matrix, int outerStart, int outerSize)
@@ -193,16 +203,21 @@ class SparseInnerVectorSet<SparseMatrix<_Scalar, _Options>, Size>
   : public SparseMatrixBase<SparseInnerVectorSet<SparseMatrix<_Scalar, _Options>, Size> >
 {
     typedef SparseMatrix<_Scalar, _Options> MatrixType;
-    enum { IsRowMajor = ei_traits<SparseInnerVectorSet>::IsRowMajor };
   public:
 
+    enum { IsRowMajor = ei_traits<SparseInnerVectorSet>::IsRowMajor };
+
     EIGEN_SPARSE_GENERIC_PUBLIC_INTERFACE(SparseInnerVectorSet)
     class InnerIterator: public MatrixType::InnerIterator
     {
       public:
         inline InnerIterator(const SparseInnerVectorSet& xpr, int outer)
-          : MatrixType::InnerIterator(xpr.m_matrix, xpr.m_outerStart + outer)
+          : MatrixType::InnerIterator(xpr.m_matrix, xpr.m_outerStart + outer), m_outer(outer)
         {}
+        inline int row() const { return IsRowMajor ? m_outer : this->index(); }
+        inline int col() const { return IsRowMajor ? this->index() : m_outer; }
+      protected:
+        int m_outer;
     };
 
     inline SparseInnerVectorSet(const MatrixType& matrix, int outerStart, int outerSize)

Eigen/src/Sparse/SparseCwiseBinaryOp.h

@@ -186,8 +186,8 @@ class ei_sparse_cwise_binary_op_inner_iterator_selector<BinaryOp, Lhs, Rhs, Deri
     EIGEN_STRONG_INLINE Scalar value() const { return m_value; }
 
     EIGEN_STRONG_INLINE int index() const { return m_id; }
-    EIGEN_STRONG_INLINE int row() const { return m_lhsIter.row(); }
+    EIGEN_STRONG_INLINE int row() const { return Lhs::IsRowMajor ? m_lhsIter.row() : index(); }
-    EIGEN_STRONG_INLINE int col() const { return m_lhsIter.col(); }
+    EIGEN_STRONG_INLINE int col() const { return Lhs::IsRowMajor ? index() : m_lhsIter.col(); }
 
     EIGEN_STRONG_INLINE operator bool() const { return m_id>=0; }
 

Eigen/src/Sparse/SparseMatrix.h

@@ -25,17 +25,24 @@
 #ifndef EIGEN_SPARSEMATRIX_H
 #define EIGEN_SPARSEMATRIX_H
 
-/** \class SparseMatrix
+/** \ingroup Sparse_Module
   *
-  * \brief Sparse matrix
+  * \class SparseMatrix
+  *
+  * \brief The main sparse matrix class
+  *
+  * This class implements a sparse matrix using the very common compressed row/column storage
+  * scheme.
   *
   * \param _Scalar the scalar type, i.e. the type of the coefficients
+  * \param _Options Union of bit flags controlling the storage scheme. Currently the only possibility
+  *                 is RowMajor. The default is 0 which means column-major.
   *
   * See http://www.netlib.org/linalg/html_templates/node91.html for details on the storage scheme.
   *
   */
-template<typename _Scalar, int _Flags>
+template<typename _Scalar, int _Options>
-struct ei_traits<SparseMatrix<_Scalar, _Flags> >
+struct ei_traits<SparseMatrix<_Scalar, _Options> >
 {
   typedef _Scalar Scalar;
   enum {
@@ -43,17 +50,15 @@ struct ei_traits<SparseMatrix<_Scalar, _Flags> >
     ColsAtCompileTime = Dynamic,
     MaxRowsAtCompileTime = Dynamic,
     MaxColsAtCompileTime = Dynamic,
-    Flags = SparseBit | _Flags,
+    Flags = SparseBit | _Options,
     CoeffReadCost = NumTraits<Scalar>::ReadCost,
     SupportedAccessPatterns = InnerRandomAccessPattern
   };
 };
 
-
+template<typename _Scalar, int _Options>
-
-template<typename _Scalar, int _Flags>
 class SparseMatrix
-  : public SparseMatrixBase<SparseMatrix<_Scalar, _Flags> >
+  : public SparseMatrixBase<SparseMatrix<_Scalar, _Options> >
 {
   public:
     EIGEN_SPARSE_GENERIC_PUBLIC_INTERFACE(SparseMatrix)
@@ -64,10 +69,10 @@ class SparseMatrix
     // EIGEN_SPARSE_INHERIT_SCALAR_ASSIGNMENT_OPERATOR(SparseMatrix, /=)
 
     typedef MappedSparseMatrix<Scalar,Flags> Map;
+    using Base::IsRowMajor;
 
   protected:
 
-    enum { IsRowMajor = Base::IsRowMajor };
     typedef SparseMatrix<Scalar,(Flags&~RowMajorBit)|(IsRowMajor?RowMajorBit:0)> TransposedSparseMatrix;
 
     int m_outerSize;
@@ -508,10 +513,13 @@ class SparseMatrix
     {
       delete[] m_outerIndex;
     }
+
+    /** Overloaded for performance */
+    Scalar sum() const;
 };
 
-template<typename Scalar, int _Flags>
+template<typename Scalar, int _Options>
-class SparseMatrix<Scalar,_Flags>::InnerIterator
+class SparseMatrix<Scalar,_Options>::InnerIterator
 {
   public:
     InnerIterator(const SparseMatrix& mat, int outer)

Eigen/src/Sparse/SparseMatrixBase.h

@@ -25,6 +25,17 @@
 #ifndef EIGEN_SPARSEMATRIXBASE_H
 #define EIGEN_SPARSEMATRIXBASE_H
 
+/** \ingroup Sparse_Module
+  *
+  * \class SparseMatrixBase
+  *
+  * \brief Base class of any sparse matrices or sparse expressions
+  *
+  * \param Derived
+  *
+  *
+  *
+  */
 template<typename Derived> class SparseMatrixBase
 {
   public:
@@ -432,10 +443,16 @@ template<typename Derived> class SparseMatrixBase
         for (int j=0; j<outerSize(); ++j)
         {
           for (typename Derived::InnerIterator i(derived(),j); i; ++i)
+          {
             if(IsRowMajor)
-              res.coeffRef(j,i.index()) = i.value();
+              std::cerr << i.row() << "," << i.col() << " == " << j << "," << i.index() << "\n";
             else
-              res.coeffRef(i.index(),j) = i.value();
+              std::cerr << i.row() << "," << i.col() << " == " << i.index() << "," << j << "\n";
+//             if(IsRowMajor)
+              res.coeffRef(i.row(),i.col()) = i.value();
+//             else
+//               res.coeffRef(i.index(),j) = i.value();
+          }
         }
         return res;
       }

Eigen/src/Sparse/SparseRedux.h

@@ -37,4 +37,20 @@ SparseMatrixBase<Derived>::sum() const
   return res;
 }
 
+template<typename _Scalar, int _Options>
+typename ei_traits<SparseMatrix<_Scalar,_Options> >::Scalar
+SparseMatrix<_Scalar,_Options>::sum() const
+{
+  ei_assert(rows()>0 && cols()>0 && "you are using a non initialized matrix");
+  return Matrix<Scalar,1,Dynamic>::Map(m_data.value(0), m_data.size()).sum();
+}
+
+template<typename _Scalar, int _Options>
+typename ei_traits<SparseVector<_Scalar,_Options> >::Scalar
+SparseVector<_Scalar,_Options>::sum() const
+{
+  ei_assert(rows()>0 && cols()>0 && "you are using a non initialized matrix");
+  return Matrix<Scalar,1,Dynamic>::Map(m_data.value(0), m_data.size()).sum();
+}
+
 #endif // EIGEN_SPARSEREDUX_H

Eigen/src/Sparse/SparseTranspose.h

@@ -66,20 +66,26 @@ template<typename MatrixType> class SparseTranspose
 
 template<typename MatrixType> class SparseTranspose<MatrixType>::InnerIterator : public MatrixType::InnerIterator
 {
+    typedef typename MatrixType::InnerIterator Base;
   public:
 
     EIGEN_STRONG_INLINE InnerIterator(const SparseTranspose& trans, int outer)
-      : MatrixType::InnerIterator(trans.m_matrix, outer)
+      : Base(trans.m_matrix, outer)
     {}
+    inline int row() const { return Base::col(); }
+    inline int col() const { return Base::row(); }
 };
 
 template<typename MatrixType> class SparseTranspose<MatrixType>::ReverseInnerIterator : public MatrixType::ReverseInnerIterator
 {
+    typedef typename MatrixType::ReverseInnerIterator Base;
   public:
 
     EIGEN_STRONG_INLINE ReverseInnerIterator(const SparseTranspose& xpr, int outer)
-      : MatrixType::ReverseInnerIterator(xpr.m_matrix, outer)
+      : Base(xpr.m_matrix, outer)
     {}
+    inline int row() const { return Base::col(); }
+    inline int col() const { return Base::row(); }
 };
 
 #endif // EIGEN_SPARSETRANSPOSE_H

Eigen/src/Sparse/SparseVector.h

@@ -34,26 +34,26 @@
   * See http://www.netlib.org/linalg/html_templates/node91.html for details on the storage scheme.
   *
   */
-template<typename _Scalar, int _Flags>
+template<typename _Scalar, int _Options>
-struct ei_traits<SparseVector<_Scalar, _Flags> >
+struct ei_traits<SparseVector<_Scalar, _Options> >
 {
   typedef _Scalar Scalar;
   enum {
-    IsColVector = _Flags & RowMajorBit ? 0 : 1,
+    IsColVector = _Options & RowMajorBit ? 0 : 1,
 
     RowsAtCompileTime = IsColVector ? Dynamic : 1,
     ColsAtCompileTime = IsColVector ? 1 : Dynamic,
     MaxRowsAtCompileTime = RowsAtCompileTime,
     MaxColsAtCompileTime = ColsAtCompileTime,
-    Flags = SparseBit | _Flags,
+    Flags = SparseBit | _Options,
     CoeffReadCost = NumTraits<Scalar>::ReadCost,
     SupportedAccessPatterns = InnerRandomAccessPattern
   };
 };
 
-template<typename _Scalar, int _Flags>
+template<typename _Scalar, int _Options>
 class SparseVector
-  : public SparseMatrixBase<SparseVector<_Scalar, _Flags> >
+  : public SparseMatrixBase<SparseVector<_Scalar, _Options> >
 {
   public:
     EIGEN_SPARSE_GENERIC_PUBLIC_INTERFACE(SparseVector)
@@ -357,10 +357,13 @@ class SparseVector
 
     /** Destructor */
     inline ~SparseVector() {}
+
+    /** Overloaded for performance */
+    Scalar sum() const;
 };
 
-template<typename Scalar, int _Flags>
+template<typename Scalar, int _Options>
-class SparseVector<Scalar,_Flags>::InnerIterator
+class SparseVector<Scalar,_Options>::InnerIterator
 {
   public:
     InnerIterator(const SparseVector& vec, int outer=0)