3#ifndef DUNE_ISTL_COLCOMPMATRIX_HH
4#define DUNE_ISTL_COLCOMPMATRIX_HH
58 template<
class M,
class S>
76 const Matrix& matrix()
const
92 typename RowIndexSet::const_iterator pos)
93 : firstRow_(firstRow), pos_(pos)
99 return *(firstRow_+ *pos_);
109 typename RowIndexSet::value_type index()
const
116 typename Matrix::const_iterator firstRow_;
118 typename RowIndexSet::const_iterator pos_;
156 template<
class M,
class X,
class TM,
class TD,
class T1>
159 template<
class T,
bool flag>
160 struct SeqOverlappingSchwarzAssemblerHelper;
166 template<
class B,
class TA,
int n,
int m>
198 size_type nnz()
const
217 int* getRowIndex()
const
222 int* getColStart()
const
227 ColCompMatrix& operator=(
const Matrix& mat);
228 ColCompMatrix& operator=(
const ColCompMatrix& mat);
249 template<
class T,
class A,
int n,
int m>
252 template<
class I,
class S,
class D>
266 template<
typename Iter>
267 void addRowNnz(
const Iter& row)
const;
269 template<
typename Iter,
typename Set>
270 void addRowNnz(
const Iter& row,
const Set& s)
const;
274 template<
typename Iter>
275 void countEntries(
const Iter& row,
const CIter& col)
const;
277 void countEntries(size_type colidx)
const;
279 void calcColstart()
const;
281 template<
typename Iter>
282 void copyValue(
const Iter& row,
const CIter& col)
const;
284 void copyValue(
const CIter& col, size_type rowindex, size_type colidx)
const;
286 virtual void createMatrix()
const;
290 void allocateMatrixStorage()
const;
292 void allocateMarker();
296 mutable size_type *marker;
299 template<
class T,
class A,
int n,
int m>
300 ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::ColCompMatrixInitializer(ColCompMatrix& mat_)
301 : mat(&mat_), cols(mat_.M()), marker(0)
306 template<
class T,
class A,
int n,
int m>
307 ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::ColCompMatrixInitializer()
308 : mat(0), cols(0), marker(0)
311 template<
class T,
class A,
int n,
int m>
312 ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::~ColCompMatrixInitializer()
318 template<
class T,
class A,
int n,
int m>
319 template<
typename Iter>
320 void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::addRowNnz(
const Iter& row)
const
322 mat->Nnz_+=row->getsize();
325 template<
class T,
class A,
int n,
int m>
326 template<
typename Iter,
typename Map>
327 void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::addRowNnz(
const Iter& row,
328 const Map& indices)
const
330 typedef typename Iter::value_type::const_iterator RIter;
331 typedef typename Map::const_iterator MIter;
332 MIter siter =indices.
begin();
333 for(RIter entry=row->begin(); entry!=row->end(); ++entry)
335 for(; siter!=indices.end() && *siter<entry.index(); ++siter) ;
336 if(siter==indices.end())
338 if(*siter==entry.index())
344 template<
class T,
class A,
int n,
int m>
345 void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::allocate()
347 allocateMatrixStorage();
351 template<
class T,
class A,
int n,
int m>
352 void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::allocateMatrixStorage()
const
356 mat->values=
new T[mat->Nnz_];
357 mat->rowindex=
new int[mat->Nnz_];
358 mat->colstart=
new int[cols+1];
361 template<
class T,
class A,
int n,
int m>
362 void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::allocateMarker()
366 for(size_type i=0; i < cols; ++i)
370 template<
class T,
class A,
int n,
int m>
371 template<
typename Iter>
372 void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::countEntries(
const Iter& row,
const CIter& col)
const
375 countEntries(col.index());
378 template<
class T,
class A,
int n,
int m>
379 void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::countEntries(size_type colindex)
const
381 for(size_type i=0; i < m; ++i)
383 assert(colindex*m+i<cols);
384 marker[colindex*m+i]+=n;
388 template<
class T,
class A,
int n,
int m>
389 void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::calcColstart()
const
392 for(size_type i=0; i < cols; ++i) {
394 mat->colstart[i+1]=mat->colstart[i]+marker[i];
395 marker[i]=mat->colstart[i];
399 template<
class T,
class A,
int n,
int m>
400 template<
typename Iter>
401 void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::copyValue(
const Iter& row,
const CIter& col)
const
403 copyValue(col, row.index(), col.index());
406 template<
class T,
class A,
int n,
int m>
407 void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::copyValue(
const CIter& col, size_type rowindex, size_type colindex)
const
409 for(size_type i=0; i<n; i++) {
410 for(size_type j=0; j<m; j++) {
411 assert(colindex*m+j<cols-1 || (
int)marker[colindex*m+j]<mat->colstart[colindex*m+j+1]);
412 assert((
int)marker[colindex*m+j]<mat->Nnz_);
413 mat->rowindex[marker[colindex*m+j]]=rowindex*n+i;
414 mat->values[marker[colindex*m+j]]=(*col)[i][j];
415 ++marker[colindex*m+j];
420 template<
class T,
class A,
int n,
int m>
421 void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::createMatrix()
const
427 template<
class F,
class MRS>
428 void copyToColCompMatrix(F& initializer,
const MRS& mrs)
430 typedef typename MRS::const_iterator Iter;
431 typedef typename std::iterator_traits<Iter>::value_type::const_iterator CIter;
432 for(Iter row=mrs.begin(); row!= mrs.end(); ++row)
433 initializer.addRowNnz(row);
435 initializer.allocate();
437 for(Iter row=mrs.begin(); row!= mrs.end(); ++row) {
439 for(CIter col=row->begin(); col != row->end(); ++col)
440 initializer.countEntries(row, col);
443 initializer.calcColstart();
445 for(Iter row=mrs.begin(); row!= mrs.end(); ++row) {
446 for(CIter col=row->begin(); col != row->end(); ++col) {
447 initializer.copyValue(row, col);
451 initializer.createMatrix();
454 template<
class F,
class M,
class S>
455 void copyToColCompMatrix(F& initializer,
const MatrixRowSubset<M,S>& mrs)
457 typedef MatrixRowSubset<M,S> MRS;
458 typedef typename MRS::RowIndexSet SIS;
459 typedef typename SIS::const_iterator SIter;
460 typedef typename MRS::const_iterator Iter;
461 typedef typename std::iterator_traits<Iter>::value_type row_type;
462 typedef typename row_type::const_iterator CIter;
465 for(Iter row=mrs.begin(); row!= mrs.end(); ++row)
466 initializer.addRowNnz(row, mrs.rowIndexSet());
468 initializer.allocate();
470 typedef typename MRS::Matrix::size_type size_type;
476 std::vector<size_type> subMatrixIndex(mrs.matrix().N(),
477 std::numeric_limits<size_type>::max());
479 for(SIter index = mrs.rowIndexSet().begin(); index!=mrs.rowIndexSet().end(); ++index)
480 subMatrixIndex[*index]=s++;
482 for(Iter row=mrs.begin(); row!= mrs.end(); ++row)
483 for(CIter col=row->begin(); col != row->end(); ++col) {
484 if(subMatrixIndex[col.index()]!=std::numeric_limits<size_type>::max())
486 initializer.countEntries(subMatrixIndex[col.index()]);
489 initializer.calcColstart();
491 for(Iter row=mrs.begin(); row!= mrs.end(); ++row)
492 for(CIter col=row->begin(); col != row->end(); ++col) {
493 if(subMatrixIndex[col.index()]!=std::numeric_limits<size_type>::max())
495 initializer.copyValue(col, subMatrixIndex[row.index()], subMatrixIndex[col.index()]);
497 initializer.createMatrix();
502 template<
class B,
class TA,
int n,
int m>
503 ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::ColCompMatrix()
504 : N_(0), M_(0), Nnz_(0), values(0), rowindex(0), colstart(0)
507 template<
class B,
class TA,
int n,
int m>
508 ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >
509 ::ColCompMatrix(
const Matrix& mat)
510 : N_(n*mat.N()), M_(m*mat.M()), Nnz_(n*m*mat.nonzeroes())
513 template<
class B,
class TA,
int n,
int m>
514 ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >&
515 ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::operator=(
const Matrix& mat)
523 template<
class B,
class TA,
int n,
int m>
524 ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >&
525 ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::operator=(
const ColCompMatrix& mat)
533 colstart=
new int[M_+1];
534 for(
int i=0; i<=M_; ++i)
535 colstart[i]=mat.colstart[i];
539 values =
new B[Nnz_];
540 rowindex =
new int[Nnz_];
542 for(
int i=0; i<Nnz_; ++i)
543 values[i]=mat.values[i];
545 for(
int i=0; i<Nnz_; ++i)
546 rowindex[i]=mat.rowindex[i];
551 template<
class B,
class TA,
int n,
int m>
552 void ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >
553 ::setMatrix(
const Matrix& mat)
557 ColCompMatrixInitializer<Matrix> initializer(*
this);
559 copyToColCompMatrix(initializer, MatrixRowSet<Matrix>(mat));
562 template<
class B,
class TA,
int n,
int m>
563 void ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >
564 ::setMatrix(
const Matrix& mat,
const std::set<std::size_t>& mrs)
570 ColCompMatrixInitializer<Matrix> initializer(*
this);
572 copyToColCompMatrix(initializer, MatrixRowSubset<Matrix,std::set<std::size_t> >(mat,mrs));
575 template<
class B,
class TA,
int n,
int m>
576 ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::~ColCompMatrix()
582 template<
class B,
class TA,
int n,
int m>
583 void ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::free()
Implementation of the BCRSMatrix class.
This file implements a vector space as a tensor product of a given vector space. The number of compon...
A sparse block matrix with compressed row storage.
Definition: bcrsmatrix.hh:423
A::size_type size_type
The type for the index access and the size.
Definition: bcrsmatrix.hh:457
Definition: bvector.hh:650
BCRSMatrix< FieldMatrix< B, n, m >, TA > Matrix
The type of the matrix to convert.
Definition: colcompmatrix.hh:173
ColCompMatrix(const Matrix &mat)
Constructor that initializes the data.
virtual void setMatrix(const Matrix &mat, const std::set< std::size_t > &mrs)
Initialize data from a given set of matrix rows and columns.
size_type N() const
Get the number of rows.
Definition: colcompmatrix.hh:193
virtual void free()
free allocated space.
virtual ~ColCompMatrix()
Destructor.
virtual void setMatrix(const Matrix &mat)
Initialize data from given matrix.
size_type M() const
Get the number of columns.
Definition: colcompmatrix.hh:207
A dense n x m matrix.
Definition: fmatrix.hh:68
Base class for stl conformant forward iterators.
Definition: iteratorfacades.hh:144
ConstIterator class for sequential access.
Definition: matrix.hh:398
Provides access to an iterator over all matrix rows.
Definition: colcompmatrix.hh:22
Matrix::ConstRowIterator const_iterator
The matrix row iterator type.
Definition: colcompmatrix.hh:27
const_iterator begin() const
Get the row iterator at the first row.
Definition: colcompmatrix.hh:38
M Matrix
The type of the matrix.
Definition: colcompmatrix.hh:25
MatrixRowSet(const Matrix &m)
Construct an row set over all matrix rows.
Definition: colcompmatrix.hh:33
const_iterator end() const
Get the row iterator at the end of all rows.
Definition: colcompmatrix.hh:43
The matrix row iterator type.
Definition: colcompmatrix.hh:89
Provides access to an iterator over an arbitrary subset of matrix rows.
Definition: colcompmatrix.hh:60
S RowIndexSet
the type of the set of valid row indices.
Definition: colcompmatrix.hh:65
const_iterator begin() const
Get the row iterator at the first row.
Definition: colcompmatrix.hh:122
M Matrix
the type of the matrix class.
Definition: colcompmatrix.hh:63
MatrixRowSubset(const Matrix &m, const RowIndexSet &s)
Construct an row set over all matrix rows.
Definition: colcompmatrix.hh:72
const_iterator end() const
Get the row iterator at the end of all rows.
Definition: colcompmatrix.hh:127
A generic dynamic dense matrix.
Definition: matrix.hh:555
A::size_type size_type
Type for indices and sizes.
Definition: matrix.hh:571
RowIterator begin()
Get iterator to first row.
Definition: matrix.hh:609
size_type M() const
Return the number of columns.
Definition: matrix.hh:695
size_type N() const
Return the number of rows.
Definition: matrix.hh:690
Initializer for SuperLU Matrices representing the subdomains.
Definition: overlappingschwarz.hh:43
Sequential overlapping Schwarz preconditioner.
Definition: overlappingschwarz.hh:742
RealIterator< const B > const_iterator
iterator class for sequential access
Definition: basearray.hh:204
Implements a matrix constructed from a given type representing a field and compile-time given number ...
Implements a vector constructed from a given type representing a field and a compile-time given size.
Dune namespace.
Definition: alignment.hh:11
Inititializer for the ColCompMatrix as needed by OverlappingSchwarz.
Definition: colcompmatrix.hh:154
Utility class for converting an ISTL Matrix into a column-compressed matrix.
Definition: colcompmatrix.hh:145
Traits for type conversions and type information.
Definition of the DUNE_UNUSED macro for the case that config.h is not available.
#define DUNE_UNUSED_PARAMETER(parm)
A macro to mark intentionally unused function parameters with.
Definition: unused.hh:18