$treeview $search $mathjax
Eigen
3.2.5
$projectbrief
|
$projectbrief
|
$searchbox |
00001 00002 // This file is part of Eigen, a lightweight C++ template library 00003 // for linear algebra. 00004 // 00005 // Copyright (C) 2012 Désiré Nuentsa-Wakam <desire.nuentsa_wakam@inria.fr> 00006 // 00007 // This Source Code Form is subject to the terms of the Mozilla 00008 // Public License v. 2.0. If a copy of the MPL was not distributed 00009 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 00010 00011 #ifndef EIGEN_ORDERING_H 00012 #define EIGEN_ORDERING_H 00013 00014 namespace Eigen { 00015 00016 #include "Eigen_Colamd.h" 00017 00018 namespace internal { 00019 00025 template<typename MatrixType> 00026 void ordering_helper_at_plus_a(const MatrixType& mat, MatrixType& symmat) 00027 { 00028 MatrixType C; 00029 C = mat.transpose(); // NOTE: Could be costly 00030 for (int i = 0; i < C.rows(); i++) 00031 { 00032 for (typename MatrixType::InnerIterator it(C, i); it; ++it) 00033 it.valueRef() = 0.0; 00034 } 00035 symmat = C + mat; 00036 } 00037 00038 } 00039 00040 #ifndef EIGEN_MPL2_ONLY 00041 00050 template <typename Index> 00051 class AMDOrdering 00052 { 00053 public: 00054 typedef PermutationMatrix<Dynamic, Dynamic, Index> PermutationType; 00055 00059 template <typename MatrixType> 00060 void operator()(const MatrixType& mat, PermutationType& perm) 00061 { 00062 // Compute the symmetric pattern 00063 SparseMatrix<typename MatrixType::Scalar, ColMajor, Index> symm; 00064 internal::ordering_helper_at_plus_a(mat,symm); 00065 00066 // Call the AMD routine 00067 //m_mat.prune(keep_diag()); 00068 internal::minimum_degree_ordering(symm, perm); 00069 } 00070 00072 template <typename SrcType, unsigned int SrcUpLo> 00073 void operator()(const SparseSelfAdjointView<SrcType, SrcUpLo>& mat, PermutationType& perm) 00074 { 00075 SparseMatrix<typename SrcType::Scalar, ColMajor, Index> C; C = mat; 00076 00077 // Call the AMD routine 00078 // m_mat.prune(keep_diag()); //Remove the diagonal elements 00079 internal::minimum_degree_ordering(C, perm); 00080 } 00081 }; 00082 00083 #endif // EIGEN_MPL2_ONLY 00084 00093 template <typename Index> 00094 class NaturalOrdering 00095 { 00096 public: 00097 typedef PermutationMatrix<Dynamic, Dynamic, Index> PermutationType; 00098 00100 template <typename MatrixType> 00101 void operator()(const MatrixType& /*mat*/, PermutationType& perm) 00102 { 00103 perm.resize(0); 00104 } 00105 00106 }; 00107 00114 template<typename Index> 00115 class COLAMDOrdering 00116 { 00117 public: 00118 typedef PermutationMatrix<Dynamic, Dynamic, Index> PermutationType; 00119 typedef Matrix<Index, Dynamic, 1> IndexVector; 00120 00124 template <typename MatrixType> 00125 void operator() (const MatrixType& mat, PermutationType& perm) 00126 { 00127 eigen_assert(mat.isCompressed() && "COLAMDOrdering requires a sparse matrix in compressed mode. Call .makeCompressed() before passing it to COLAMDOrdering"); 00128 00129 Index m = mat.rows(); 00130 Index n = mat.cols(); 00131 Index nnz = mat.nonZeros(); 00132 // Get the recommended value of Alen to be used by colamd 00133 Index Alen = internal::colamd_recommended(nnz, m, n); 00134 // Set the default parameters 00135 double knobs [COLAMD_KNOBS]; 00136 Index stats [COLAMD_STATS]; 00137 internal::colamd_set_defaults(knobs); 00138 00139 IndexVector p(n+1), A(Alen); 00140 for(Index i=0; i <= n; i++) p(i) = mat.outerIndexPtr()[i]; 00141 for(Index i=0; i < nnz; i++) A(i) = mat.innerIndexPtr()[i]; 00142 // Call Colamd routine to compute the ordering 00143 Index info = internal::colamd(m, n, Alen, A.data(), p.data(), knobs, stats); 00144 EIGEN_UNUSED_VARIABLE(info); 00145 eigen_assert( info && "COLAMD failed " ); 00146 00147 perm.resize(n); 00148 for (Index i = 0; i < n; i++) perm.indices()(p(i)) = i; 00149 } 00150 }; 00151 00152 } // end namespace Eigen 00153 00154 #endif