From 07cd78ce1090a2da46c6142d984b91533fed6db3 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 13 Apr 2015 08:50:46 -0500 Subject: [PATCH] Introduce iterators over DynamicSparsityPattern objects. --- doc/news/changes.h | 7 + .../deal.II/lac/dynamic_sparsity_pattern.h | 534 +++++++++++++++++- 2 files changed, 540 insertions(+), 1 deletion(-) diff --git a/doc/news/changes.h b/doc/news/changes.h index 8be699e962..cc55557c02 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -440,6 +440,13 @@ inconvenience this causes.

Specific improvements

    +
  1. New: The DynamicSparsityPattern class (formerly called + CompressedSparsityPattern) now has an iterator class that allows to + walk over the nonzero elements of a matrix represented by this class. +
    + (Wolfgang Bangerth, 2015/04/13) +
  2. +
  3. New: The GridGenerator::subdivided_hyper_cube() and GridGenerator::subdivided_hyper_rectangle() now work also for codimension one and two Triangulation; diff --git a/include/deal.II/lac/dynamic_sparsity_pattern.h b/include/deal.II/lac/dynamic_sparsity_pattern.h index 086493e33d..88af4d84c4 100644 --- a/include/deal.II/lac/dynamic_sparsity_pattern.h +++ b/include/deal.II/lac/dynamic_sparsity_pattern.h @@ -30,6 +30,7 @@ DEAL_II_NAMESPACE_OPEN template class SparseMatrix; +class DynamicSparsityPattern; /*! @addtogroup Sparsity @@ -37,6 +38,203 @@ template class SparseMatrix; */ +/** + * Iterators on objects of type DynamicSparsityPattern. + */ +namespace DynamicSparsityPatternIterators +{ + // forward declaration + class Iterator; + + /** + * Declare type for container size. + */ + typedef types::global_dof_index size_type; + + /** + * Accessor class for iterators into objects of type + * DynamicSparsityPattern. + * + * Note that this class only allows read access to elements, providing their + * row and column number (or alternatively the index within the complete + * sparsity pattern). It does not allow modifying the sparsity pattern + * itself. + * + * @author Wolfgang Bangerth + * @date 2015 + */ + class Accessor + { + public: + /** + * Constructor. + */ + Accessor (const DynamicSparsityPattern *sparsity_pattern, + const unsigned int row, + const unsigned int index_within_row); + + /** + * Constructor. Construct the end accessor for the given sparsity pattern. + */ + Accessor (const DynamicSparsityPattern *sparsity_pattern); + + /** + * Row number of the element represented by this object. + */ + size_type row () const; + + /** + * Index within the current row of the element represented by this object. + */ + size_type index () const; + + /** + * Column number of the element represented by this object. + */ + size_type column () const; + + /** + * Comparison. True, if both iterators point to the same matrix position. + */ + bool operator == (const Accessor &) const; + + /** + * Comparison operator. Result is true if either the first row number is + * smaller or if the row numbers are equal and the first index is smaller. + * + * This function is only valid if both iterators point into the same + * sparsity pattern. + */ + bool operator < (const Accessor &) const; + + protected: + /** + * The sparsity pattern we operate on accessed. + */ + const DynamicSparsityPattern *sparsity_pattern; + + /** + * The row we currently point into. + */ + unsigned int current_row; + + /** + * The index of the element we currently point to. + */ + std::size_t index_within_row; + + /** + * Move the accessor to the next nonzero entry in the matrix. + */ + void advance (); + + /** + * Grant access to iterator class. + */ + friend class Iterator; + }; + + + + /** + * An iterator class for walking over the elements of a sparsity pattern. + * + * The typical use for these iterators is to iterate over the elements of a + * sparsity pattern (or, since they also serve as the basis for iterating + * over the elements of an associated matrix, over the elements of a sparse + * matrix), or over the elements of individual rows. There is no + * guarantee that the elements of a row are actually traversed in an order + * in which column numbers monotonically increase. See the documentation of the + * SparsityPattern class for more information. + * + * @note This class operates directly on the internal data + * structures of the SparsityPattern class. As a consequence, some + * operations are cheap and some are not. In particular, it is cheap + * to access the column index of the sparsity pattern entry pointed + * to. On the other hand, it is expensive to compute the distance + * between two iterators. As a consequence, when you design + * algorithms that use these iterators, it is common practice to not + * loop over all elements of a sparsity pattern at once, but + * to have an outer loop over all rows and within this loop iterate + * over the elements of this row. This way, you only ever need to + * dereference the iterator to obtain the column indices whereas the + * (expensive) lookup of the row index can be avoided by using the + * loop index instead. + */ + class Iterator + { + public: + /** + * Constructor. Create an iterator into the sparsity pattern @p sp for the + * given global index (i.e., the index of the given element counting from + * the zeroth row). + */ + Iterator (const DynamicSparsityPattern *sp, + const unsigned int row, + const unsigned int index_within_row); + + /** + * Constructor. Create an invalid (end) iterator into the sparsity + * pattern @p sp. + */ + Iterator (const DynamicSparsityPattern *sp); + + /** + * Prefix increment. + */ + Iterator &operator++ (); + + /** + * Postfix increment. + */ + Iterator operator++ (int); + + /** + * Dereferencing operator. + */ + const Accessor &operator* () const; + + /** + * Dereferencing operator. + */ + const Accessor *operator-> () const; + + /** + * Comparison. True, if both iterators point to the same matrix position. + */ + bool operator == (const Iterator &) const; + + /** + * Inverse of ==. + */ + bool operator != (const Iterator &) const; + + /** + * Comparison operator. Result is true if either the first row number is + * smaller or if the row numbers are equal and the first index is smaller. + * + * This function is only valid if both iterators point into the same + * matrix. + */ + bool operator < (const Iterator &) const; + + /** + * Return the distance between the current iterator and the argument. The + * distance is given by how many times one has to apply operator++ to the + * current iterator to get the argument (for a positive return value), or + * operator-- (for a negative return value). + */ + int operator - (const Iterator &p) const; + + private: + /** + * Store an object of the accessor class. + */ + Accessor accessor; + }; +} + + /** * This class acts as an intermediate form of the SparsityPattern class. From * the interface it mostly represents a SparsityPattern object that is kept @@ -56,7 +254,7 @@ template class SparseMatrix; * module. * * This class is an example of the "dynamic" type of - * @ref Sparsity. + * @ref Sparsity. It is used in most tutorial programs in one way or another. * *

    Interface

    * @@ -89,6 +287,25 @@ public: */ typedef types::global_dof_index size_type; + /** + * Typedef an for iterator class that allows to walk over all nonzero elements + * of a sparsity pattern. + * + * Since the iterator does not allow to modify the sparsity pattern, this + * type is the same as that for @p const_iterator. + */ + typedef + DynamicSparsityPatternIterators::Iterator + iterator; + + /** + * Typedef for an iterator class that allows to walk over all nonzero elements + * of a sparsity pattern. + */ + typedef + DynamicSparsityPatternIterators::Iterator + const_iterator; + /** * An iterator that can be used to iterate over the elements of a single * row. The result of dereferencing such an iterator is a column index. @@ -246,6 +463,49 @@ public: size_type column_number (const size_type row, const size_type index) const; + /** + * @name Iterators + */ +// @{ + + /** + * STL-like iterator with the first entry of the matrix. The resulting + * iterator can be used to walk over all nonzero entries of the sparsity + * pattern. + * + * Note the discussion in the general documentation of this class about the + * order in which elements are accessed. + */ + iterator begin () const; + + /** + * Final iterator. + */ + iterator end () const; + + /** + * STL-like iterator with the first entry of row r. + * + * Note that if the given row is empty, i.e. does not contain any nonzero + * entries, then the iterator returned by this function equals + * end(r). Note also that the iterator may not be dereferencable in + * that case. + * + * Note also the discussion in the general documentation of this class about + * the order in which elements are accessed. + */ + iterator begin (const size_type r) const; + + /** + * Final iterator of row r. It points to the first element past the + * end of line @p r, or past the end of the entire sparsity pattern. + * + * Note that the end iterator is not necessarily dereferencable. This is in + * particular the case if it is the end iterator for the last row of a + * matrix. + */ + iterator end (const size_type r) const; + /** * Return an iterator that can loop over all entries in the given row. * Dereferencing the iterator yields a column index. @@ -256,6 +516,9 @@ public: * Returns the end of the current row. */ row_iterator row_end (const size_type row) const; + +// @} + /** * Compute the bandwidth of the matrix represented by this structure. The * bandwidth is the maximum of $|i-j|$ for which the index pair $(i,j)$ @@ -357,12 +620,224 @@ private: * Actual data: store for each row the set of nonzero entries. */ std::vector lines; + + // make the accessor class a friend + friend class DynamicSparsityPatternIterators::Accessor; }; /*@}*/ /*---------------------- Inline functions -----------------------------------*/ +namespace DynamicSparsityPatternIterators +{ + inline + Accessor:: + Accessor (const DynamicSparsityPattern *sparsity_pattern, + const unsigned int row, + const unsigned int index_within_row) + : + sparsity_pattern(sparsity_pattern), + current_row (row), + index_within_row(index_within_row) + {} + + + inline + Accessor:: + Accessor (const DynamicSparsityPattern *sparsity_pattern) + : + sparsity_pattern(sparsity_pattern), + current_row(numbers::invalid_unsigned_int), + index_within_row(numbers::invalid_unsigned_int) + {} + + + + inline + size_type + Accessor::row() const + { + Assert (current_row < sparsity_pattern->n_rows(), + ExcInternalError()); + + return current_row; + } + + + inline + size_type + Accessor::column() const + { + Assert (current_row < sparsity_pattern->n_rows(), + ExcInternalError()); + + return sparsity_pattern->lines[current_row].entries[index_within_row]; + } + + + inline + size_type + Accessor::index() const + { + Assert (current_row < sparsity_pattern->n_rows(), + ExcInternalError()); + + return index_within_row; + } + + + + + inline + bool + Accessor::operator == (const Accessor &other) const + { + return (sparsity_pattern == other.sparsity_pattern && + current_row == other.current_row && + index_within_row == other.index_within_row); + } + + + + inline + bool + Accessor::operator < (const Accessor &other) const + { + Assert (sparsity_pattern == other.sparsity_pattern, + ExcInternalError()); + + return ((current_row < other.current_row) || + ((current_row == other.current_row) && + (index_within_row < other.index_within_row))); + } + + + inline + void + Accessor::advance () + { + Assert (current_row < sparsity_pattern->n_rows(), + ExcInternalError()); + + // move to the next element in this row + ++index_within_row; + + // if this moves us beyond the end of the row, go to the next row + // if possible, or set the iterator to an invalid state if not. + // + // use a while loop to ensure that we properly skip over lines + // that are empty + while (index_within_row == + sparsity_pattern->lines[current_row].entries.size()) + { + ++current_row; + index_within_row = 0; + + if (current_row == sparsity_pattern->n_rows()) + { + *this = Accessor(sparsity_pattern); // invalid object + return; + } + } + } + + + + inline + Iterator::Iterator (const DynamicSparsityPattern *sparsity_pattern, + const unsigned int row, + const unsigned int index_within_row) + : + accessor(sparsity_pattern, row, index_within_row) + {} + + + + inline + Iterator::Iterator (const DynamicSparsityPattern *sparsity_pattern) + : + accessor(sparsity_pattern) + {} + + + + inline + Iterator & + Iterator::operator++ () + { + accessor.advance (); + return *this; + } + + + + inline + Iterator + Iterator::operator++ (int) + { + const Iterator iter = *this; + accessor.advance (); + return iter; + } + + + + inline + const Accessor & + Iterator::operator* () const + { + return accessor; + } + + + + inline + const Accessor * + Iterator::operator-> () const + { + return &accessor; + } + + + inline + bool + Iterator::operator == (const Iterator &other) const + { + return (accessor == other.accessor); + } + + + + inline + bool + Iterator::operator != (const Iterator &other) const + { + return ! (*this == other); + } + + + inline + bool + Iterator::operator < (const Iterator &other) const + { + return accessor < other.accessor; + } + + + inline + int + Iterator::operator - (const Iterator &other) const + { + Assert (accessor.sparsity_pattern == other.accessor.sparsity_pattern, + ExcInternalError()); + Assert (false, ExcNotImplemented()); + + return 0; + } +} + + inline void DynamicSparsityPattern::Line::add (const size_type j) @@ -484,6 +959,63 @@ DynamicSparsityPattern::column_number (const size_type row, +inline +DynamicSparsityPattern::iterator +DynamicSparsityPattern::begin () const +{ + return begin(0); +} + + +inline +DynamicSparsityPattern::iterator +DynamicSparsityPattern::end () const +{ + return iterator(this); +} + + + +inline +DynamicSparsityPattern::iterator +DynamicSparsityPattern::begin (const size_type r) const +{ + Assert (r(r,0,n_rows())); + + // find the first row starting at r that has entries and return the + // begin iterator to it + unsigned int row = r; + while ((row(r,0,n_rows())); + + // find the first row after r that has entries and return the begin + // iterator to it + unsigned int row = r+1; + while ((row