dune-istl
2.2.0
|
00001 #ifndef DUNE_MATRIX_INDEX_SET_HH 00002 #define DUNE_MATRIX_INDEX_SET_HH 00003 00004 #include <vector> 00005 #include <set> 00006 00007 namespace Dune { 00008 00009 00011 class MatrixIndexSet 00012 { 00013 00014 public: 00015 typedef std::size_t size_type; 00016 00018 MatrixIndexSet() : rows_(0), cols_(0) 00019 {} 00020 00022 MatrixIndexSet(size_type rows, size_type cols) : rows_(rows), cols_(cols) { 00023 indices_.resize(rows_); 00024 } 00025 00027 void resize(size_type rows, size_type cols) { 00028 rows_ = rows; 00029 cols_ = cols; 00030 indices_.resize(rows_); 00031 } 00032 00034 void add(size_type i, size_type j) { 00035 indices_[i].insert(j); 00036 } 00037 00039 size_type size() const { 00040 size_type entries = 0; 00041 for (size_type i=0; i<rows_; i++) 00042 entries += indices_[i].size(); 00043 00044 return entries; 00045 } 00046 00048 size_type rows() const {return rows_;} 00049 00050 00052 size_type rowsize(size_type row) const {return indices_[row].size();} 00053 00060 template <class MatrixType> 00061 void import(const MatrixType& m, size_type rowOffset=0, size_type colOffset=0) { 00062 00063 typedef typename MatrixType::row_type RowType; 00064 typedef typename RowType::ConstIterator ColumnIterator; 00065 00066 for (size_type rowIdx=0; rowIdx<m.N(); rowIdx++) { 00067 00068 const RowType& row = m[rowIdx]; 00069 00070 ColumnIterator cIt = row.begin(); 00071 ColumnIterator cEndIt = row.end(); 00072 00073 for(; cIt!=cEndIt; ++cIt) 00074 add(rowIdx+rowOffset, cIt.index()+colOffset); 00075 00076 } 00077 00078 } 00079 00085 template <class MatrixType> 00086 void exportIdx(MatrixType& matrix) const { 00087 00088 matrix.setSize(rows_, cols_); 00089 matrix.setBuildMode(MatrixType::random); 00090 00091 for (size_type i=0; i<rows_; i++) 00092 matrix.setrowsize(i, indices_[i].size()); 00093 00094 matrix.endrowsizes(); 00095 00096 for (size_type i=0; i<rows_; i++) { 00097 00098 typename std::set<size_type>::iterator it = indices_[i].begin(); 00099 for (; it!=indices_[i].end(); ++it) 00100 matrix.addindex(i, *it); 00101 00102 } 00103 00104 matrix.endindices(); 00105 00106 } 00107 00108 private: 00109 00110 std::vector<std::set<size_type> > indices_; 00111 00112 size_type rows_, cols_; 00113 00114 }; 00115 00116 00117 } // end namespace Dune 00118 00119 #endif