diff options
Diffstat (limited to 'Eigen/src/SparseCore/SparseMatrixBase.h')
| -rw-r--r-- | Eigen/src/SparseCore/SparseMatrixBase.h | 285 |
1 files changed, 118 insertions, 167 deletions
diff --git a/Eigen/src/SparseCore/SparseMatrixBase.h b/Eigen/src/SparseCore/SparseMatrixBase.h index 6b2169a7b..c6b548f11 100644 --- a/Eigen/src/SparseCore/SparseMatrixBase.h +++ b/Eigen/src/SparseCore/SparseMatrixBase.h @@ -1,7 +1,7 @@ // This file is part of Eigen, a lightweight C++ template library // for linear algebra. // -// Copyright (C) 2008-2011 Gael Guennebaud <gael.guennebaud@inria.fr> +// Copyright (C) 2008-2014 Gael Guennebaud <gael.guennebaud@inria.fr> // // This Source Code Form is subject to the terms of the Mozilla // Public License v. 2.0. If a copy of the MPL was not distributed @@ -18,32 +18,41 @@ namespace Eigen { * * \brief Base class of any sparse matrices or sparse expressions * - * \tparam Derived + * \tparam Derived is the derived type, e.g. a sparse matrix type, or an expression, etc. * * This class can be extended with the help of the plugin mechanism described on the page - * \ref TopicCustomizingEigen by defining the preprocessor symbol \c EIGEN_SPARSEMATRIXBASE_PLUGIN. + * \ref TopicCustomizing_Plugins by defining the preprocessor symbol \c EIGEN_SPARSEMATRIXBASE_PLUGIN. */ -template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> +template<typename Derived> class SparseMatrixBase + : public EigenBase<Derived> { public: typedef typename internal::traits<Derived>::Scalar Scalar; + + /** The numeric type of the expression' coefficients, e.g. float, double, int or std::complex<float>, etc. + * + * It is an alias for the Scalar type */ + typedef Scalar value_type; + typedef typename internal::packet_traits<Scalar>::type PacketScalar; typedef typename internal::traits<Derived>::StorageKind StorageKind; - typedef typename internal::traits<Derived>::Index Index; + + /** The integer type used to \b store indices within a SparseMatrix. + * For a \c SparseMatrix<Scalar,Options,IndexType> it an alias of the third template parameter \c IndexType. */ + typedef typename internal::traits<Derived>::StorageIndex StorageIndex; + typedef typename internal::add_const_on_value_type_if_arithmetic< typename internal::packet_traits<Scalar>::type >::type PacketReturnType; typedef SparseMatrixBase StorageBaseType; - typedef EigenBase<Derived> Base; + + typedef Matrix<StorageIndex,Dynamic,1> IndexVector; + typedef Matrix<Scalar,Dynamic,1> ScalarVector; template<typename OtherDerived> - Derived& operator=(const EigenBase<OtherDerived> &other) - { - other.derived().evalTo(derived()); - return derived(); - } + Derived& operator=(const EigenBase<OtherDerived> &other); enum { @@ -83,11 +92,6 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> * constructed from this one. See the \ref flags "list of flags". */ - CoeffReadCost = internal::traits<Derived>::CoeffReadCost, - /**< This is a rough measure of how expensive it is to read one coefficient from - * this expression. - */ - IsRowMajor = Flags&RowMajorBit ? 1 : 0, InnerSizeAtCompileTime = int(IsVectorAtCompileTime) ? int(SizeAtCompileTime) @@ -103,10 +107,11 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> CwiseUnaryOp<internal::scalar_conjugate_op<Scalar>, Eigen::Transpose<const Derived> >, Transpose<const Derived> >::type AdjointReturnType; + typedef Transpose<Derived> TransposeReturnType; + typedef typename internal::add_const<Transpose<const Derived> >::type ConstTransposeReturnType; - - typedef SparseMatrix<Scalar, Flags&RowMajorBit ? RowMajor : ColMajor, Index> PlainObject; - + // FIXME storage order do not match evaluator storage order + typedef SparseMatrix<Scalar, Flags&RowMajorBit ? RowMajor : ColMajor, StorageIndex> PlainObject; #ifndef EIGEN_PARSED_BY_DOXYGEN /** This is the "real scalar" type; if the \a Scalar type is already real numbers @@ -124,6 +129,8 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> /** \internal Represents a matrix with all coefficients equal to one another*/ typedef CwiseNullaryOp<internal::scalar_constant_op<Scalar>,Matrix<Scalar,Dynamic,Dynamic> > ConstantReturnType; + /** type of the equivalent dense matrix */ + typedef Matrix<Scalar,RowsAtCompileTime,ColsAtCompileTime> DenseMatrixType; /** type of the equivalent square matrix */ typedef Matrix<Scalar,EIGEN_SIZE_MAX(RowsAtCompileTime,ColsAtCompileTime), EIGEN_SIZE_MAX(RowsAtCompileTime,ColsAtCompileTime)> SquareMatrixType; @@ -132,9 +139,21 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> inline Derived& derived() { return *static_cast<Derived*>(this); } inline Derived& const_cast_derived() const { return *static_cast<Derived*>(const_cast<SparseMatrixBase*>(this)); } + + typedef EigenBase<Derived> Base; + #endif // not EIGEN_PARSED_BY_DOXYGEN #define EIGEN_CURRENT_STORAGE_BASE_CLASS Eigen::SparseMatrixBase +#ifdef EIGEN_PARSED_BY_DOXYGEN +#define EIGEN_DOC_UNARY_ADDONS(METHOD,OP) /** <p>This method does not change the sparsity of \c *this: the OP is applied to explicitly stored coefficients only. \sa SparseCompressedBase::coeffs() </p> */ +#define EIGEN_DOC_BLOCK_ADDONS_NOT_INNER_PANEL /** <p> \warning This method returns a read-only expression for any sparse matrices. \sa \ref TutorialSparse_SubMatrices "Sparse block operations" </p> */ +#define EIGEN_DOC_BLOCK_ADDONS_INNER_PANEL_IF(COND) /** <p> \warning This method returns a read-write expression for COND sparse matrices only. Otherwise, the returned expression is read-only. \sa \ref TutorialSparse_SubMatrices "Sparse block operations" </p> */ +#else +#define EIGEN_DOC_UNARY_ADDONS(X,Y) +#define EIGEN_DOC_BLOCK_ADDONS_NOT_INNER_PANEL +#define EIGEN_DOC_BLOCK_ADDONS_INNER_PANEL_IF(COND) +#endif # include "../plugins/CommonCwiseUnaryOps.h" # include "../plugins/CommonCwiseBinaryOps.h" # include "../plugins/MatrixCwiseUnaryOps.h" @@ -143,8 +162,10 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> # ifdef EIGEN_SPARSEMATRIXBASE_PLUGIN # include EIGEN_SPARSEMATRIXBASE_PLUGIN # endif -# undef EIGEN_CURRENT_STORAGE_BASE_CLASS #undef EIGEN_CURRENT_STORAGE_BASE_CLASS +#undef EIGEN_DOC_UNARY_ADDONS +#undef EIGEN_DOC_BLOCK_ADDONS_NOT_INNER_PANEL +#undef EIGEN_DOC_BLOCK_ADDONS_INNER_PANEL_IF /** \returns the number of rows. \sa cols() */ inline Index rows() const { return derived().rows(); } @@ -153,9 +174,6 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> /** \returns the number of coefficients, which is \a rows()*cols(). * \sa rows(), cols(). */ inline Index size() const { return rows() * cols(); } - /** \returns the number of nonzero coefficients which is in practice the number - * of stored coefficients. */ - inline Index nonZeros() const { return derived().nonZeros(); } /** \returns true if either the number of rows or the number of columns is equal to 1. * In other words, this function returns * \code rows()==1 || cols()==1 \endcode @@ -175,93 +193,23 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> template<typename OtherDerived> - Derived& operator=(const ReturnByValue<OtherDerived>& other) - { - other.evalTo(derived()); - return derived(); - } - + Derived& operator=(const ReturnByValue<OtherDerived>& other); template<typename OtherDerived> - inline Derived& operator=(const SparseMatrixBase<OtherDerived>& other) - { - return assign(other.derived()); - } + inline Derived& operator=(const SparseMatrixBase<OtherDerived>& other); - inline Derived& operator=(const Derived& other) - { -// if (other.isRValue()) -// derived().swap(other.const_cast_derived()); -// else - return assign(other.derived()); - } + inline Derived& operator=(const Derived& other); protected: template<typename OtherDerived> - inline Derived& assign(const OtherDerived& other) - { - const bool transpose = (Flags & RowMajorBit) != (OtherDerived::Flags & RowMajorBit); - const Index outerSize = (int(OtherDerived::Flags) & RowMajorBit) ? other.rows() : other.cols(); - if ((!transpose) && other.isRValue()) - { - // eval without temporary - derived().resize(other.rows(), other.cols()); - derived().setZero(); - derived().reserve((std::max)(this->rows(),this->cols())*2); - for (Index j=0; j<outerSize; ++j) - { - derived().startVec(j); - for (typename OtherDerived::InnerIterator it(other, j); it; ++it) - { - Scalar v = it.value(); - derived().insertBackByOuterInner(j,it.index()) = v; - } - } - derived().finalize(); - } - else - { - assignGeneric(other); - } - return derived(); - } + inline Derived& assign(const OtherDerived& other); template<typename OtherDerived> - inline void assignGeneric(const OtherDerived& other) - { - //const bool transpose = (Flags & RowMajorBit) != (OtherDerived::Flags & RowMajorBit); - eigen_assert(( ((internal::traits<Derived>::SupportedAccessPatterns&OuterRandomAccessPattern)==OuterRandomAccessPattern) || - (!((Flags & RowMajorBit) != (OtherDerived::Flags & RowMajorBit)))) && - "the transpose operation is supposed to be handled in SparseMatrix::operator="); - - enum { Flip = (Flags & RowMajorBit) != (OtherDerived::Flags & RowMajorBit) }; - - const Index outerSize = other.outerSize(); - //typedef typename internal::conditional<transpose, LinkedVectorMatrix<Scalar,Flags&RowMajorBit>, Derived>::type TempType; - // thanks to shallow copies, we always eval to a tempary - Derived temp(other.rows(), other.cols()); - - temp.reserve((std::max)(this->rows(),this->cols())*2); - for (Index j=0; j<outerSize; ++j) - { - temp.startVec(j); - for (typename OtherDerived::InnerIterator it(other.derived(), j); it; ++it) - { - Scalar v = it.value(); - temp.insertBackByOuterInner(Flip?it.index():j,Flip?j:it.index()) = v; - } - } - temp.finalize(); - - derived() = temp.markAsRValue(); - } + inline void assignGeneric(const OtherDerived& other); public: - template<typename Lhs, typename Rhs> - inline Derived& operator=(const SparseSparseProduct<Lhs,Rhs>& product); - friend std::ostream & operator << (std::ostream & s, const SparseMatrixBase& m) { typedef typename Derived::Nested Nested; @@ -269,11 +217,12 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> if (Flags&RowMajorBit) { - const Nested nm(m.derived()); + Nested nm(m.derived()); + internal::evaluator<NestedCleaned> thisEval(nm); for (Index row=0; row<nm.outerSize(); ++row) { Index col = 0; - for (typename NestedCleaned::InnerIterator it(nm.derived(), row); it; ++it) + for (typename internal::evaluator<NestedCleaned>::InnerIterator it(thisEval, row); it; ++it) { for ( ; col<it.index(); ++col) s << "0 "; @@ -287,10 +236,11 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> } else { - const Nested nm(m.derived()); + Nested nm(m.derived()); + internal::evaluator<NestedCleaned> thisEval(nm); if (m.cols() == 1) { Index row = 0; - for (typename NestedCleaned::InnerIterator it(nm.derived(), 0); it; ++it) + for (typename internal::evaluator<NestedCleaned>::InnerIterator it(thisEval, 0); it; ++it) { for ( ; row<it.index(); ++row) s << "0" << std::endl; @@ -302,8 +252,8 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> } else { - SparseMatrix<Scalar, RowMajorBit, Index> trans = m; - s << static_cast<const SparseMatrixBase<SparseMatrix<Scalar, RowMajorBit, Index> >&>(trans); + SparseMatrix<Scalar, RowMajorBit, StorageIndex> trans = m; + s << static_cast<const SparseMatrixBase<SparseMatrix<Scalar, RowMajorBit, StorageIndex> >&>(trans); } } return s; @@ -313,56 +263,65 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> Derived& operator+=(const SparseMatrixBase<OtherDerived>& other); template<typename OtherDerived> Derived& operator-=(const SparseMatrixBase<OtherDerived>& other); + + template<typename OtherDerived> + Derived& operator+=(const DiagonalBase<OtherDerived>& other); + template<typename OtherDerived> + Derived& operator-=(const DiagonalBase<OtherDerived>& other); + + template<typename OtherDerived> + Derived& operator+=(const EigenBase<OtherDerived> &other); + template<typename OtherDerived> + Derived& operator-=(const EigenBase<OtherDerived> &other); Derived& operator*=(const Scalar& other); Derived& operator/=(const Scalar& other); - #define EIGEN_SPARSE_CWISE_PRODUCT_RETURN_TYPE \ - CwiseBinaryOp< \ - internal::scalar_product_op< \ - typename internal::scalar_product_traits< \ - typename internal::traits<Derived>::Scalar, \ - typename internal::traits<OtherDerived>::Scalar \ - >::ReturnType \ - >, \ - const Derived, \ - const OtherDerived \ - > + template<typename OtherDerived> struct CwiseProductDenseReturnType { + typedef CwiseBinaryOp<internal::scalar_product_op<typename ScalarBinaryOpTraits< + typename internal::traits<Derived>::Scalar, + typename internal::traits<OtherDerived>::Scalar + >::ReturnType>, + const Derived, + const OtherDerived + > Type; + }; template<typename OtherDerived> - EIGEN_STRONG_INLINE const EIGEN_SPARSE_CWISE_PRODUCT_RETURN_TYPE + EIGEN_STRONG_INLINE const typename CwiseProductDenseReturnType<OtherDerived>::Type cwiseProduct(const MatrixBase<OtherDerived> &other) const; - // sparse * sparse - template<typename OtherDerived> - const typename SparseSparseProductReturnType<Derived,OtherDerived>::Type - operator*(const SparseMatrixBase<OtherDerived> &other) const; - // sparse * diagonal template<typename OtherDerived> - const SparseDiagonalProduct<Derived,OtherDerived> - operator*(const DiagonalBase<OtherDerived> &other) const; + const Product<Derived,OtherDerived> + operator*(const DiagonalBase<OtherDerived> &other) const + { return Product<Derived,OtherDerived>(derived(), other.derived()); } // diagonal * sparse template<typename OtherDerived> friend - const SparseDiagonalProduct<OtherDerived,Derived> + const Product<OtherDerived,Derived> operator*(const DiagonalBase<OtherDerived> &lhs, const SparseMatrixBase& rhs) - { return SparseDiagonalProduct<OtherDerived,Derived>(lhs.derived(), rhs.derived()); } - - /** dense * sparse (return a dense object unless it is an outer product) */ - template<typename OtherDerived> friend - const typename DenseSparseProductReturnType<OtherDerived,Derived>::Type - operator*(const MatrixBase<OtherDerived>& lhs, const Derived& rhs) - { return typename DenseSparseProductReturnType<OtherDerived,Derived>::Type(lhs.derived(),rhs); } - - /** sparse * dense (returns a dense object unless it is an outer product) */ + { return Product<OtherDerived,Derived>(lhs.derived(), rhs.derived()); } + + // sparse * sparse + template<typename OtherDerived> + const Product<Derived,OtherDerived,AliasFreeProduct> + operator*(const SparseMatrixBase<OtherDerived> &other) const; + + // sparse * dense template<typename OtherDerived> - const typename SparseDenseProductReturnType<Derived,OtherDerived>::Type + const Product<Derived,OtherDerived> operator*(const MatrixBase<OtherDerived> &other) const - { return typename SparseDenseProductReturnType<Derived,OtherDerived>::Type(derived(), other.derived()); } + { return Product<Derived,OtherDerived>(derived(), other.derived()); } + + // dense * sparse + template<typename OtherDerived> friend + const Product<OtherDerived,Derived> + operator*(const MatrixBase<OtherDerived> &lhs, const SparseMatrixBase& rhs) + { return Product<OtherDerived,Derived>(lhs.derived(), rhs.derived()); } /** \returns an expression of P H P^-1 where H is the matrix represented by \c *this */ - SparseSymmetricPermutationProduct<Derived,Upper|Lower> twistedBy(const PermutationMatrix<Dynamic,Dynamic,Index>& perm) const + SparseSymmetricPermutationProduct<Derived,Upper|Lower> twistedBy(const PermutationMatrix<Dynamic,Dynamic,StorageIndex>& perm) const { return SparseSymmetricPermutationProduct<Derived,Upper|Lower>(derived(), perm); } @@ -370,22 +329,16 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> template<typename OtherDerived> Derived& operator*=(const SparseMatrixBase<OtherDerived>& other); - #ifdef EIGEN2_SUPPORT - // deprecated - template<typename OtherDerived> - typename internal::plain_matrix_type_column_major<OtherDerived>::type - solveTriangular(const MatrixBase<OtherDerived>& other) const; - - // deprecated - template<typename OtherDerived> - void solveTriangularInPlace(MatrixBase<OtherDerived>& other) const; - #endif // EIGEN2_SUPPORT - template<int Mode> - inline const SparseTriangularView<Derived, Mode> triangularView() const; + inline const TriangularView<const Derived, Mode> triangularView() const; + + template<unsigned int UpLo> struct SelfAdjointViewReturnType { typedef SparseSelfAdjointView<Derived, UpLo> Type; }; + template<unsigned int UpLo> struct ConstSelfAdjointViewReturnType { typedef const SparseSelfAdjointView<const Derived, UpLo> Type; }; - template<unsigned int UpLo> inline const SparseSelfAdjointView<Derived, UpLo> selfadjointView() const; - template<unsigned int UpLo> inline SparseSelfAdjointView<Derived, UpLo> selfadjointView(); + template<unsigned int UpLo> inline + typename ConstSelfAdjointViewReturnType<UpLo>::Type selfadjointView() const; + template<unsigned int UpLo> inline + typename SelfAdjointViewReturnType<UpLo>::Type selfadjointView(); template<typename OtherDerived> Scalar dot(const MatrixBase<OtherDerived>& other) const; template<typename OtherDerived> Scalar dot(const SparseMatrixBase<OtherDerived>& other) const; @@ -393,9 +346,9 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> RealScalar norm() const; RealScalar blueNorm() const; - Transpose<Derived> transpose() { return derived(); } - const Transpose<const Derived> transpose() const { return derived(); } - const AdjointReturnType adjoint() const { return transpose(); } + TransposeReturnType transpose() { return TransposeReturnType(derived()); } + const ConstTransposeReturnType transpose() const { return ConstTransposeReturnType(derived()); } + const AdjointReturnType adjoint() const { return AdjointReturnType(transpose()); } // inner-vector typedef Block<Derived,IsRowMajor?1:Dynamic,IsRowMajor?Dynamic:1,true> InnerVectorReturnType; @@ -409,25 +362,14 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> InnerVectorsReturnType innerVectors(Index outerStart, Index outerSize); const ConstInnerVectorsReturnType innerVectors(Index outerStart, Index outerSize) const; - /** \internal use operator= */ - template<typename DenseDerived> - void evalTo(MatrixBase<DenseDerived>& dst) const - { - dst.setZero(); - for (Index j=0; j<outerSize(); ++j) - for (typename Derived::InnerIterator i(derived(),j); i; ++i) - dst.coeffRef(i.row(),i.col()) = i.value(); - } - - Matrix<Scalar,RowsAtCompileTime,ColsAtCompileTime> toDense() const + DenseMatrixType toDense() const { - return derived(); + return DenseMatrixType(derived()); } template<typename OtherDerived> bool isApprox(const SparseMatrixBase<OtherDerived>& other, - const RealScalar& prec = NumTraits<Scalar>::dummy_precision()) const - { return toDense().isApprox(other.toDense(),prec); } + const RealScalar& prec = NumTraits<Scalar>::dummy_precision()) const; template<typename OtherDerived> bool isApprox(const MatrixBase<OtherDerived>& other, @@ -443,10 +385,19 @@ template<typename Derived> class SparseMatrixBase : public EigenBase<Derived> { return typename internal::eval<Derived>::type(derived()); } Scalar sum() const; + + inline const SparseView<Derived> + pruned(const Scalar& reference = Scalar(0), const RealScalar& epsilon = NumTraits<Scalar>::dummy_precision()) const; protected: bool m_isRValue; + + static inline StorageIndex convert_index(const Index idx) { + return internal::convert_index<StorageIndex>(idx); + } + private: + template<typename Dest> void evalTo(Dest &) const; }; } // end namespace Eigen |