From 390b6efb1f3ab1039dede50c3c642f7e02bbf9f0 Mon Sep 17 00:00:00 2001 From: wolf Date: Mon, 3 May 2004 13:43:06 +0000 Subject: [PATCH] Pull iterator/accessor classes out of the mother class. Factorize some parts to make generalization later on easier. git-svn-id: https://svn.dealii.org/trunk@9137 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/block_sparse_matrix.h | 583 ++++++++++-------- .../lac/block_sparse_matrix.templates.h | 49 +- deal.II/lac/include/lac/block_vector.h | 2 +- 3 files changed, 356 insertions(+), 278 deletions(-) diff --git a/deal.II/lac/include/lac/block_sparse_matrix.h b/deal.II/lac/include/lac/block_sparse_matrix.h index 9cc422c38d..5547148f01 100644 --- a/deal.II/lac/include/lac/block_sparse_matrix.h +++ b/deal.II/lac/include/lac/block_sparse_matrix.h @@ -23,62 +23,41 @@ template class Vector; template class BlockVector; +template class BlockSparseMatrix; /*! @addtogroup Matrix1 *@{ */ + +namespace internal +{ + /** - * Blocked sparse matrix. The behaviour of objects of this type is - * almost as for the SparseMatrix objects, with most of the - * functions being implemented in both classes. The main difference is - * that the matrix represented by this object is composed of an array - * of sparse matrices (i.e. of type SparseMatrix) and all - * accesses to the elements of this object are relayed to accesses of - * the base matrices. - * - * In addition to the usual matrix access and linear algebra - * functions, there are functions block() which allow access to the - * different blocks of the matrix. This may, for example, be of help - * when you want to implement Schur complement methods, or block - * preconditioners, where each block belongs to a specific component - * of the equation you are presently discretizing. - * - * Note that the numbers of blocks and rows are implicitly determined - * by the sparsity pattern objects used. - * - * Objects of this type are frequently used when a system of differential - * equations has solutions with variables that fall into different - * classes. For example, solutions of the Stokes or Navier-Stokes equations - * have @p dim velocity components and one pressure component. In this case, - * it may make sense to consider the linear system of equations as a system of - * 2x2 blocks, and one can construct preconditioners or solvers based on this - * 2x2 block structure. This class can help you in these cases, as it allows - * to view the matrix alternatively as one big matrix, or as a number of - * individual blocks. - * - * @ref Instantiations: some (@ @) + * Namespace in which iterators in block matrices are implemented. * - * @author Wolfgang Bangerth, 2000 + * @author Wolfgang Bangerth, 2004 */ -template -class BlockSparseMatrix : public Subscriptor -{ - public: - /** - * Type of matrix entries. In analogy to - * the STL container classes. - */ - typedef number value_type; - - class const_iterator; + namespace BlockMatrixIterators + { + template class ConstIterator; + /** * Accessor class for iterators */ + template class Accessor { public: + /** + * Typedef the value type of the + * matrix we point into. + */ + typedef + typename BlockSparseMatrix::value_type + value_type; + /** * Constructor. Since we use * accessors only for read @@ -87,7 +66,7 @@ class BlockSparseMatrix : public Subscriptor */ Accessor (const BlockSparseMatrix *m, const unsigned int row, - const unsigned short index); + const unsigned int index); /** * Row number of the element @@ -101,7 +80,7 @@ class BlockSparseMatrix : public Subscriptor * represented by this * object. */ - unsigned short index() const; + unsigned int index() const; /** * Column number of the @@ -113,7 +92,7 @@ class BlockSparseMatrix : public Subscriptor /** * Value of this matrix entry. */ - number value() const; + value_type value() const; /** * Block row of the @@ -139,7 +118,9 @@ class BlockSparseMatrix : public Subscriptor * Iterator of the underlying matrix * class. */ - typename SparseMatrix::const_iterator base_iterator; + typename + BlockSparseMatrix::BlockType::const_iterator + base_iterator; /** * Number of block where row lies in. @@ -167,41 +148,47 @@ class BlockSparseMatrix : public Subscriptor */ unsigned int a_index; - friend class const_iterator; + /** + * Let the iterator class be a + * friend. + */ + template + friend class ConstIterator; }; /** * STL conforming iterator. */ - class const_iterator : private Accessor + template + class ConstIterator : private Accessor { public: /** * Constructor. */ - const_iterator(const BlockSparseMatrix*, - unsigned int row, - unsigned short index); + ConstIterator(const BlockSparseMatrix*, + const unsigned int row, + const unsigned int index); /** * Prefix increment. */ - const_iterator& operator++ (); + ConstIterator& operator++ (); /** * Postfix increment. */ - const_iterator& operator++ (int); + ConstIterator& operator++ (int); /** * Dereferencing operator. */ - const Accessor& operator* () const; + const Accessor & operator* () const; /** * Dereferencing operator. */ - const Accessor* operator-> () const; + const Accessor * operator-> () const; /** * Comparison. True, if @@ -209,11 +196,11 @@ class BlockSparseMatrix : public Subscriptor * the same matrix * position. */ - bool operator == (const const_iterator&) const; + bool operator == (const ConstIterator&) const; /** * Inverse of operator==(). */ - bool operator != (const const_iterator&) const; + bool operator != (const ConstIterator&) const; /** * Comparison @@ -224,8 +211,67 @@ class BlockSparseMatrix : public Subscriptor * equal and the first * index is smaller. */ - bool operator < (const const_iterator&) const; + bool operator < (const ConstIterator&) const; }; + } +} + + +/** + * Blocked sparse matrix. The behaviour of objects of this type is + * almost as for the SparseMatrix objects, with most of the + * functions being implemented in both classes. The main difference is + * that the matrix represented by this object is composed of an array + * of sparse matrices (i.e. of type SparseMatrix) and all + * accesses to the elements of this object are relayed to accesses of + * the base matrices. + * + * In addition to the usual matrix access and linear algebra + * functions, there are functions block() which allow access to the + * different blocks of the matrix. This may, for example, be of help + * when you want to implement Schur complement methods, or block + * preconditioners, where each block belongs to a specific component + * of the equation you are presently discretizing. + * + * Note that the numbers of blocks and rows are implicitly determined + * by the sparsity pattern objects used. + * + * Objects of this type are frequently used when a system of differential + * equations has solutions with variables that fall into different + * classes. For example, solutions of the Stokes or Navier-Stokes equations + * have @p dim velocity components and one pressure component. In this case, + * it may make sense to consider the linear system of equations as a system of + * 2x2 blocks, and one can construct preconditioners or solvers based on this + * 2x2 block structure. This class can help you in these cases, as it allows + * to view the matrix alternatively as one big matrix, or as a number of + * individual blocks. + * + * @ref Instantiations: some (@ @) + * + * @author Wolfgang Bangerth, 2000 + */ +template +class BlockSparseMatrix : public Subscriptor +{ + public: + /** + * Typedef the type of the underlying + * matrix. + */ + typedef SparseMatrix BlockType; + + /** + * Type of matrix entries. In analogy to + * the STL container classes. + */ + typedef typename BlockType::value_type value_type; + typedef value_type *pointer; + typedef const value_type *const_pointer; + typedef internal::BlockMatrixIterators::ConstIterator iterator; + typedef internal::BlockMatrixIterators::ConstIterator const_iterator; + typedef value_type &reference; + typedef const value_type &const_reference; + typedef size_t size_type; /** * Constructor; initializes the @@ -333,7 +379,7 @@ class BlockSparseMatrix : public Subscriptor * Access the block with the * given coordinates. */ - SparseMatrix & + BlockType & block (const unsigned int row, const unsigned int column); @@ -343,7 +389,7 @@ class BlockSparseMatrix : public Subscriptor * given coordinates. Version for * constant objects. */ - const SparseMatrix & + const BlockType & block (const unsigned int row, const unsigned int column) const; @@ -825,12 +871,12 @@ class BlockSparseMatrix : public Subscriptor * STL-like iterator with the * first entry of row r. */ - const_iterator begin (unsigned int r) const; + const_iterator begin (const unsigned int r) const; /** * Final iterator of row r. */ - const_iterator end (unsigned int r) const; + const_iterator end (const unsigned int r) const; /** * Determine an estimate for the @@ -885,10 +931,25 @@ class BlockSparseMatrix : public Subscriptor /** * Array of sub-matrices. */ - Table<2,SmartPointer > > sub_objects; + Table<2,SmartPointer > sub_objects; + /** + * Make the iterator class a + * friend. We have to work around + * a compiler bug here again. + */ +#ifndef DEAL_II_NAMESP_TEMPL_FRIEND_BUG + template + friend class internal::BlockMatrixIterators::Accessor; + + template + friend class internal::BlockMatrixIterators::ConstIterator; +#else + typedef internal::BlockMatrixIterators::Accessor Accessor; friend class Accessor; + friend class const_iterator; +#endif }; @@ -896,179 +957,184 @@ class BlockSparseMatrix : public Subscriptor /*@}*/ /* ------------------------- Template functions ---------------------- */ -template -inline -BlockSparseMatrix::Accessor:: -Accessor (const BlockSparseMatrix *matrix, - const unsigned int r, - const unsigned short i) - : - matrix(matrix), - base_iterator(matrix->block(0,0).begin()), - row_block(0), - row_start(0), - col_block(0), - col_start(0), - a_index(0) -{ - Assert (i==0, ExcNotImplemented()); - if (r < matrix->m()) - { - std::pair indices - = matrix->sparsity_pattern->get_row_indices().global_to_local(r); - row_block = indices.first; - base_iterator = matrix->block(indices.first, 0).begin(indices.second); - row_start = matrix->sparsity_pattern - ->get_row_indices().local_to_global(row_block, 0); - } - else +namespace internal +{ + namespace BlockMatrixIterators + { + template + inline + Accessor:: + Accessor (const BlockSparseMatrix *matrix, + const unsigned int r, + const unsigned int i) + : + matrix(matrix), + base_iterator(matrix->block(0,0).begin()), + row_block(0), + row_start(0), + col_block(0), + col_start(0), + a_index(0) { - row_block = matrix->n_block_rows(); - base_iterator = matrix->block(0, 0).begin(); + Assert (i==0, ExcNotImplemented()); + + if (r < matrix->m()) + { + std::pair indices + = matrix->sparsity_pattern->get_row_indices().global_to_local(r); + row_block = indices.first; + base_iterator = matrix->block(indices.first, 0).begin(indices.second); + row_start = matrix->sparsity_pattern + ->get_row_indices().local_to_global(row_block, 0); + } + else + { + row_block = matrix->n_block_rows(); + base_iterator = matrix->block(0, 0).begin(); + } } -} -template -inline -unsigned int -BlockSparseMatrix::Accessor::row() const -{ - return row_start + base_iterator->row(); -} + template + inline + unsigned int + Accessor::row() const + { + return row_start + base_iterator->row(); + } -template -inline -short unsigned int -BlockSparseMatrix::Accessor::index() const -{ - return a_index; -} + template + inline + unsigned int + Accessor::index() const + { + return a_index; + } -template -inline -unsigned int -BlockSparseMatrix::Accessor::column() const -{ - return col_start + base_iterator->column(); -} + template + inline + unsigned int + Accessor::column() const + { + return col_start + base_iterator->column(); + } -template -inline -unsigned int -BlockSparseMatrix::Accessor::block_row() const -{ - return row_block; -} + template + inline + unsigned int + Accessor::block_row() const + { + return row_block; + } -template -inline -unsigned int -BlockSparseMatrix::Accessor::block_column() const -{ - return col_block; -} + template + inline + unsigned int + Accessor::block_column() const + { + return col_block; + } -template -inline -number -BlockSparseMatrix::Accessor::value () const -{ - return base_iterator->value(); -} + template + inline + typename Accessor::value_type + Accessor::value () const + { + return base_iterator->value(); + } //----------------------------------------------------------------------// -template -inline -BlockSparseMatrix::const_iterator:: -const_iterator(const BlockSparseMatrix* m, - unsigned int r, - unsigned short i) - : - BlockSparseMatrix::Accessor(m, r, i) -{} + template + inline + ConstIterator:: + ConstIterator (const BlockSparseMatrix* m, + const unsigned int r, + const unsigned int i) + : + Accessor(m, r, i) + {} -template -inline -typename BlockSparseMatrix::const_iterator& -BlockSparseMatrix::const_iterator::operator++ () -{ - Assert (this->row_blockmatrix->n_block_rows(), ExcIteratorPastEnd()); + template + inline + ConstIterator & + ConstIterator::operator++ () + { + Assert (this->row_blockmatrix->n_block_rows(), ExcIteratorPastEnd()); - // Remember current row inside block - unsigned int local_row = this->base_iterator->row(); + // Remember current row inside block + unsigned int local_row = this->base_iterator->row(); - // Advance inside block - ++this->base_iterator; - ++this->a_index; + // Advance inside block + ++this->base_iterator; + ++this->a_index; - // If end of row inside block, - // advance to next block - if (this->base_iterator == - this->matrix->block(this->row_block, this->col_block).end(local_row)) - { - if (this->col_block < this->matrix->n_block_cols()-1) - { - // Advance to next block in - // row + // If end of row inside block, + // advance to next block + if (this->base_iterator == + this->matrix->block(this->row_block, this->col_block).end(local_row)) + { + if (this->col_block < this->matrix->n_block_cols()-1) + { + // Advance to next block in + // row //TODO: what if this row in that block is empty? - ++this->col_block; - this->col_start = this->matrix->sparsity_pattern - ->get_column_indices().local_to_global(this->col_block, 0); - } - else - { - // Advance to first block - // in next row - this->col_block = 0; - this->col_start = 0; - this->a_index = 0; - ++local_row; - if (local_row >= this->matrix->block(this->row_block,0).m()) - { - // If final row in - // block, go to next - // block row - local_row = 0; - ++this->row_block; - if (this->row_block < this->matrix->n_block_rows()) - this->row_start = this->matrix->sparsity_pattern - ->get_row_indices().local_to_global(this->row_block, 0); - } - } - // Finally, set base_iterator - // to start of row determined - // above - if (this->row_block < this->matrix->n_block_rows()) - this->base_iterator = this->matrix->block(this->row_block, this->col_block).begin(local_row); - else - // Set base_iterator to a - // defined state for - // comparison. This is the - // end() state. + ++this->col_block; + this->col_start = this->matrix->sparsity_pattern + ->get_column_indices().local_to_global(this->col_block, 0); + } + else + { + // Advance to first block + // in next row + this->col_block = 0; + this->col_start = 0; + this->a_index = 0; + ++local_row; + if (local_row >= this->matrix->block(this->row_block,0).m()) + { + // If final row in + // block, go to next + // block row + local_row = 0; + ++this->row_block; + if (this->row_block < this->matrix->n_block_rows()) + this->row_start = this->matrix->sparsity_pattern + ->get_row_indices().local_to_global(this->row_block, 0); + } + } + // Finally, set base_iterator + // to start of row determined + // above + if (this->row_block < this->matrix->n_block_rows()) + this->base_iterator = this->matrix->block(this->row_block, this->col_block).begin(local_row); + else + // Set base_iterator to a + // defined state for + // comparison. This is the + // end() state. //TODO: this is a particularly bad choice, since it is actually a valid iterator. it should rather be something like the end iterator of the last block! - this->base_iterator = this->matrix->block(0, 0).begin(); + this->base_iterator = this->matrix->block(0, 0).begin(); + } + return *this; } - return *this; -} // template // inline // const_iterator& -// BlockSparseMatrix::const_iterator::operator++ (int) +// ConstIterator::operator++ (int) // { // Assert (false, ExcNotImplemented()); // } @@ -1076,73 +1142,76 @@ BlockSparseMatrix::const_iterator::operator++ () -template -inline -const typename BlockSparseMatrix::Accessor& -BlockSparseMatrix::const_iterator::operator* () const -{ - return *this; -} + template + inline + const Accessor & + ConstIterator::operator* () const + { + return *this; + } -template -inline -const typename BlockSparseMatrix::Accessor* -BlockSparseMatrix::const_iterator::operator-> () const -{ - return this; -} + template + inline + const Accessor * + ConstIterator::operator-> () const + { + return this; + } -template -inline -bool -BlockSparseMatrix::const_iterator:: -operator == (const const_iterator& i) const -{ - if (this->matrix != i->matrix) - return false; + template + inline + bool + ConstIterator:: + operator == (const ConstIterator& i) const + { + if (this->matrix != i->matrix) + return false; - if (this->row_block == i->row_block - && this->col_block == i->col_block - && this->base_iterator == i->base_iterator) - return true; - return false; -} + if (this->row_block == i->row_block + && this->col_block == i->col_block + && this->base_iterator == i->base_iterator) + return true; + return false; + } -template -inline -bool -BlockSparseMatrix::const_iterator:: -operator != (const const_iterator& i) const -{ - return !(*this == i); -} + template + inline + bool + ConstIterator:: + operator != (const ConstIterator& i) const + { + return !(*this == i); + } -template -inline -bool -BlockSparseMatrix::const_iterator:: -operator < (const const_iterator& i) const -{ - if (this->row_block < i->row_block) - return true; - if (this->row_block == i->row_block) + template + inline + bool + ConstIterator:: + operator < (const ConstIterator& i) const { - if (this->base_iterator->row() < i->base_iterator->row()) - return true; - if (this->base_iterator->row() == i->base_iterator->row()) - { - if (this->a_index < i->a_index) - return true; - } + if (this->row_block < i->row_block) + return true; + if (this->row_block == i->row_block) + { + if (this->base_iterator->row() < i->base_iterator->row()) + return true; + if (this->base_iterator->row() == i->base_iterator->row()) + { + if (this->a_index < i->a_index) + return true; + } + } + return false; } - return false; + + } } //----------------------------------------------------------------------// @@ -1169,7 +1238,7 @@ BlockSparseMatrix::n_block_rows () const template inline -SparseMatrix & +typename BlockSparseMatrix::BlockType & BlockSparseMatrix::block (const unsigned int row, const unsigned int column) { @@ -1183,7 +1252,7 @@ BlockSparseMatrix::block (const unsigned int row, template inline -const SparseMatrix & +const typename BlockSparseMatrix::BlockType & BlockSparseMatrix::block (const unsigned int row, const unsigned int column) const { @@ -1705,7 +1774,7 @@ BlockSparseMatrix::end () const template inline typename BlockSparseMatrix::const_iterator -BlockSparseMatrix::begin (unsigned int r) const +BlockSparseMatrix::begin (const unsigned int r) const { Assert (r::begin (unsigned int r) const template inline typename BlockSparseMatrix::const_iterator -BlockSparseMatrix::end (unsigned int r) const +BlockSparseMatrix::end (const unsigned int r) const { Assert (r -BlockSparseMatrix::BlockSparseMatrix () : +BlockSparseMatrix::BlockSparseMatrix () + : rows (0), columns (0), sparsity_pattern (0) @@ -31,7 +32,8 @@ BlockSparseMatrix::BlockSparseMatrix () : template BlockSparseMatrix:: -BlockSparseMatrix (const BlockSparsityPattern &sparsity) : +BlockSparseMatrix (const BlockSparsityPattern &sparsity) + : rows (0), columns (0) { @@ -48,10 +50,10 @@ BlockSparseMatrix::~BlockSparseMatrix () for (unsigned int r=0; r *p = sub_objects[r][c]; + BlockType *p = sub_objects[r][c]; sub_objects[r][c] = 0; delete p; - }; + } } @@ -70,6 +72,7 @@ operator = (const BlockSparseMatrix &m) for (unsigned int r=0; r *p = sub_objects[r][c]; + BlockType *p = sub_objects[r][c]; sub_objects[r][c] = 0; delete p; - }; + } + sub_objects.clear (); // then associate new sparsity @@ -113,7 +117,7 @@ reinit (const BlockSparsityPattern &sparsity) for (unsigned int r=0; r *p = new SparseMatrix(); + BlockType *p = new SparseMatrix(); p->reinit (sparsity.block(r,c)); sub_objects[r][c] = p; } @@ -141,6 +145,7 @@ BlockSparseMatrix::empty () const for (unsigned int c=0; c::n_actually_nonzero_elements () const for (unsigned int i=0; in_actually_nonzero_elements (); + return count; } @@ -177,26 +183,28 @@ BlockSparseMatrix::get_sparsity_pattern () const } + template void -BlockSparseMatrix::print_formatted (std::ostream &out, - const unsigned int precision, - const bool scientific, - const unsigned int width, - const char *zero_string, - const double denominator) const +BlockSparseMatrix:: +print_formatted (std::ostream &out, + const unsigned int precision, + const bool scientific, + const unsigned int width, + const char *zero_string, + const double denominator) const { for (unsigned int r=0;r unsigned int BlockSparseMatrix::memory_consumption () const @@ -206,6 +214,7 @@ BlockSparseMatrix::memory_consumption () const for (unsigned int r=0; r const_iterator; typedef value_type &reference; typedef const value_type &const_reference; - typedef size_t size_type; + typedef std::size_t size_type; /** * Constructor. There are three -- 2.39.5