dune-istl  2.1.1
matrixindexset.hh
Go to the documentation of this file.
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 
00017     MatrixIndexSet() : rows_(0), cols_(0)
00018     {}
00019 
00021     MatrixIndexSet(int rows, int cols) : rows_(rows), cols_(cols) {
00022         indices_.resize(rows_);
00023     }
00024 
00026     void resize(int rows, int cols) {
00027         rows_ = rows;
00028         cols_ = cols;
00029         indices_.resize(rows_);
00030     }
00031 
00033     void add(int i, int j) {
00034         indices_[i].insert(j);
00035     }
00036 
00038     int size() const {
00039         int entries = 0;
00040         for (int i=0; i<rows_; i++)
00041             entries += indices_[i].size();
00042 
00043         return entries;
00044     }
00045 
00047     int rows() const {return rows_;}
00048 
00049 
00051     int rowsize(int row) const {return indices_[row].size();}
00052 
00059     template <class MatrixType>
00060     void import(const MatrixType& m, int rowOffset=0, int colOffset=0) {
00061 
00062         typedef typename MatrixType::row_type RowType;
00063         typedef typename RowType::ConstIterator ColumnIterator;
00064 
00065         for (int rowIdx=0; rowIdx<m.N(); rowIdx++) {
00066 
00067             const RowType& row = m[rowIdx];
00068             
00069             ColumnIterator cIt    = row.begin();
00070             ColumnIterator cEndIt = row.end();
00071             
00072             for(; cIt!=cEndIt; ++cIt)
00073                 add(rowIdx+rowOffset, cIt.index()+colOffset);
00074 
00075         }
00076 
00077     }
00078 
00084     template <class MatrixType>
00085     void exportIdx(MatrixType& matrix) const {
00086 
00087         matrix.setSize(rows_, cols_);
00088         matrix.setBuildMode(MatrixType::random);
00089     
00090         for (int i=0; i<rows_; i++)
00091             matrix.setrowsize(i, indices_[i].size());
00092     
00093         matrix.endrowsizes();
00094     
00095         for (int i=0; i<rows_; i++) {
00096 
00097             typename std::set<unsigned int>::iterator it = indices_[i].begin();
00098             for (; it!=indices_[i].end(); ++it)
00099                 matrix.addindex(i, *it);
00100 
00101         }
00102 
00103         matrix.endindices();
00104 
00105     }
00106 
00107 private:
00108 
00109     std::vector<std::set<unsigned int> > indices_;
00110 
00111     int rows_, cols_;
00112 
00113 };
00114 
00115 
00116 } // end namespace Dune
00117 
00118 #endif