Big refactoring/cleaning in the spasre module with

in particular the addition of a selfadjointView, and the
extension of triangularView. The rest is cleaning and does not
change/extend the API.

18 files changed, 448 insertions, 271 deletions

diff --git a/Eigen/Sparse b/Eigen/Sparse
index 552dbde4a..f4dcad07e 100644
--- a/Eigen/Sparse
+++ b/Eigen/Sparse
@@ -98,7 +98,6 @@ struct Sparse {};
 #include "src/Sparse/SparseVector.h"
 #include "src/Sparse/CoreIterators.h"
 #include "src/Sparse/SparseTranspose.h"
-#include "src/Sparse/SparseCwise.h"
 #include "src/Sparse/SparseCwiseUnaryOp.h"
 #include "src/Sparse/SparseCwiseBinaryOp.h"
 #include "src/Sparse/SparseDot.h"
@@ -107,12 +106,12 @@ struct Sparse {};
 #include "src/Sparse/SparseFuzzy.h"
 #include "src/Sparse/SparseProduct.h"
 #include "src/Sparse/SparseDiagonalProduct.h"
-#include "src/Sparse/SparseTriangular.h"
+#include "src/Sparse/SparseTriangularView.h"
+#include "src/Sparse/SparseSelfAdjointView.h"
 #include "src/Sparse/TriangularSolver.h"
 #include "src/Sparse/SparseLLT.h"
 #include "src/Sparse/SparseLDLT.h"
 #include "src/Sparse/SparseLU.h"
-#include "src/Sparse/SparseExpressionMaker.h"
 
 #ifdef EIGEN_CHOLMOD_SUPPORT
 # include "src/Sparse/CholmodSupport.h"
diff --git a/Eigen/src/Core/MatrixBase.h b/Eigen/src/Core/MatrixBase.h
index b752c7821..67e270af7 100644
--- a/Eigen/src/Core/MatrixBase.h
+++ b/Eigen/src/Core/MatrixBase.h
@@ -710,15 +710,6 @@ template<typename Derived> class MatrixBase
     typedef Homogeneous<Derived,MatrixBase<Derived>::ColsAtCompileTime==1?Vertical:Horizontal> HomogeneousReturnType;
     const HomogeneousReturnType homogeneous() const;
 
-/////////// Sparse module ///////////
-
-    // dense = sparse * dense
-    template<typename Derived1, typename Derived2>
-    Derived& lazyAssign(const SparseProduct<Derived1,Derived2,SparseTimeDenseProduct>& product);
-    // dense = dense * sparse
-    template<typename Derived1, typename Derived2>
-    Derived& lazyAssign(const SparseProduct<Derived1,Derived2,DenseTimeSparseProduct>& product);
-
 ////////// Householder module ///////////
 
     void makeHouseholderInPlace(Scalar& tau, RealScalar& beta);
diff --git a/Eigen/src/Core/SelfAdjointView.h b/Eigen/src/Core/SelfAdjointView.h
index 95ce666c9..3e435853c 100644
--- a/Eigen/src/Core/SelfAdjointView.h
+++ b/Eigen/src/Core/SelfAdjointView.h
@@ -201,15 +201,15 @@ struct ei_triangular_assignment_selector<Derived1, Derived2, SelfAdjoint, Dynami
 ***************************************************************************/
 
 template<typename Derived>
-template<unsigned int Mode>
-const SelfAdjointView<Derived, Mode> MatrixBase<Derived>::selfadjointView() const
+template<unsigned int UpLo>
+const SelfAdjointView<Derived, UpLo> MatrixBase<Derived>::selfadjointView() const
 {
   return derived();
 }
 
 template<typename Derived>
-template<unsigned int Mode>
-SelfAdjointView<Derived, Mode> MatrixBase<Derived>::selfadjointView()
+template<unsigned int UpLo>
+SelfAdjointView<Derived, UpLo> MatrixBase<Derived>::selfadjointView()
 {
   return derived();
 }
diff --git a/Eigen/src/Core/util/Constants.h b/Eigen/src/Core/util/Constants.h
index ed1cfa35c..c801d8049 100644
--- a/Eigen/src/Core/util/Constants.h
+++ b/Eigen/src/Core/util/Constants.h
@@ -193,7 +193,7 @@ enum { AsRequested=0, EnforceAlignedAccess=2 };
 enum { ConditionalJumpCost = 5 };
 enum CornerType { TopLeft, TopRight, BottomLeft, BottomRight };
 enum DirectionType { Vertical, Horizontal, BothDirections };
-enum ProductEvaluationMode { NormalProduct, CacheFriendlyProduct, SparseTimeSparseProduct, SparseTimeDenseProduct, DenseTimeSparseProduct };
+enum ProductEvaluationMode { NormalProduct, CacheFriendlyProduct };
 
 enum {
   /** \internal Equivalent to a slice vectorization for fixed-size matrices having good alignment
diff --git a/Eigen/src/Core/util/ForwardDeclarations.h b/Eigen/src/Core/util/ForwardDeclarations.h
index 5591c754d..2fedbbc07 100644
--- a/Eigen/src/Core/util/ForwardDeclarations.h
+++ b/Eigen/src/Core/util/ForwardDeclarations.h
@@ -146,8 +146,4 @@ template<typename Scalar,int Dim> class Translation;
 template<typename Scalar> class UniformScaling;
 template<typename MatrixType,int Direction> class Homogeneous;
 
-// Sparse module:
-template<typename Lhs, typename Rhs, int ProductMode> class SparseProduct;
-
-
 #endif // EIGEN_FORWARDDECLARATIONS_H
diff --git a/Eigen/src/Sparse/SparseCwiseBinaryOp.h b/Eigen/src/Sparse/SparseCwiseBinaryOp.h
index 9c354b4d3..5dd7edc00 100644
--- a/Eigen/src/Sparse/SparseCwiseBinaryOp.h
+++ b/Eigen/src/Sparse/SparseCwiseBinaryOp.h
@@ -42,35 +42,6 @@
 //  4 - dense op dense     product      dense
 //                         generic      dense
 
-// template<typename BinaryOp, typename Lhs, typename Rhs>
-// struct ei_traits<SparseCwiseBinaryOp<BinaryOp, Lhs, Rhs> >
-// {
-//   typedef typename ei_result_of<
-//                      BinaryOp(
-//                        typename Lhs::Scalar,
-//                        typename Rhs::Scalar
-//                      )
-//                    >::type Scalar;
-//   typedef typename ei_promote_storage_type<typename ei_traits<Lhs>::StorageType,
-//                                            typename ei_traits<Rhs>::StorageType>::ret StorageType;
-//   typedef typename Lhs::Nested LhsNested;
-//   typedef typename Rhs::Nested RhsNested;
-//   typedef typename ei_unref<LhsNested>::type _LhsNested;
-//   typedef typename ei_unref<RhsNested>::type _RhsNested;
-//   enum {
-//     LhsCoeffReadCost = _LhsNested::CoeffReadCost,
-//     RhsCoeffReadCost = _RhsNested::CoeffReadCost,
-//     LhsFlags = _LhsNested::Flags,
-//     RhsFlags = _RhsNested::Flags,
-//     RowsAtCompileTime = Lhs::RowsAtCompileTime,
-//     ColsAtCompileTime = Lhs::ColsAtCompileTime,
-//     MaxRowsAtCompileTime = Lhs::MaxRowsAtCompileTime,
-//     MaxColsAtCompileTime = Lhs::MaxColsAtCompileTime,
-//     Flags = (int(LhsFlags) | int(RhsFlags)) & HereditaryBits,
-//     CoeffReadCost = LhsCoeffReadCost + RhsCoeffReadCost + ei_functor_traits<BinaryOp>::Cost
-//   };
-// };
-
 template<> struct ei_promote_storage_type<Dense,Sparse>
 { typedef Sparse ret; };
 
@@ -82,16 +53,9 @@ class CwiseBinaryOpImpl<BinaryOp, Lhs, Rhs, Sparse>
   : public SparseMatrixBase<CwiseBinaryOp<BinaryOp, Lhs, Rhs> >
 {
   public:
-
     class InnerIterator;
-
     typedef CwiseBinaryOp<BinaryOp, Lhs, Rhs> Derived;
     EIGEN_SPARSE_PUBLIC_INTERFACE(Derived)
-//     typedef typename ei_traits<SparseCwiseBinaryOp>::LhsNested LhsNested;
-//     typedef typename ei_traits<SparseCwiseBinaryOp>::RhsNested RhsNested;
-//     typedef typename ei_unref<LhsNested>::type _LhsNested;
-//     typedef typename ei_unref<RhsNested>::type _RhsNested;
-
 };
 
 template<typename BinaryOp, typename Lhs, typename Rhs, typename Derived,
diff --git a/Eigen/src/Sparse/SparseLDLT.h b/Eigen/src/Sparse/SparseLDLT.h
index 631bff04b..ea5d345f5 100644
--- a/Eigen/src/Sparse/SparseLDLT.h
+++ b/Eigen/src/Sparse/SparseLDLT.h
@@ -333,12 +333,12 @@ bool SparseLDLT<MatrixType, Backend>::solveInPlace(MatrixBase<Derived> &b) const
     return false;
 
   if (m_matrix.nonZeros()>0) // otherwise L==I
-    m_matrix.template triangular<UnitLowerTriangular>().solveInPlace(b);
+    m_matrix.template triangularView<UnitLowerTriangular>().solveInPlace(b);
   b = b.cwise() / m_diag;
   // FIXME should be .adjoint() but it fails to compile...
 
   if (m_matrix.nonZeros()>0) // otherwise L==I
-    m_matrix.transpose().template triangular<UnitUpperTriangular>().solveInPlace(b);
+    m_matrix.transpose().template triangularView<UnitUpperTriangular>().solveInPlace(b);
 
   return true;
 }
diff --git a/Eigen/src/Sparse/SparseLLT.h b/Eigen/src/Sparse/SparseLLT.h
index b2f65b944..6307b2493 100644
--- a/Eigen/src/Sparse/SparseLLT.h
+++ b/Eigen/src/Sparse/SparseLLT.h
@@ -193,15 +193,15 @@ bool SparseLLT<MatrixType, Backend>::solveInPlace(MatrixBase<Derived> &b) const
   const int size = m_matrix.rows();
   ei_assert(size==b.rows());
 
-  m_matrix.template triangular<LowerTriangular>().solveInPlace(b);
+  m_matrix.template triangularView<LowerTriangular>().solveInPlace(b);
   // FIXME should be simply .adjoint() but it fails to compile...
   if (NumTraits<Scalar>::IsComplex)
   {
     CholMatrixType aux = m_matrix.conjugate();
-    aux.transpose().template triangular<UpperTriangular>().solveInPlace(b);
+    aux.transpose().template triangularView<UpperTriangular>().solveInPlace(b);
   }
   else
-    m_matrix.transpose().template triangular<UpperTriangular>().solveInPlace(b);
+    m_matrix.transpose().template triangularView<UpperTriangular>().solveInPlace(b);
 
   return true;
 }
diff --git a/Eigen/src/Sparse/SparseMatrix.h b/Eigen/src/Sparse/SparseMatrix.h
index f3915af70..5e89d3f53 100644
--- a/Eigen/src/Sparse/SparseMatrix.h
+++ b/Eigen/src/Sparse/SparseMatrix.h
@@ -540,6 +540,7 @@ class SparseMatrix<Scalar,_Options>::InnerIterator
     inline Scalar& valueRef() { return const_cast<Scalar&>(m_matrix.m_data.value(m_id)); }
 
     inline int index() const { return m_matrix.m_data.index(m_id); }
+    inline int outer() const { return m_outer; }
     inline int row() const { return IsRowMajor ? m_outer : index(); }
     inline int col() const { return IsRowMajor ? index() : m_outer; }
 
diff --git a/Eigen/src/Sparse/SparseMatrixBase.h b/Eigen/src/Sparse/SparseMatrixBase.h
index 2eda425a7..2b6970818 100644
--- a/Eigen/src/Sparse/SparseMatrixBase.h
+++ b/Eigen/src/Sparse/SparseMatrixBase.h
@@ -244,7 +244,7 @@ template<typename Derived> class SparseMatrixBase : public AnyMatrixBase<Derived
     }
 
     template<typename Lhs, typename Rhs>
-    inline Derived& operator=(const SparseProduct<Lhs,Rhs,SparseTimeSparseProduct>& product);
+    inline Derived& operator=(const SparseProduct<Lhs,Rhs>& product);
 
     friend std::ostream & operator << (std::ostream & s, const SparseMatrixBase& m)
     {
@@ -337,12 +337,13 @@ template<typename Derived> class SparseMatrixBase : public AnyMatrixBase<Derived
 
     // dense * sparse (return a dense object)
     template<typename OtherDerived> friend
-    const typename SparseProductReturnType<OtherDerived,Derived>::Type
+    const DenseTimeSparseProduct<OtherDerived,Derived>
     operator*(const MatrixBase<OtherDerived>& lhs, const Derived& rhs)
-    { return typename SparseProductReturnType<OtherDerived,Derived>::Type(lhs.derived(),rhs); }
+    { return DenseTimeSparseProduct<OtherDerived,Derived>(lhs.derived(),rhs); }
 
+    // sparse * dense (returns a dense object)
     template<typename OtherDerived>
-    const typename SparseProductReturnType<Derived,OtherDerived>::Type
+    const SparseTimeDenseProduct<Derived,OtherDerived>
     operator*(const MatrixBase<OtherDerived> &other) const;
 
     template<typename OtherDerived>
@@ -360,7 +361,10 @@ template<typename Derived> class SparseMatrixBase : public AnyMatrixBase<Derived
 //     void solveTriangularInPlace(SparseMatrixBase<OtherDerived>& other) const;
 
     template<int Mode>
-    inline const SparseTriangular<Derived, Mode> triangular() const;
+    inline const SparseTriangularView<Derived, Mode> triangularView() const;
+
+    template<unsigned int UpLo> inline const SparseSelfAdjointView<Derived, UpLo> selfadjointView() const;
+    template<unsigned int UpLo> inline SparseSelfAdjointView<Derived, UpLo> selfadjointView();
 
     template<typename OtherDerived> Scalar dot(const MatrixBase<OtherDerived>& other) const;
     template<typename OtherDerived> Scalar dot(const SparseMatrixBase<OtherDerived>& other) const;
diff --git a/Eigen/src/Sparse/SparseProduct.h b/Eigen/src/Sparse/SparseProduct.h
index 16f36d640..2aaa8713c 100644
--- a/Eigen/src/Sparse/SparseProduct.h
+++ b/Eigen/src/Sparse/SparseProduct.h
@@ -25,22 +25,8 @@
 #ifndef EIGEN_SPARSEPRODUCT_H
 #define EIGEN_SPARSEPRODUCT_H
 
-template<typename Lhs, typename Rhs> struct ei_sparse_product_mode<Lhs,Rhs,Sparse,Sparse> { enum { value = SparseTimeSparseProduct }; };
-template<typename Lhs, typename Rhs> struct ei_sparse_product_mode<Lhs,Rhs,Dense, Sparse> { enum { value = DenseTimeSparseProduct }; };
-template<typename Lhs, typename Rhs> struct ei_sparse_product_mode<Lhs,Rhs,Sparse,Dense>  { enum { value = SparseTimeDenseProduct }; };
-
-template<typename Lhs, typename Rhs, int ProductMode>
-struct SparseProductReturnType
-{
-  typedef const typename ei_nested<Lhs,Rhs::RowsAtCompileTime>::type LhsNested;
-  typedef const typename ei_nested<Rhs,Lhs::RowsAtCompileTime>::type RhsNested;
-
-  typedef SparseProduct<LhsNested, RhsNested, ProductMode> Type;
-};
-
-// sparse product return type specialization
 template<typename Lhs, typename Rhs>
-struct SparseProductReturnType<Lhs,Rhs,SparseTimeSparseProduct>
+struct SparseProductReturnType
 {
   typedef typename ei_traits<Lhs>::Scalar Scalar;
   enum {
@@ -60,11 +46,11 @@ struct SparseProductReturnType<Lhs,Rhs,SparseTimeSparseProduct>
     SparseMatrix<Scalar,0>,
     const typename ei_nested<Rhs,Lhs::RowsAtCompileTime>::type>::ret RhsNested;
 
-  typedef SparseProduct<LhsNested, RhsNested, SparseTimeSparseProduct> Type;
+  typedef SparseProduct<LhsNested, RhsNested> Type;
 };
 
-template<typename LhsNested, typename RhsNested, int ProductMode>
-struct ei_traits<SparseProduct<LhsNested, RhsNested, ProductMode> >
+template<typename LhsNested, typename RhsNested>
+struct ei_traits<SparseProduct<LhsNested, RhsNested> >
 {
   // clean the nested types:
   typedef typename ei_cleantype<LhsNested>::type _LhsNested;
@@ -85,7 +71,6 @@ struct ei_traits<SparseProduct<LhsNested, RhsNested, ProductMode> >
     MaxColsAtCompileTime = _RhsNested::MaxColsAtCompileTime,
 
     EvalToRowMajor = (RhsFlags & LhsFlags & RowMajorBit),
-    ResultIsSparse = ProductMode==SparseTimeSparseProduct,
 
     RemovedBits = ~(EvalToRowMajor ? 0 : RowMajorBit),
 
@@ -96,16 +81,14 @@ struct ei_traits<SparseProduct<LhsNested, RhsNested, ProductMode> >
     CoeffReadCost = Dynamic
   };
 
-  typedef typename ei_meta_if<ResultIsSparse, Sparse, Dense>::ret StorageType;
+  typedef Sparse StorageType;
 
-  typedef typename ei_meta_if<ResultIsSparse,
-    SparseMatrixBase<SparseProduct<LhsNested, RhsNested, ProductMode> >,
-    MatrixBase<SparseProduct<LhsNested, RhsNested, ProductMode> > >::ret Base;
+  typedef SparseMatrixBase<SparseProduct<LhsNested, RhsNested> > Base;
 };
 
-template<typename LhsNested, typename RhsNested, int ProductMode>
+template<typename LhsNested, typename RhsNested>
 class SparseProduct : ei_no_assignment_operator,
-  public ei_traits<SparseProduct<LhsNested, RhsNested, ProductMode> >::Base
+  public ei_traits<SparseProduct<LhsNested, RhsNested> >::Base
 {
   public:
 
@@ -256,38 +239,14 @@ struct ei_sparse_product_selector<Lhs,Rhs,ResultType,RowMajor,RowMajor,ColMajor>
   }
 };
 
-// NOTE eventually let's transpose one argument even in this case since it might be expensive if
-// the result is not dense.
-// template<typename Lhs, typename Rhs, typename ResultType, int ResStorageOrder>
-// struct ei_sparse_product_selector<Lhs,Rhs,ResultType,RowMajor,ColMajor,ResStorageOrder>
-// {
-//   static void run(const Lhs& lhs, const Rhs& rhs, ResultType& res)
-//   {
-//     // trivial product as lhs.row/rhs.col dot products
-//     // loop over the preferred order of the result
-//   }
-// };
-
 // NOTE the 2 others cases (col row *) must never occurs since they are caught
 // by ProductReturnType which transform it to (col col *) by evaluating rhs.
 
 
-// template<typename Derived>
-// template<typename Lhs, typename Rhs>
-// inline Derived& SparseMatrixBase<Derived>::lazyAssign(const SparseProduct<Lhs,Rhs>& product)
-// {
-// //   std::cout << "sparse product to dense\n";
-//   ei_sparse_product_selector<
-//     typename ei_cleantype<Lhs>::type,
-//     typename ei_cleantype<Rhs>::type,
-//     typename ei_cleantype<Derived>::type>::run(product.lhs(),product.rhs(),derived());
-//   return derived();
-// }
-
 // sparse = sparse * sparse
 template<typename Derived>
 template<typename Lhs, typename Rhs>
-inline Derived& SparseMatrixBase<Derived>::operator=(const SparseProduct<Lhs,Rhs,SparseTimeSparseProduct>& product)
+inline Derived& SparseMatrixBase<Derived>::operator=(const SparseProduct<Lhs,Rhs>& product)
 {
   ei_sparse_product_selector<
     typename ei_cleantype<Lhs>::type,
@@ -296,78 +255,79 @@ inline Derived& SparseMatrixBase<Derived>::operator=(const SparseProduct<Lhs,Rhs
   return derived();
 }
 
-// dense = sparse * dense
-// Note that here we force no inlining and separate the setZero() because GCC messes up otherwise
-template<typename Lhs, typename Rhs, typename Dest>
-EIGEN_DONT_INLINE void ei_sparse_time_dense_product(const Lhs& lhs, const Rhs& rhs, Dest& dst)
+
+template<typename Lhs, typename Rhs>
+struct ei_traits<SparseTimeDenseProduct<Lhs,Rhs> >
+ : ei_traits<ProductBase<SparseTimeDenseProduct<Lhs,Rhs>, Lhs, Rhs> >
 {
-  typedef typename ei_cleantype<Lhs>::type _Lhs;
-  typedef typename ei_cleantype<Rhs>::type _Rhs;
-  typedef typename _Lhs::InnerIterator LhsInnerIterator;
-  enum {
-    LhsIsRowMajor = (_Lhs::Flags&RowMajorBit)==RowMajorBit,
-    LhsIsSelfAdjoint = (_Lhs::Flags&SelfAdjointBit)==SelfAdjointBit,
-    ProcessFirstHalf = LhsIsSelfAdjoint
-      && (   ((_Lhs::Flags&(UpperTriangularBit|LowerTriangularBit))==0)
-          || ( (_Lhs::Flags&UpperTriangularBit) && !LhsIsRowMajor)
-          || ( (_Lhs::Flags&LowerTriangularBit) && LhsIsRowMajor) ),
-    ProcessSecondHalf = LhsIsSelfAdjoint && (!ProcessFirstHalf)
-  };
-  for (int j=0; j<lhs.outerSize(); ++j)
-  {
-    LhsInnerIterator i(lhs,j);
-    if (ProcessSecondHalf && i && (i.index()==j))
-    {
-      dst.row(j) += i.value() * rhs.row(j);
-      ++i;
-    }
-    typename Rhs::Scalar rhs_j = rhs.coeff(j,0);
-    Block<Dest,1,Dest::ColsAtCompileTime> dst_j(dst.row(LhsIsRowMajor ? j : 0));
-    for(; (ProcessFirstHalf ? i && i.index() < j : i) ; ++i)
+  typedef Dense StorageType;
+};
+
+template<typename Lhs, typename Rhs>
+class SparseTimeDenseProduct
+  : public ProductBase<SparseTimeDenseProduct<Lhs,Rhs>, Lhs, Rhs>
+{
+  public:
+    EIGEN_PRODUCT_PUBLIC_INTERFACE(SparseTimeDenseProduct)
+
+    SparseTimeDenseProduct(const Lhs& lhs, const Rhs& rhs) : Base(lhs,rhs)
+    {}
+
+    template<typename Dest> void scaleAndAddTo(Dest& dest, Scalar alpha) const
     {
-      if(LhsIsSelfAdjoint)
+      typedef typename ei_cleantype<Lhs>::type _Lhs;
+      typedef typename ei_cleantype<Rhs>::type _Rhs;
+      typedef typename _Lhs::InnerIterator LhsInnerIterator;
+      enum { LhsIsRowMajor = (_Lhs::Flags&RowMajorBit)==RowMajorBit };
+      for(int j=0; j<m_lhs.outerSize(); ++j)
       {
-        int a = LhsIsRowMajor ? j : i.index();
-        int b = LhsIsRowMajor ? i.index() : j;
-        typename Lhs::Scalar v = i.value();
-        dst.row(a) += (v) * rhs.row(b);
-        dst.row(b) += ei_conj(v) * rhs.row(a);
+        typename Rhs::Scalar rhs_j = alpha * m_rhs.coeff(j,0);
+        Block<Dest,1,Dest::ColsAtCompileTime> dest_j(dest.row(LhsIsRowMajor ? j : 0));
+        for(LhsInnerIterator it(m_lhs,j); it ;++it)
+        {
+          if(LhsIsRowMajor)                   dest_j += (alpha*it.value()) * m_rhs.row(it.index());
+          else if(Rhs::ColsAtCompileTime==1)  dest.coeffRef(it.index()) += it.value() * rhs_j;
+          else                                dest.row(it.index()) += (alpha*it.value()) * m_rhs.row(j);
+        }
       }
-      else if(LhsIsRowMajor)
-        dst_j += i.value() * rhs.row(i.index());
-      else if(Rhs::ColsAtCompileTime==1)
-        dst.coeffRef(i.index()) += i.value() * rhs_j;
-      else
-        dst.row(i.index()) += i.value() * rhs.row(j);
     }
-    if (ProcessFirstHalf && i && (i.index()==j))
-      dst.row(j) += i.value() * rhs.row(j);
-  }
-}
 
-template<typename Derived>
+  private:
+    SparseTimeDenseProduct& operator=(const SparseTimeDenseProduct&);
+};
+
+
+// dense = dense * sparse
 template<typename Lhs, typename Rhs>
-Derived& MatrixBase<Derived>::lazyAssign(const SparseProduct<Lhs,Rhs,SparseTimeDenseProduct>& product)
+struct ei_traits<DenseTimeSparseProduct<Lhs,Rhs> >
+ : ei_traits<ProductBase<DenseTimeSparseProduct<Lhs,Rhs>, Lhs, Rhs> >
 {
-  derived().setZero();
-  ei_sparse_time_dense_product(product.lhs(), product.rhs(), derived());
-  return derived();
-}
+  typedef Dense StorageType;
+};
 
-// dense = dense * sparse
-template<typename Derived>
 template<typename Lhs, typename Rhs>
-Derived& MatrixBase<Derived>::lazyAssign(const SparseProduct<Lhs,Rhs,DenseTimeSparseProduct>& product)
+class DenseTimeSparseProduct
+  : public ProductBase<DenseTimeSparseProduct<Lhs,Rhs>, Lhs, Rhs>
 {
-  typedef typename ei_cleantype<Rhs>::type _Rhs;
-  typedef typename _Rhs::InnerIterator RhsInnerIterator;
-  enum { RhsIsRowMajor = (_Rhs::Flags&RowMajorBit)==RowMajorBit };
-  derived().setZero();
-  for (int j=0; j<product.rhs().outerSize(); ++j)
-    for (RhsInnerIterator i(product.rhs(),j); i; ++i)
-      derived().col(RhsIsRowMajor ? i.index() : j) += i.value() * product.lhs().col(RhsIsRowMajor ? j : i.index());
-  return derived();
-}
+  public:
+    EIGEN_PRODUCT_PUBLIC_INTERFACE(DenseTimeSparseProduct)
+
+    DenseTimeSparseProduct(const Lhs& lhs, const Rhs& rhs) : Base(lhs,rhs)
+    {}
+
+    template<typename Dest> void scaleAndAddTo(Dest& dest, Scalar alpha) const
+    {
+      typedef typename ei_cleantype<Rhs>::type _Rhs;
+      typedef typename _Rhs::InnerIterator RhsInnerIterator;
+      enum { RhsIsRowMajor = (_Rhs::Flags&RowMajorBit)==RowMajorBit };
+      for(int j=0; j<m_rhs.outerSize(); ++j)
+        for(RhsInnerIterator i(m_rhs,j); i; ++i)
+          dest.col(RhsIsRowMajor ? i.index() : j) += (alpha*i.value()) * m_lhs.col(RhsIsRowMajor ? j : i.index());
+    }
+
+  private:
+    DenseTimeSparseProduct& operator=(const DenseTimeSparseProduct&);
+};
 
 // sparse * sparse
 template<typename Derived>
@@ -381,10 +341,10 @@ SparseMatrixBase<Derived>::operator*(const SparseMatrixBase<OtherDerived> &other
 // sparse * dense
 template<typename Derived>
 template<typename OtherDerived>
-EIGEN_STRONG_INLINE const typename SparseProductReturnType<Derived,OtherDerived>::Type
+EIGEN_STRONG_INLINE const SparseTimeDenseProduct<Derived,OtherDerived>
 SparseMatrixBase<Derived>::operator*(const MatrixBase<OtherDerived> &other) const
 {
-  return typename SparseProductReturnType<Derived,OtherDerived>::Type(derived(), other.derived());
+  return SparseTimeDenseProduct<Derived,OtherDerived>(derived(), other.derived());
 }
 
 #endif // EIGEN_SPARSEPRODUCT_H
diff --git a/Eigen/src/Sparse/SparseSelfAdjointView.h b/Eigen/src/Sparse/SparseSelfAdjointView.h
new file mode 100644
index 000000000..f5296accf
--- /dev/null
+++ b/Eigen/src/Sparse/SparseSelfAdjointView.h
@@ -0,0 +1,223 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra.
+//
+// Copyright (C) 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
+// License as published by the Free Software Foundation; either
+// version 3 of the License, or (at your option) any later version.
+//
+// Alternatively, you can redistribute it and/or
+// modify it under the terms of the GNU General Public License as
+// published by the Free Software Foundation; either version 2 of
+// the License, or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License and a copy of the GNU General Public License along with
+// Eigen. If not, see <http://www.gnu.org/licenses/>.
+
+#ifndef EIGEN_SPARSE_SELFADJOINTVIEW_H
+#define EIGEN_SPARSE_SELFADJOINTVIEW_H
+
+/** \class SparseSelfAdjointView
+  * \nonstableyet
+  *
+  * \brief Pseudo expression to manipulate a triangular sparse matrix as a selfadjoint matrix.
+  *
+  * \param MatrixType the type of the dense matrix storing the coefficients
+  * \param UpLo can be either \c LowerTriangular or \c UpperTriangular
+  *
+  * This class is an expression of a sefladjoint matrix from a triangular part of a matrix
+  * with given dense storage of the coefficients. It is the return type of MatrixBase::selfadjointView()
+  * and most of the time this is the only way that it is used.
+  *
+  * \sa SparseMatrixBase::selfAdjointView()
+  */
+template<typename Lhs, typename Rhs, int UpLo>
+class SparseSelfAdjointTimeDenseProduct;
+
+template<typename Lhs, typename Rhs, int UpLo>
+class DenseTimeSparseSelfAdjointProduct;
+
+template<typename MatrixType, unsigned int UpLo> class SparseSelfAdjointView
+{
+  public:
+
+    typedef typename ei_traits<MatrixType>::Scalar Scalar;
+
+    inline SparseSelfAdjointView(const MatrixType& matrix) : m_matrix(matrix)
+    {
+      ei_assert(ei_are_flags_consistent<UpLo>::ret);
+      ei_assert(rows()==cols() && "SelfAdjointView is only for squared matrices");
+    }
+
+    inline int rows() const { return m_matrix.rows(); }
+    inline int cols() const { return m_matrix.cols(); }
+
+    /** \internal \returns a reference to the nested matrix */
+    const MatrixType& matrix() const { return m_matrix; }
+
+    /** Efficient sparse self-adjoint matrix times dense vector/matrix product */
+    template<typename OtherDerived>
+    SparseSelfAdjointTimeDenseProduct<MatrixType,OtherDerived,UpLo>
+    operator*(const MatrixBase<OtherDerived>& rhs) const
+    {
+      return SparseSelfAdjointTimeDenseProduct<MatrixType,OtherDerived,UpLo>(m_matrix, rhs.derived());
+    }
+
+    /** Efficient dense vector/matrix times sparse self-adjoint matrix product */
+    template<typename OtherDerived> friend
+    DenseTimeSparseSelfAdjointProduct<OtherDerived,MatrixType,UpLo>
+    operator*(const MatrixBase<OtherDerived>& lhs, const SparseSelfAdjointView& rhs)
+    {
+      return DenseTimeSparseSelfAdjointProduct<OtherDerived,MatrixType,UpLo>(lhs.derived(), rhs.m_matrix);
+    }
+
+    /** Perform a symmetric rank K update of the selfadjoint matrix \c *this:
+      * \f$ this = this + \alpha ( u u^* ) \f$ where \a u is a vector or matrix.
+      *
+      * \returns a reference to \c *this
+      *
+      * Note that it is faster to set alpha=0 than initializing the matrix to zero
+      * and then keep the default value alpha=1.
+      *
+      * To perform \f$ this = this + \alpha ( u^* u ) \f$ you can simply
+      * call this function with u.adjoint().
+      */
+    template<typename DerivedU>
+    SparseSelfAdjointView& rankUpdate(const MatrixBase<DerivedU>& u, Scalar alpha = Scalar(1));
+
+    // const SparseLLT<PlainMatrixType, UpLo> llt() const;
+    // const SparseLDLT<PlainMatrixType, UpLo> ldlt() const;
+
+  protected:
+
+    const typename MatrixType::Nested m_matrix;
+};
+
+/***************************************************************************
+* Implementation of SparseMatrixBase methods
+***************************************************************************/
+
+template<typename Derived>
+template<unsigned int UpLo>
+const SparseSelfAdjointView<Derived, UpLo> SparseMatrixBase<Derived>::selfadjointView() const
+{
+  return derived();
+}
+
+template<typename Derived>
+template<unsigned int UpLo>
+SparseSelfAdjointView<Derived, UpLo> SparseMatrixBase<Derived>::selfadjointView()
+{
+  return derived();
+}
+
+/***************************************************************************
+* Implementation of SparseSelfAdjointView methods
+***************************************************************************/
+
+template<typename MatrixType, unsigned int UpLo>
+template<typename DerivedU>
+SparseSelfAdjointView<MatrixType,UpLo>&
+SparseSelfAdjointView<MatrixType,UpLo>::rankUpdate(const MatrixBase<DerivedU>& u, Scalar alpha)
+{
+  SparseMatrix<Scalar,MatrixType::Flags&RowMajorBit?RowMajor:ColMajor> tmp = u * u.adjoint();
+  if(alpha==Scalar(0))
+    m_matrix = tmp.template triangularView<UpLo>();
+  else
+    m_matrix += alpha * tmp.template triangularView<UpLo>();
+
+  return this;
+}
+
+/***************************************************************************
+* Implementation of sparse self-adjoint time dense matrix
+***************************************************************************/
+
+template<typename Lhs, typename Rhs, int UpLo>
+struct ei_traits<SparseSelfAdjointTimeDenseProduct<Lhs,Rhs,UpLo> >
+ : ei_traits<ProductBase<SparseSelfAdjointTimeDenseProduct<Lhs,Rhs,UpLo>, Lhs, Rhs> >
+{};
+
+template<typename Lhs, typename Rhs, int UpLo>
+class SparseSelfAdjointTimeDenseProduct
+  : public ProductBase<SparseSelfAdjointTimeDenseProduct<Lhs,Rhs,UpLo>, Lhs, Rhs>
+{
+  public:
+    EIGEN_PRODUCT_PUBLIC_INTERFACE(SparseSelfAdjointTimeDenseProduct)
+
+    SparseSelfAdjointTimeDenseProduct(const Lhs& lhs, const Rhs& rhs) : Base(lhs,rhs)
+    {}
+
+    template<typename Dest> void scaleAndAddTo(Dest& dest, Scalar alpha) const
+    {
+      // TODO use alpha
+      ei_assert(alpha==Scalar(1) && "alpha != 1 is not implemented yet, sorry");
+      typedef typename ei_cleantype<Lhs>::type _Lhs;
+      typedef typename ei_cleantype<Rhs>::type _Rhs;
+      typedef typename _Lhs::InnerIterator LhsInnerIterator;
+      enum {
+        LhsIsRowMajor = (_Lhs::Flags&RowMajorBit)==RowMajorBit,
+        ProcessFirstHalf =
+                 ((UpLo&(UpperTriangularBit|LowerTriangularBit))==(UpperTriangularBit|LowerTriangularBit))
+              || ( (UpLo&UpperTriangularBit) && !LhsIsRowMajor)
+              || ( (UpLo&LowerTriangularBit) && LhsIsRowMajor),
+        ProcessSecondHalf = !ProcessFirstHalf
+      };
+      for (int j=0; j<m_lhs.outerSize(); ++j)
+      {
+        LhsInnerIterator i(m_lhs,j);
+        if (ProcessSecondHalf && i && (i.index()==j))
+        {
+          dest.row(j) += i.value() * m_rhs.row(j);
+          ++i;
+        }
+        Block<Dest,1,Dest::ColsAtCompileTime> dest_j(dest.row(LhsIsRowMajor ? j : 0));
+        for(; (ProcessFirstHalf ? i && i.index() < j : i) ; ++i)
+        {
+          int a = LhsIsRowMajor ? j : i.index();
+          int b = LhsIsRowMajor ? i.index() : j;
+          typename Lhs::Scalar v = i.value();
+          dest.row(a) += (v) * m_rhs.row(b);
+          dest.row(b) += ei_conj(v) * m_rhs.row(a);
+        }
+        if (ProcessFirstHalf && i && (i.index()==j))
+          dest.row(j) += i.value() * m_rhs.row(j);
+      }
+    }
+
+  private:
+    SparseSelfAdjointTimeDenseProduct& operator=(const SparseSelfAdjointTimeDenseProduct&);
+};
+
+template<typename Lhs, typename Rhs, int UpLo>
+struct ei_traits<DenseTimeSparseSelfAdjointProduct<Lhs,Rhs,UpLo> >
+ : ei_traits<ProductBase<DenseTimeSparseSelfAdjointProduct<Lhs,Rhs,UpLo>, Lhs, Rhs> >
+{};
+
+template<typename Lhs, typename Rhs, int UpLo>
+class DenseTimeSparseSelfAdjointProduct
+  : public ProductBase<DenseTimeSparseSelfAdjointProduct<Lhs,Rhs,UpLo>, Lhs, Rhs>
+{
+  public:
+    EIGEN_PRODUCT_PUBLIC_INTERFACE(DenseTimeSparseSelfAdjointProduct)
+
+    DenseTimeSparseSelfAdjointProduct(const Lhs& lhs, const Rhs& rhs) : Base(lhs,rhs)
+    {}
+
+    template<typename Dest> void scaleAndAddTo(Dest& dest, Scalar alpha) const
+    {
+      // TODO
+    }
+
+  private:
+    DenseTimeSparseSelfAdjointProduct& operator=(const DenseTimeSparseSelfAdjointProduct&);
+};
+#endif // EIGEN_SPARSE_SELFADJOINTVIEW_H
diff --git a/Eigen/src/Sparse/SparseTriangular.h b/Eigen/src/Sparse/SparseTriangular.h
deleted file mode 100644
index 42e7ff02a..000000000
--- a/Eigen/src/Sparse/SparseTriangular.h
+++ /dev/null
@@ -1,60 +0,0 @@
-// This file is part of Eigen, a lightweight C++ template library
-// for linear algebra.
-//
-// Copyright (C) 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
-// License as published by the Free Software Foundation; either
-// version 3 of the License, or (at your option) any later version.
-//
-// Alternatively, you can redistribute it and/or
-// modify it under the terms of the GNU General Public License as
-// published by the Free Software Foundation; either version 2 of
-// the License, or (at your option) any later version.
-//
-// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
-// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
-// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
-// GNU General Public License for more details.
-//
-// You should have received a copy of the GNU Lesser General Public
-// License and a copy of the GNU General Public License along with
-// Eigen. If not, see <http://www.gnu.org/licenses/>.
-
-#ifndef EIGEN_SPARSE_TRIANGULAR_H
-#define EIGEN_SPARSE_TRIANGULAR_H
-
-template<typename ExpressionType, int Mode> class SparseTriangular
-{
-  public:
-
-    typedef typename ei_traits<ExpressionType>::Scalar Scalar;
-    typedef typename ei_meta_if<ei_must_nest_by_value<ExpressionType>::ret,
-        ExpressionType, const ExpressionType&>::ret ExpressionTypeNested;
-    typedef CwiseUnaryOp<ei_scalar_add_op<Scalar>, ExpressionType> ScalarAddReturnType;
-
-    inline SparseTriangular(const ExpressionType& matrix) : m_matrix(matrix) {}
-
-    /** \internal */
-    inline const ExpressionType& _expression() const { return m_matrix; }
-
-    template<typename OtherDerived>
-    typename ei_plain_matrix_type_column_major<OtherDerived>::type
-    solve(const MatrixBase<OtherDerived>& other) const;
-    template<typename OtherDerived> void solveInPlace(MatrixBase<OtherDerived>& other) const;
-    template<typename OtherDerived> void solveInPlace(SparseMatrixBase<OtherDerived>& other) const;
-    
-  protected:
-    ExpressionTypeNested m_matrix;
-};
-
-template<typename Derived>
-template<int Mode>
-inline const SparseTriangular<Derived, Mode>
-SparseMatrixBase<Derived>::triangular() const
-{
-  return derived();
-}
-
-#endif // EIGEN_SPARSE_TRIANGULAR_H
diff --git a/Eigen/src/Sparse/SparseTriangularView.h b/Eigen/src/Sparse/SparseTriangularView.h
new file mode 100644
index 000000000..6a9461528
--- /dev/null
+++ b/Eigen/src/Sparse/SparseTriangularView.h
@@ -0,0 +1,95 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra.
+//
+// Copyright (C) 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
+// License as published by the Free Software Foundation; either
+// version 3 of the License, or (at your option) any later version.
+//
+// Alternatively, you can redistribute it and/or
+// modify it under the terms of the GNU General Public License as
+// published by the Free Software Foundation; either version 2 of
+// the License, or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License and a copy of the GNU General Public License along with
+// Eigen. If not, see <http://www.gnu.org/licenses/>.
+
+#ifndef EIGEN_SPARSE_TRIANGULARVIEW_H
+#define EIGEN_SPARSE_TRIANGULARVIEW_H
+
+template<typename MatrixType, int Mode>
+struct ei_traits<SparseTriangularView<MatrixType,Mode> >
+: public ei_traits<MatrixType>
+{};
+
+template<typename MatrixType, int Mode> class SparseTriangularView
+  : public SparseMatrixBase<SparseTriangularView<MatrixType,Mode> >
+{
+    enum { SkipFirst = (Mode==LowerTriangular && (!MatrixType::Flags&RowMajorBit))
+                    || (Mode==UpperTriangular && ( MatrixType::Flags&RowMajorBit)) };
+  public:
+
+    class InnerIterator;
+
+    inline int rows() { return m_matrix.rows(); }
+    inline int cols() { return m_matrix.cols(); }
+
+    typedef typename ei_traits<MatrixType>::Scalar Scalar;
+    typedef typename ei_meta_if<ei_must_nest_by_value<MatrixType>::ret,
+        MatrixType, const MatrixType&>::ret MatrixTypeNested;
+
+    inline SparseTriangularView(const MatrixType& matrix) : m_matrix(matrix) {}
+
+    /** \internal */
+    inline const MatrixType& nestedExpression() const { return m_matrix; }
+
+    template<typename OtherDerived>
+    typename ei_plain_matrix_type_column_major<OtherDerived>::type
+    solve(const MatrixBase<OtherDerived>& other) const;
+
+    template<typename OtherDerived> void solveInPlace(MatrixBase<OtherDerived>& other) const;
+    template<typename OtherDerived> void solveInPlace(SparseMatrixBase<OtherDerived>& other) const;
+
+  protected:
+    MatrixTypeNested m_matrix;
+};
+
+template<typename MatrixType, int Mode>
+class SparseTriangularView<MatrixType,Mode>::InnerIterator : public MatrixType::InnerIterator
+{
+    typedef typename MatrixType::InnerIterator Base;
+  public:
+
+    EIGEN_STRONG_INLINE InnerIterator(const SparseTriangularView& view, int outer)
+      : Base(view.nestedExpression(), outer)
+    {
+      if(SkipFirst)
+        while((*this) && this->index()<outer)
+          ++(*this);
+    }
+    inline int row() const { return Base::row(); }
+    inline int col() const { return Base::col(); }
+
+    EIGEN_STRONG_INLINE operator bool() const
+    {
+      return SkipFirst ? Base::operator bool() : (Base::operator bool() && this->index() < this->outer());
+    }
+};
+
+template<typename Derived>
+template<int Mode>
+inline const SparseTriangularView<Derived, Mode>
+SparseMatrixBase<Derived>::triangularView() const
+{
+  return derived();
+}
+
+#endif // EIGEN_SPARSE_TRIANGULARVIEW_H
diff --git a/Eigen/src/Sparse/SparseUtil.h b/Eigen/src/Sparse/SparseUtil.h
index b99be580c..6da1bee25 100644
--- a/Eigen/src/Sparse/SparseUtil.h
+++ b/Eigen/src/Sparse/SparseUtil.h
@@ -118,25 +118,29 @@ enum {
 };
 
 template<typename Derived> class SparseMatrixBase;
-template<typename _Scalar, int _Flags = 0> class SparseMatrix;
-template<typename _Scalar, int _Flags = 0> class DynamicSparseMatrix;
-template<typename _Scalar, int _Flags = 0> class SparseVector;
-template<typename _Scalar, int _Flags = 0> class MappedSparseMatrix;
-
-template<typename MatrixType>                            class SparseNestByValue;
-template<typename MatrixType, int Size>                  class SparseInnerVectorSet;
-template<typename ExpressionType,
-         unsigned int Added, unsigned int Removed>       class SparseFlagged;
-template<typename ExpressionType, int Mode>              class SparseTriangular;
-template<typename Lhs, typename Rhs>                     class SparseDiagonalProduct;
+template<typename _Scalar, int _Flags = 0>  class SparseMatrix;
+template<typename _Scalar, int _Flags = 0>  class DynamicSparseMatrix;
+template<typename _Scalar, int _Flags = 0>  class SparseVector;
+template<typename _Scalar, int _Flags = 0>  class MappedSparseMatrix;
+
+template<typename MatrixType>               class SparseNestByValue;
+template<typename MatrixType, int Size>     class SparseInnerVectorSet;
+template<typename MatrixType, int Mode>     class SparseTriangularView;
+template<typename MatrixType, unsigned int UpLo>     class SparseSelfAdjointView;
+template<typename Lhs, typename Rhs>        class SparseDiagonalProduct;
+
+
+template<typename Lhs, typename Rhs>        class SparseProduct;
+template<typename Lhs, typename Rhs>        class SparseTimeDenseProduct;
+template<typename Lhs, typename Rhs>        class DenseTimeSparseProduct;
 
 template<typename Lhs, typename Rhs,
          typename LhsStorage = typename ei_traits<Lhs>::StorageType,
          typename RhsStorage = typename ei_traits<Rhs>::StorageType> struct ei_sparse_product_mode;
 
-template<typename Lhs, typename Rhs, int ProductMode = ei_sparse_product_mode<Lhs,Rhs>::value> struct SparseProductReturnType;
+template<typename Lhs, typename Rhs> struct SparseProductReturnType;
 
-const int CoherentAccessPattern  = 0x1;
+const int CoherentAccessPattern     = 0x1;
 const int InnerRandomAccessPattern  = 0x2 | CoherentAccessPattern;
 const int OuterRandomAccessPattern  = 0x4 | CoherentAccessPattern;
 const int RandomAccessPattern       = 0x8 | OuterRandomAccessPattern | InnerRandomAccessPattern;
diff --git a/Eigen/src/Sparse/TriangularSolver.h b/Eigen/src/Sparse/TriangularSolver.h
index cb2bc4d58..2c5b485b0 100644
--- a/Eigen/src/Sparse/TriangularSolver.h
+++ b/Eigen/src/Sparse/TriangularSolver.h
@@ -160,7 +160,7 @@ struct ei_sparse_solve_triangular_selector<Lhs,Rhs,Mode,UpperTriangular,ColMajor
 
 template<typename ExpressionType,int Mode>
 template<typename OtherDerived>
-void SparseTriangular<ExpressionType,Mode>::solveInPlace(MatrixBase<OtherDerived>& other) const
+void SparseTriangularView<ExpressionType,Mode>::solveInPlace(MatrixBase<OtherDerived>& other) const
 {
   ei_assert(m_matrix.cols() == m_matrix.rows());
   ei_assert(m_matrix.cols() == other.rows());
@@ -182,7 +182,7 @@ void SparseTriangular<ExpressionType,Mode>::solveInPlace(MatrixBase<OtherDerived
 template<typename ExpressionType,int Mode>
 template<typename OtherDerived>
 typename ei_plain_matrix_type_column_major<OtherDerived>::type
-SparseTriangular<ExpressionType,Mode>::solve(const MatrixBase<OtherDerived>& other) const
+SparseTriangularView<ExpressionType,Mode>::solve(const MatrixBase<OtherDerived>& other) const
 {
   typename ei_plain_matrix_type_column_major<OtherDerived>::type res(other);
   solveInPlace(res);
@@ -210,10 +210,10 @@ struct ei_sparse_solve_triangular_sparse_selector<Lhs,Rhs,Mode,UpLo,ColMajor>
     const bool IsLowerTriangular = (UpLo==LowerTriangular);
     AmbiVector<Scalar> tempVector(other.rows()*2);
     tempVector.setBounds(0,other.rows());
-    
+
     Rhs res(other.rows(), other.cols());
     res.reserve(other.nonZeros());
-    
+
     for(int col=0 ; col<other.cols() ; ++col)
     {
       // FIXME estimate number of non zeros
@@ -224,7 +224,7 @@ struct ei_sparse_solve_triangular_sparse_selector<Lhs,Rhs,Mode,UpLo,ColMajor>
       {
         tempVector.coeffRef(rhsIt.index()) = rhsIt.value();
       }
-      
+
       for(int i=IsLowerTriangular?0:lhs.cols()-1;
           IsLowerTriangular?i<lhs.cols():i>=0;
           i+=IsLowerTriangular?1:-1)
@@ -233,7 +233,7 @@ struct ei_sparse_solve_triangular_sparse_selector<Lhs,Rhs,Mode,UpLo,ColMajor>
         Scalar& ci = tempVector.coeffRef(i);
         if (ci!=Scalar(0))
         {
-          // find 
+          // find
           typename Lhs::InnerIterator it(lhs, i);
           if(!(Mode & UnitDiagBit))
           {
@@ -260,8 +260,8 @@ struct ei_sparse_solve_triangular_sparse_selector<Lhs,Rhs,Mode,UpLo,ColMajor>
           }
         }
       }
-      
-      
+
+
       int count = 0;
       // FIXME compute a reference value to filter zeros
       for (typename AmbiVector<Scalar>::Iterator it(tempVector/*,1e-12*/); it; ++it)
@@ -281,7 +281,7 @@ struct ei_sparse_solve_triangular_sparse_selector<Lhs,Rhs,Mode,UpLo,ColMajor>
 
 template<typename ExpressionType,int Mode>
 template<typename OtherDerived>
-void SparseTriangular<ExpressionType,Mode>::solveInPlace(SparseMatrixBase<OtherDerived>& other) const
+void SparseTriangularView<ExpressionType,Mode>::solveInPlace(SparseMatrixBase<OtherDerived>& other) const
 {
   ei_assert(m_matrix.cols() == m_matrix.rows());
   ei_assert(m_matrix.cols() == other.rows());
diff --git a/test/sparse_product.cpp b/test/sparse_product.cpp
index 3f0e793d5..e944d6c53 100644
--- a/test/sparse_product.cpp
+++ b/test/sparse_product.cpp
@@ -114,10 +114,10 @@ template<typename SparseMatrixType> void sparse_product(const SparseMatrixType&
     VERIFY_IS_APPROX(mS.transpose().conjugate(), mS);
     VERIFY_IS_APPROX(mS, refS);
     VERIFY_IS_APPROX(x=mS*b, refX=refS*b);
-    // TODO properly implement triangular/selfadjoint views
-//     VERIFY_IS_APPROX(x=mUp.template marked<UpperTriangular|SelfAdjoint>()*b, refX=refS*b);
-//     VERIFY_IS_APPROX(x=mLo.template marked<LowerTriangular|SelfAdjoint>()*b, refX=refS*b);
-//     VERIFY_IS_APPROX(x=mS.template marked<SelfAdjoint>()*b, refX=refS*b);
+
+    VERIFY_IS_APPROX(x=mUp.template selfadjointView<UpperTriangular>()*b, refX=refS*b);
+    VERIFY_IS_APPROX(x=mLo.template selfadjointView<LowerTriangular>()*b, refX=refS*b);
+    VERIFY_IS_APPROX(x=mS.template selfadjointView<UpperTriangular|LowerTriangular>()*b, refX=refS*b);
   }
 }
 
diff --git a/test/sparse_solvers.cpp b/test/sparse_solvers.cpp
index a530a9515..24107977c 100644
--- a/test/sparse_solvers.cpp
+++ b/test/sparse_solvers.cpp
@@ -66,12 +66,12 @@ template<typename Scalar> void sparse_solvers(int rows, int cols)
     // lower - dense
     initSparse<Scalar>(density, refMat2, m2, ForceNonZeroDiag|MakeLowerTriangular, &zeroCoords, &nonzeroCoords);
     VERIFY_IS_APPROX(refMat2.template triangularView<LowerTriangular>().solve(vec2),
-                     m2.template triangular<LowerTriangular>().solve(vec3));
+                     m2.template triangularView<LowerTriangular>().solve(vec3));
 
     // upper - dense
     initSparse<Scalar>(density, refMat2, m2, ForceNonZeroDiag|MakeUpperTriangular, &zeroCoords, &nonzeroCoords);
     VERIFY_IS_APPROX(refMat2.template triangularView<UpperTriangular>().solve(vec2),
-                     m2.template triangular<UpperTriangular>().solve(vec3));
+                     m2.template triangularView<UpperTriangular>().solve(vec3));
 
     // TODO test row major
 
@@ -82,20 +82,20 @@ template<typename Scalar> void sparse_solvers(int rows, int cols)
     initSparse<Scalar>(density, refMat2, m2, ForceNonZeroDiag|MakeLowerTriangular);
     initSparse<Scalar>(density, refMatB, matB);
     refMat2.template triangularView<LowerTriangular>().solveInPlace(refMatB);
-    m2.template triangular<LowerTriangular>().solveInPlace(matB);
+    m2.template triangularView<LowerTriangular>().solveInPlace(matB);
     VERIFY_IS_APPROX(matB.toDense(), refMatB);
 
     // upper - sparse
     initSparse<Scalar>(density, refMat2, m2, ForceNonZeroDiag|MakeUpperTriangular);
     initSparse<Scalar>(density, refMatB, matB);
     refMat2.template triangularView<UpperTriangular>().solveInPlace(refMatB);
-    m2.template triangular<UpperTriangular>().solveInPlace(matB);
+    m2.template triangularView<UpperTriangular>().solveInPlace(matB);
     VERIFY_IS_APPROX(matB, refMatB);
 
     // test deprecated API
     initSparse<Scalar>(density, refMat2, m2, ForceNonZeroDiag|MakeLowerTriangular, &zeroCoords, &nonzeroCoords);
     VERIFY_IS_APPROX(refMat2.template triangularView<LowerTriangular>().solve(vec2),
-                     m2.template triangular<LowerTriangular>().solve(vec3));
+                     m2.template triangularView<LowerTriangular>().solve(vec3));
   }
 
   // test LLT