Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-02-22 10:34:43

0001 // This file is part of Eigen, a lightweight C++ template library
0002 // for linear algebra.
0003 //
0004 // Copyright (C) 2008-2014 Gael Guennebaud <gael.guennebaud@inria.fr>
0005 //
0006 // This Source Code Form is subject to the terms of the Mozilla
0007 // Public License v. 2.0. If a copy of the MPL was not distributed
0008 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
0009 
0010 #ifndef EIGEN_SPARSEMATRIXBASE_H
0011 #define EIGEN_SPARSEMATRIXBASE_H
0012 
0013 namespace Eigen { 
0014 
0015 /** \ingroup SparseCore_Module
0016   *
0017   * \class SparseMatrixBase
0018   *
0019   * \brief Base class of any sparse matrices or sparse expressions
0020   *
0021   * \tparam Derived is the derived type, e.g. a sparse matrix type, or an expression, etc.
0022   *
0023   * This class can be extended with the help of the plugin mechanism described on the page
0024   * \ref TopicCustomizing_Plugins by defining the preprocessor symbol \c EIGEN_SPARSEMATRIXBASE_PLUGIN.
0025   */
0026 template<typename Derived> class SparseMatrixBase
0027   : public EigenBase<Derived>
0028 {
0029   public:
0030 
0031     typedef typename internal::traits<Derived>::Scalar Scalar;
0032     
0033     /** The numeric type of the expression' coefficients, e.g. float, double, int or std::complex<float>, etc.
0034       *
0035       * It is an alias for the Scalar type */
0036     typedef Scalar value_type;
0037     
0038     typedef typename internal::packet_traits<Scalar>::type PacketScalar;
0039     typedef typename internal::traits<Derived>::StorageKind StorageKind;
0040 
0041     /** The integer type used to \b store indices within a SparseMatrix.
0042       * For a \c SparseMatrix<Scalar,Options,IndexType> it an alias of the third template parameter \c IndexType. */
0043     typedef typename internal::traits<Derived>::StorageIndex StorageIndex;
0044 
0045     typedef typename internal::add_const_on_value_type_if_arithmetic<
0046                          typename internal::packet_traits<Scalar>::type
0047                      >::type PacketReturnType;
0048 
0049     typedef SparseMatrixBase StorageBaseType;
0050 
0051     typedef Matrix<StorageIndex,Dynamic,1> IndexVector;
0052     typedef Matrix<Scalar,Dynamic,1> ScalarVector;
0053     
0054     template<typename OtherDerived>
0055     Derived& operator=(const EigenBase<OtherDerived> &other);
0056 
0057     enum {
0058 
0059       RowsAtCompileTime = internal::traits<Derived>::RowsAtCompileTime,
0060         /**< The number of rows at compile-time. This is just a copy of the value provided
0061           * by the \a Derived type. If a value is not known at compile-time,
0062           * it is set to the \a Dynamic constant.
0063           * \sa MatrixBase::rows(), MatrixBase::cols(), ColsAtCompileTime, SizeAtCompileTime */
0064 
0065       ColsAtCompileTime = internal::traits<Derived>::ColsAtCompileTime,
0066         /**< The number of columns at compile-time. This is just a copy of the value provided
0067           * by the \a Derived type. If a value is not known at compile-time,
0068           * it is set to the \a Dynamic constant.
0069           * \sa MatrixBase::rows(), MatrixBase::cols(), RowsAtCompileTime, SizeAtCompileTime */
0070 
0071 
0072       SizeAtCompileTime = (internal::size_at_compile_time<internal::traits<Derived>::RowsAtCompileTime,
0073                                                    internal::traits<Derived>::ColsAtCompileTime>::ret),
0074         /**< This is equal to the number of coefficients, i.e. the number of
0075           * rows times the number of columns, or to \a Dynamic if this is not
0076           * known at compile-time. \sa RowsAtCompileTime, ColsAtCompileTime */
0077 
0078       MaxRowsAtCompileTime = RowsAtCompileTime,
0079       MaxColsAtCompileTime = ColsAtCompileTime,
0080 
0081       MaxSizeAtCompileTime = (internal::size_at_compile_time<MaxRowsAtCompileTime,
0082                                                       MaxColsAtCompileTime>::ret),
0083 
0084       IsVectorAtCompileTime = RowsAtCompileTime == 1 || ColsAtCompileTime == 1,
0085         /**< This is set to true if either the number of rows or the number of
0086           * columns is known at compile-time to be equal to 1. Indeed, in that case,
0087           * we are dealing with a column-vector (if there is only one column) or with
0088           * a row-vector (if there is only one row). */
0089 
0090       NumDimensions = int(MaxSizeAtCompileTime) == 1 ? 0 : bool(IsVectorAtCompileTime) ? 1 : 2,
0091         /**< This value is equal to Tensor::NumDimensions, i.e. 0 for scalars, 1 for vectors,
0092          * and 2 for matrices.
0093          */
0094 
0095       Flags = internal::traits<Derived>::Flags,
0096         /**< This stores expression \ref flags flags which may or may not be inherited by new expressions
0097           * constructed from this one. See the \ref flags "list of flags".
0098           */
0099 
0100       IsRowMajor = Flags&RowMajorBit ? 1 : 0,
0101       
0102       InnerSizeAtCompileTime = int(IsVectorAtCompileTime) ? int(SizeAtCompileTime)
0103                              : int(IsRowMajor) ? int(ColsAtCompileTime) : int(RowsAtCompileTime),
0104 
0105       #ifndef EIGEN_PARSED_BY_DOXYGEN
0106       _HasDirectAccess = (int(Flags)&DirectAccessBit) ? 1 : 0 // workaround sunCC
0107       #endif
0108     };
0109 
0110     /** \internal the return type of MatrixBase::adjoint() */
0111     typedef typename internal::conditional<NumTraits<Scalar>::IsComplex,
0112                         CwiseUnaryOp<internal::scalar_conjugate_op<Scalar>, Eigen::Transpose<const Derived> >,
0113                         Transpose<const Derived>
0114                      >::type AdjointReturnType;
0115     typedef Transpose<Derived> TransposeReturnType;
0116     typedef typename internal::add_const<Transpose<const Derived> >::type ConstTransposeReturnType;
0117 
0118     // FIXME storage order do not match evaluator storage order
0119     typedef SparseMatrix<Scalar, Flags&RowMajorBit ? RowMajor : ColMajor, StorageIndex> PlainObject;
0120 
0121 #ifndef EIGEN_PARSED_BY_DOXYGEN
0122     /** This is the "real scalar" type; if the \a Scalar type is already real numbers
0123       * (e.g. int, float or double) then \a RealScalar is just the same as \a Scalar. If
0124       * \a Scalar is \a std::complex<T> then RealScalar is \a T.
0125       *
0126       * \sa class NumTraits
0127       */
0128     typedef typename NumTraits<Scalar>::Real RealScalar;
0129 
0130     /** \internal the return type of coeff()
0131       */
0132     typedef typename internal::conditional<_HasDirectAccess, const Scalar&, Scalar>::type CoeffReturnType;
0133 
0134     /** \internal Represents a matrix with all coefficients equal to one another*/
0135     typedef CwiseNullaryOp<internal::scalar_constant_op<Scalar>,Matrix<Scalar,Dynamic,Dynamic> > ConstantReturnType;
0136 
0137     /** type of the equivalent dense matrix */
0138     typedef Matrix<Scalar,RowsAtCompileTime,ColsAtCompileTime> DenseMatrixType;
0139     /** type of the equivalent square matrix */
0140     typedef Matrix<Scalar,EIGEN_SIZE_MAX(RowsAtCompileTime,ColsAtCompileTime),
0141                           EIGEN_SIZE_MAX(RowsAtCompileTime,ColsAtCompileTime)> SquareMatrixType;
0142 
0143     inline const Derived& derived() const { return *static_cast<const Derived*>(this); }
0144     inline Derived& derived() { return *static_cast<Derived*>(this); }
0145     inline Derived& const_cast_derived() const
0146     { return *static_cast<Derived*>(const_cast<SparseMatrixBase*>(this)); }
0147 
0148     typedef EigenBase<Derived> Base;
0149 
0150 #endif // not EIGEN_PARSED_BY_DOXYGEN
0151 
0152 #define EIGEN_CURRENT_STORAGE_BASE_CLASS Eigen::SparseMatrixBase
0153 #ifdef EIGEN_PARSED_BY_DOXYGEN
0154 #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> */
0155 #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> */
0156 #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> */
0157 #else
0158 #define EIGEN_DOC_UNARY_ADDONS(X,Y)
0159 #define EIGEN_DOC_BLOCK_ADDONS_NOT_INNER_PANEL
0160 #define EIGEN_DOC_BLOCK_ADDONS_INNER_PANEL_IF(COND)
0161 #endif
0162 #   include "../plugins/CommonCwiseUnaryOps.h"
0163 #   include "../plugins/CommonCwiseBinaryOps.h"
0164 #   include "../plugins/MatrixCwiseUnaryOps.h"
0165 #   include "../plugins/MatrixCwiseBinaryOps.h"
0166 #   include "../plugins/BlockMethods.h"
0167 #   ifdef EIGEN_SPARSEMATRIXBASE_PLUGIN
0168 #     include EIGEN_SPARSEMATRIXBASE_PLUGIN
0169 #   endif
0170 #undef EIGEN_CURRENT_STORAGE_BASE_CLASS
0171 #undef EIGEN_DOC_UNARY_ADDONS
0172 #undef EIGEN_DOC_BLOCK_ADDONS_NOT_INNER_PANEL
0173 #undef EIGEN_DOC_BLOCK_ADDONS_INNER_PANEL_IF
0174 
0175     /** \returns the number of rows. \sa cols() */
0176     inline Index rows() const { return derived().rows(); }
0177     /** \returns the number of columns. \sa rows() */
0178     inline Index cols() const { return derived().cols(); }
0179     /** \returns the number of coefficients, which is \a rows()*cols().
0180       * \sa rows(), cols(). */
0181     inline Index size() const { return rows() * cols(); }
0182     /** \returns true if either the number of rows or the number of columns is equal to 1.
0183       * In other words, this function returns
0184       * \code rows()==1 || cols()==1 \endcode
0185       * \sa rows(), cols(), IsVectorAtCompileTime. */
0186     inline bool isVector() const { return rows()==1 || cols()==1; }
0187     /** \returns the size of the storage major dimension,
0188       * i.e., the number of columns for a columns major matrix, and the number of rows otherwise */
0189     Index outerSize() const { return (int(Flags)&RowMajorBit) ? this->rows() : this->cols(); }
0190     /** \returns the size of the inner dimension according to the storage order,
0191       * i.e., the number of rows for a columns major matrix, and the number of cols otherwise */
0192     Index innerSize() const { return (int(Flags)&RowMajorBit) ? this->cols() : this->rows(); }
0193 
0194     bool isRValue() const { return m_isRValue; }
0195     Derived& markAsRValue() { m_isRValue = true; return derived(); }
0196 
0197     SparseMatrixBase() : m_isRValue(false) { /* TODO check flags */ }
0198 
0199     
0200     template<typename OtherDerived>
0201     Derived& operator=(const ReturnByValue<OtherDerived>& other);
0202 
0203     template<typename OtherDerived>
0204     inline Derived& operator=(const SparseMatrixBase<OtherDerived>& other);
0205 
0206     inline Derived& operator=(const Derived& other);
0207 
0208   protected:
0209 
0210     template<typename OtherDerived>
0211     inline Derived& assign(const OtherDerived& other);
0212 
0213     template<typename OtherDerived>
0214     inline void assignGeneric(const OtherDerived& other);
0215 
0216   public:
0217 
0218     friend std::ostream & operator << (std::ostream & s, const SparseMatrixBase& m)
0219     {
0220       typedef typename Derived::Nested Nested;
0221       typedef typename internal::remove_all<Nested>::type NestedCleaned;
0222 
0223       if (Flags&RowMajorBit)
0224       {
0225         Nested nm(m.derived());
0226         internal::evaluator<NestedCleaned> thisEval(nm);
0227         for (Index row=0; row<nm.outerSize(); ++row)
0228         {
0229           Index col = 0;
0230           for (typename internal::evaluator<NestedCleaned>::InnerIterator it(thisEval, row); it; ++it)
0231           {
0232             for ( ; col<it.index(); ++col)
0233               s << "0 ";
0234             s << it.value() << " ";
0235             ++col;
0236           }
0237           for ( ; col<m.cols(); ++col)
0238             s << "0 ";
0239           s << std::endl;
0240         }
0241       }
0242       else
0243       {
0244         Nested nm(m.derived());
0245         internal::evaluator<NestedCleaned> thisEval(nm);
0246         if (m.cols() == 1) {
0247           Index row = 0;
0248           for (typename internal::evaluator<NestedCleaned>::InnerIterator it(thisEval, 0); it; ++it)
0249           {
0250             for ( ; row<it.index(); ++row)
0251               s << "0" << std::endl;
0252             s << it.value() << std::endl;
0253             ++row;
0254           }
0255           for ( ; row<m.rows(); ++row)
0256             s << "0" << std::endl;
0257         }
0258         else
0259         {
0260           SparseMatrix<Scalar, RowMajorBit, StorageIndex> trans = m;
0261           s << static_cast<const SparseMatrixBase<SparseMatrix<Scalar, RowMajorBit, StorageIndex> >&>(trans);
0262         }
0263       }
0264       return s;
0265     }
0266 
0267     template<typename OtherDerived>
0268     Derived& operator+=(const SparseMatrixBase<OtherDerived>& other);
0269     template<typename OtherDerived>
0270     Derived& operator-=(const SparseMatrixBase<OtherDerived>& other);
0271     
0272     template<typename OtherDerived>
0273     Derived& operator+=(const DiagonalBase<OtherDerived>& other);
0274     template<typename OtherDerived>
0275     Derived& operator-=(const DiagonalBase<OtherDerived>& other);
0276 
0277     template<typename OtherDerived>
0278     Derived& operator+=(const EigenBase<OtherDerived> &other);
0279     template<typename OtherDerived>
0280     Derived& operator-=(const EigenBase<OtherDerived> &other);
0281 
0282     Derived& operator*=(const Scalar& other);
0283     Derived& operator/=(const Scalar& other);
0284 
0285     template<typename OtherDerived> struct CwiseProductDenseReturnType {
0286       typedef CwiseBinaryOp<internal::scalar_product_op<typename ScalarBinaryOpTraits<
0287                                                           typename internal::traits<Derived>::Scalar,
0288                                                           typename internal::traits<OtherDerived>::Scalar
0289                                                         >::ReturnType>,
0290                             const Derived,
0291                             const OtherDerived
0292                           > Type;
0293     };
0294 
0295     template<typename OtherDerived>
0296     EIGEN_STRONG_INLINE const typename CwiseProductDenseReturnType<OtherDerived>::Type
0297     cwiseProduct(const MatrixBase<OtherDerived> &other) const;
0298 
0299     // sparse * diagonal
0300     template<typename OtherDerived>
0301     const Product<Derived,OtherDerived>
0302     operator*(const DiagonalBase<OtherDerived> &other) const
0303     { return Product<Derived,OtherDerived>(derived(), other.derived()); }
0304 
0305     // diagonal * sparse
0306     template<typename OtherDerived> friend
0307     const Product<OtherDerived,Derived>
0308     operator*(const DiagonalBase<OtherDerived> &lhs, const SparseMatrixBase& rhs)
0309     { return Product<OtherDerived,Derived>(lhs.derived(), rhs.derived()); }
0310     
0311     // sparse * sparse
0312     template<typename OtherDerived>
0313     const Product<Derived,OtherDerived,AliasFreeProduct>
0314     operator*(const SparseMatrixBase<OtherDerived> &other) const;
0315     
0316     // sparse * dense
0317     template<typename OtherDerived>
0318     const Product<Derived,OtherDerived>
0319     operator*(const MatrixBase<OtherDerived> &other) const
0320     { return Product<Derived,OtherDerived>(derived(), other.derived()); }
0321     
0322     // dense * sparse
0323     template<typename OtherDerived> friend
0324     const Product<OtherDerived,Derived>
0325     operator*(const MatrixBase<OtherDerived> &lhs, const SparseMatrixBase& rhs)
0326     { return Product<OtherDerived,Derived>(lhs.derived(), rhs.derived()); }
0327     
0328      /** \returns an expression of P H P^-1 where H is the matrix represented by \c *this */
0329     SparseSymmetricPermutationProduct<Derived,Upper|Lower> twistedBy(const PermutationMatrix<Dynamic,Dynamic,StorageIndex>& perm) const
0330     {
0331       return SparseSymmetricPermutationProduct<Derived,Upper|Lower>(derived(), perm);
0332     }
0333 
0334     template<typename OtherDerived>
0335     Derived& operator*=(const SparseMatrixBase<OtherDerived>& other);
0336 
0337     template<int Mode>
0338     inline const TriangularView<const Derived, Mode> triangularView() const;
0339     
0340     template<unsigned int UpLo> struct SelfAdjointViewReturnType { typedef SparseSelfAdjointView<Derived, UpLo> Type; };
0341     template<unsigned int UpLo> struct ConstSelfAdjointViewReturnType { typedef const SparseSelfAdjointView<const Derived, UpLo> Type; };
0342 
0343     template<unsigned int UpLo> inline 
0344     typename ConstSelfAdjointViewReturnType<UpLo>::Type selfadjointView() const;
0345     template<unsigned int UpLo> inline
0346     typename SelfAdjointViewReturnType<UpLo>::Type selfadjointView();
0347 
0348     template<typename OtherDerived> Scalar dot(const MatrixBase<OtherDerived>& other) const;
0349     template<typename OtherDerived> Scalar dot(const SparseMatrixBase<OtherDerived>& other) const;
0350     RealScalar squaredNorm() const;
0351     RealScalar norm()  const;
0352     RealScalar blueNorm() const;
0353 
0354     TransposeReturnType transpose() { return TransposeReturnType(derived()); }
0355     const ConstTransposeReturnType transpose() const { return ConstTransposeReturnType(derived()); }
0356     const AdjointReturnType adjoint() const { return AdjointReturnType(transpose()); }
0357 
0358     DenseMatrixType toDense() const
0359     {
0360       return DenseMatrixType(derived());
0361     }
0362 
0363     template<typename OtherDerived>
0364     bool isApprox(const SparseMatrixBase<OtherDerived>& other,
0365                   const RealScalar& prec = NumTraits<Scalar>::dummy_precision()) const;
0366 
0367     template<typename OtherDerived>
0368     bool isApprox(const MatrixBase<OtherDerived>& other,
0369                   const RealScalar& prec = NumTraits<Scalar>::dummy_precision()) const
0370     { return toDense().isApprox(other,prec); }
0371 
0372     /** \returns the matrix or vector obtained by evaluating this expression.
0373       *
0374       * Notice that in the case of a plain matrix or vector (not an expression) this function just returns
0375       * a const reference, in order to avoid a useless copy.
0376       */
0377     inline const typename internal::eval<Derived>::type eval() const
0378     { return typename internal::eval<Derived>::type(derived()); }
0379 
0380     Scalar sum() const;
0381     
0382     inline const SparseView<Derived>
0383     pruned(const Scalar& reference = Scalar(0), const RealScalar& epsilon = NumTraits<Scalar>::dummy_precision()) const;
0384 
0385   protected:
0386 
0387     bool m_isRValue;
0388 
0389     static inline StorageIndex convert_index(const Index idx) {
0390       return internal::convert_index<StorageIndex>(idx);
0391     }
0392   private:
0393     template<typename Dest> void evalTo(Dest &) const;
0394 };
0395 
0396 } // end namespace Eigen
0397 
0398 #endif // EIGEN_SPARSEMATRIXBASE_H