00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025 #ifndef EIGEN_SPARSE_PERMUTATION_H
00026 #define EIGEN_SPARSE_PERMUTATION_H
00027
00028
00029
00030 namespace Eigen {
00031
00032 namespace internal {
00033
00034 template<typename PermutationType, typename MatrixType, int Side, bool Transposed>
00035 struct traits<permut_sparsematrix_product_retval<PermutationType, MatrixType, Side, Transposed> >
00036 {
00037 typedef typename remove_all<typename MatrixType::Nested>::type MatrixTypeNestedCleaned;
00038 typedef typename MatrixTypeNestedCleaned::Scalar Scalar;
00039 typedef typename MatrixTypeNestedCleaned::Index Index;
00040 enum {
00041 SrcStorageOrder = MatrixTypeNestedCleaned::Flags&RowMajorBit ? RowMajor : ColMajor,
00042 MoveOuter = SrcStorageOrder==RowMajor ? Side==OnTheLeft : Side==OnTheRight
00043 };
00044
00045 typedef typename internal::conditional<MoveOuter,
00046 SparseMatrix<Scalar,SrcStorageOrder,Index>,
00047 SparseMatrix<Scalar,int(SrcStorageOrder)==RowMajor?ColMajor:RowMajor,Index> >::type ReturnType;
00048 };
00049
00050 template<typename PermutationType, typename MatrixType, int Side, bool Transposed>
00051 struct permut_sparsematrix_product_retval
00052 : public ReturnByValue<permut_sparsematrix_product_retval<PermutationType, MatrixType, Side, Transposed> >
00053 {
00054 typedef typename remove_all<typename MatrixType::Nested>::type MatrixTypeNestedCleaned;
00055 typedef typename MatrixTypeNestedCleaned::Scalar Scalar;
00056 typedef typename MatrixTypeNestedCleaned::Index Index;
00057
00058 enum {
00059 SrcStorageOrder = MatrixTypeNestedCleaned::Flags&RowMajorBit ? RowMajor : ColMajor,
00060 MoveOuter = SrcStorageOrder==RowMajor ? Side==OnTheLeft : Side==OnTheRight
00061 };
00062
00063 permut_sparsematrix_product_retval(const PermutationType& perm, const MatrixType& matrix)
00064 : m_permutation(perm), m_matrix(matrix)
00065 {}
00066
00067 inline int rows() const { return m_matrix.rows(); }
00068 inline int cols() const { return m_matrix.cols(); }
00069
00070 template<typename Dest> inline void evalTo(Dest& dst) const
00071 {
00072 if(MoveOuter)
00073 {
00074 SparseMatrix<Scalar,SrcStorageOrder,Index> tmp(m_matrix.rows(), m_matrix.cols());
00075 VectorXi sizes(m_matrix.outerSize());
00076 for(Index j=0; j<m_matrix.outerSize(); ++j)
00077 {
00078 Index jp = m_permutation.indices().coeff(j);
00079 sizes[((Side==OnTheLeft) ^ Transposed) ? jp : j] = m_matrix.innerVector(((Side==OnTheRight) ^ Transposed) ? jp : j).size();
00080 }
00081 tmp.reserve(sizes);
00082 for(Index j=0; j<m_matrix.outerSize(); ++j)
00083 {
00084 Index jp = m_permutation.indices().coeff(j);
00085 Index jsrc = ((Side==OnTheRight) ^ Transposed) ? jp : j;
00086 Index jdst = ((Side==OnTheLeft) ^ Transposed) ? jp : j;
00087 for(typename MatrixTypeNestedCleaned::InnerIterator it(m_matrix,jsrc); it; ++it)
00088 tmp.insertByOuterInner(jdst,it.index()) = it.value();
00089 }
00090 dst = tmp;
00091 }
00092 else
00093 {
00094 SparseMatrix<Scalar,int(SrcStorageOrder)==RowMajor?ColMajor:RowMajor,Index> tmp(m_matrix.rows(), m_matrix.cols());
00095 VectorXi sizes(tmp.outerSize());
00096 sizes.setZero();
00097 PermutationMatrix<Dynamic,Dynamic,Index> perm;
00098 if((Side==OnTheLeft) ^ Transposed)
00099 perm = m_permutation;
00100 else
00101 perm = m_permutation.transpose();
00102
00103 for(Index j=0; j<m_matrix.outerSize(); ++j)
00104 for(typename MatrixTypeNestedCleaned::InnerIterator it(m_matrix,j); it; ++it)
00105 sizes[perm.indices().coeff(it.index())]++;
00106 tmp.reserve(sizes);
00107 for(Index j=0; j<m_matrix.outerSize(); ++j)
00108 for(typename MatrixTypeNestedCleaned::InnerIterator it(m_matrix,j); it; ++it)
00109 tmp.insertByOuterInner(perm.indices().coeff(it.index()),j) = it.value();
00110 dst = tmp;
00111 }
00112 }
00113
00114 protected:
00115 const PermutationType& m_permutation;
00116 typename MatrixType::Nested m_matrix;
00117 };
00118
00119 }
00120
00121
00122
00125 template<typename SparseDerived, typename PermDerived>
00126 inline const internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheRight, false>
00127 operator*(const SparseMatrixBase<SparseDerived>& matrix, const PermutationBase<PermDerived>& perm)
00128 {
00129 return internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheRight, false>(perm, matrix.derived());
00130 }
00131
00134 template<typename SparseDerived, typename PermDerived>
00135 inline const internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheLeft, false>
00136 operator*( const PermutationBase<PermDerived>& perm, const SparseMatrixBase<SparseDerived>& matrix)
00137 {
00138 return internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheLeft, false>(perm, matrix.derived());
00139 }
00140
00141
00142
00145 template<typename SparseDerived, typename PermDerived>
00146 inline const internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheRight, true>
00147 operator*(const SparseMatrixBase<SparseDerived>& matrix, const Transpose<PermutationBase<PermDerived> >& tperm)
00148 {
00149 return internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheRight, true>(tperm.nestedPermutation(), matrix.derived());
00150 }
00151
00154 template<typename SparseDerived, typename PermDerived>
00155 inline const internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheLeft, true>
00156 operator*(const Transpose<PermutationBase<PermDerived> >& tperm, const SparseMatrixBase<SparseDerived>& matrix)
00157 {
00158 return internal::permut_sparsematrix_product_retval<PermutationBase<PermDerived>, SparseDerived, OnTheLeft, true>(tperm.nestedPermutation(), matrix.derived());
00159 }
00160
00161 }
00162
00163 #endif // EIGEN_SPARSE_SELFADJOINTVIEW_H