$treeview $search $mathjax
Eigen
3.2.5
$projectbrief
|
$projectbrief
|
$searchbox |
00001 // This file is part of Eigen, a lightweight C++ template library 00002 // for linear algebra. 00003 // 00004 // Copyright (C) 2012 Gael Guennebaud <gael.guennebaud@inria.fr> 00005 // 00006 // This Source Code Form is subject to the terms of the Mozilla 00007 // Public License v. 2.0. If a copy of the MPL was not distributed 00008 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 00009 00010 #ifndef EIGEN_SPARSE_PERMUTATION_H 00011 #define EIGEN_SPARSE_PERMUTATION_H 00012 00013 // This file implements sparse * permutation products 00014 00015 namespace Eigen { 00016 00017 namespace internal { 00018 00019 template<typename PermutationType, typename MatrixType, int Side, bool Transposed> 00020 struct traits<permut_sparsematrix_product_retval<PermutationType, MatrixType, Side, Transposed> > 00021 { 00022 typedef typename remove_all<typename MatrixType::Nested>::type MatrixTypeNestedCleaned; 00023 typedef typename MatrixTypeNestedCleaned::Scalar Scalar; 00024 typedef typename MatrixTypeNestedCleaned::Index Index; 00025 enum { 00026 SrcStorageOrder = MatrixTypeNestedCleaned::Flags&RowMajorBit ? RowMajor : ColMajor, 00027 MoveOuter = SrcStorageOrder==RowMajor ? Side==OnTheLeft : Side==OnTheRight 00028 }; 00029 00030 typedef typename internal::conditional<MoveOuter, 00031 SparseMatrix<Scalar,SrcStorageOrder,Index>, 00032 SparseMatrix<Scalar,int(SrcStorageOrder)==RowMajor?ColMajor:RowMajor,Index> >::type ReturnType; 00033 }; 00034 00035 template<typename PermutationType, typename MatrixType, int Side, bool Transposed> 00036 struct permut_sparsematrix_product_retval 00037 : public ReturnByValue<permut_sparsematrix_product_retval<PermutationType, MatrixType, Side, Transposed> > 00038 { 00039 typedef typename remove_all<typename MatrixType::Nested>::type MatrixTypeNestedCleaned; 00040 typedef typename MatrixTypeNestedCleaned::Scalar Scalar; 00041 typedef typename MatrixTypeNestedCleaned::Index Index; 00042 00043 enum { 00044 SrcStorageOrder = MatrixTypeNestedCleaned::Flags&RowMajorBit ? RowMajor : ColMajor, 00045 MoveOuter = SrcStorageOrder==RowMajor ? Side==OnTheLeft : Side==OnTheRight 00046 }; 00047 00048 permut_sparsematrix_product_retval(const PermutationType& perm, const MatrixType& matrix) 00049 : m_permutation(perm), m_matrix(matrix) 00050 {} 00051 00052 inline int rows() const { return m_matrix.rows(); } 00053 inline int cols() const { return m_matrix.cols(); } 00054 00055 template<typename Dest> inline void evalTo(Dest& dst) const 00056 { 00057 if(MoveOuter) 00058 { 00059 SparseMatrix<Scalar,SrcStorageOrder,Index> tmp(m_matrix.rows(), m_matrix.cols()); 00060 Matrix<Index,Dynamic,1> sizes(m_matrix.outerSize()); 00061 for(Index j=0; j<m_matrix.outerSize(); ++j) 00062 { 00063 Index jp = m_permutation.indices().coeff(j); 00064 sizes[((Side==OnTheLeft) ^ Transposed) ? jp : j] = m_matrix.innerVector(((Side==OnTheRight) ^ Transposed) ? jp : j).nonZeros(); 00065 } 00066 tmp.reserve(sizes); 00067 for(Index j=0; j<m_matrix.outerSize(); ++j) 00068 { 00069 Index jp = m_permutation.indices().coeff(j); 00070 Index jsrc = ((Side==OnTheRight) ^ Transposed) ? jp : j; 00071 Index jdst = ((Side==OnTheLeft) ^ Transposed) ? jp : j; 00072 for(typename MatrixTypeNestedCleaned::InnerIterator it(m_matrix,jsrc); it; ++it) 00073 tmp.insertByOuterInner(jdst,it.index()) = it.value(); 00074 } 00075 dst = tmp; 00076 } 00077 else 00078 { 00079 SparseMatrix<Scalar,int(SrcStorageOrder)==RowMajor?ColMajor:RowMajor,Index> tmp(m_matrix.rows(), m_matrix.cols()); 00080 Matrix<Index,Dynamic,1> sizes(tmp.outerSize()); 00081 sizes.setZero(); 00082 PermutationMatrix<Dynamic,Dynamic,Index> perm; 00083 if((Side==OnTheLeft) ^ Transposed) 00084 perm = m_permutation; 00085 else 00086 perm = m_permutation.transpose(); 00087 00088 for(Index j=0; j<m_matrix.outerSize(); ++j) 00089 for(typename MatrixTypeNestedCleaned::InnerIterator it(m_matrix,j); it; ++it) 00090 sizes[perm.indices().coeff(it.index())]++; 00091 tmp.reserve(sizes); 00092 for(Index j=0; j<m_matrix.outerSize(); ++j) 00093 for(typename MatrixTypeNestedCleaned::InnerIterator it(m_matrix,j); it; ++it) 00094 tmp.insertByOuterInner(perm.indices().coeff(it.index()),j) = it.value(); 00095 dst = tmp; 00096 } 00097 } 00098 00099 protected: 00100 const PermutationType& m_permutation; 00101 typename MatrixType::Nested m_matrix; 00102 }; 00103 00104 } 00105 00106 00107 00110 template<typename SparseDerived, typename PermDerived> 00111 inline const internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheRight, false> 00112 operator*(const SparseMatrixBase<SparseDerived>& matrix, const PermutationBase<PermDerived>& perm) 00113 { 00114 return internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheRight, false>(perm, matrix.derived()); 00115 } 00116 00119 template<typename SparseDerived, typename PermDerived> 00120 inline const internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheLeft, false> 00121 operator*( const PermutationBase<PermDerived>& perm, const SparseMatrixBase<SparseDerived>& matrix) 00122 { 00123 return internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheLeft, false>(perm, matrix.derived()); 00124 } 00125 00126 00127 00130 template<typename SparseDerived, typename PermDerived> 00131 inline const internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheRight, true> 00132 operator*(const SparseMatrixBase<SparseDerived>& matrix, const Transpose<PermutationBase<PermDerived> >& tperm) 00133 { 00134 return internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheRight, true>(tperm.nestedPermutation(), matrix.derived()); 00135 } 00136 00139 template<typename SparseDerived, typename PermDerived> 00140 inline const internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheLeft, true> 00141 operator*(const Transpose<PermutationBase<PermDerived> >& tperm, const SparseMatrixBase<SparseDerived>& matrix) 00142 { 00143 return internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheLeft, true>(tperm.nestedPermutation(), matrix.derived()); 00144 } 00145 00146 } // end namespace Eigen 00147 00148 #endif // EIGEN_SPARSE_SELFADJOINTVIEW_H