From a1c1a402b7c269ee2296f02e24132f9f1b3b4181 Mon Sep 17 00:00:00 2001 From: David Wells Date: Mon, 20 Aug 2018 00:22:11 -0400 Subject: [PATCH] Add matrix iterators to Table (and FullMatrix). This generalizes the work done in b4b71dc9f60 by adding more template parameters to the iterator and accessor classes for TransposeTable: the new version is now used in Table (and thus FullMatrix). The old FullMatrix::const_iterator has been completely removed (since the new version provides both const and non-const iterators). --- doc/news/changes/minor/20180820DavidWells | 4 + include/deal.II/base/linear_index_iterator.h | 2 +- include/deal.II/base/table.h | 818 +++++++++++++------ include/deal.II/lac/full_matrix.h | 298 +------ tests/full_matrix/full_matrix_iterator_01.cc | 17 +- 5 files changed, 604 insertions(+), 535 deletions(-) create mode 100644 doc/news/changes/minor/20180820DavidWells diff --git a/doc/news/changes/minor/20180820DavidWells b/doc/news/changes/minor/20180820DavidWells new file mode 100644 index 0000000000..73b873eb93 --- /dev/null +++ b/doc/news/changes/minor/20180820DavidWells @@ -0,0 +1,4 @@ +New: Table and FullMatrix now have complete (both mutable and const) +random-access iterators. +
+(David Wells, 2018/08/20) diff --git a/include/deal.II/base/linear_index_iterator.h b/include/deal.II/base/linear_index_iterator.h index f7be30427b..363fb6f940 100644 --- a/include/deal.II/base/linear_index_iterator.h +++ b/include/deal.II/base/linear_index_iterator.h @@ -504,7 +504,7 @@ inline bool LinearIndexIterator:: operator>(const DerivedIterator &other) const { - return other < *this; + return other < static_cast(*this); } diff --git a/include/deal.II/base/table.h b/include/deal.II/base/table.h index bf0ad0a00b..52094a8e8d 100644 --- a/include/deal.II/base/table.h +++ b/include/deal.II/base/table.h @@ -828,6 +828,294 @@ public: +/** + * A namespace for iterators and accessors for Table<2, T> and + * TransposeTable. These classes have special accessors (that is, compared to + * Table<3, T>) since they have a matrix-like structure; i.e., the accessors + * also provide row and column information and are designed to be compatible + * with the SparseMatrix and SparsityPattern iterator classes. + */ +namespace MatrixTableIterators +{ + /** + * @brief Enumeration describing the storage order (i.e., the in-memory + * layout) of a table class. + */ + enum class Storage + { + /** + * The data are organized in row-major (i.e., C-style) order. + */ + row_major, + + /** + * The data are organized in column-major (i.e., Fortran-style) order. + */ + column_major + }; + + // Forward declaration of the iterator class. + template + class Iterator; + + /** + * @brief %Accessor class template. This class is partially specialized for + * both values of Constness. + */ + template + class Accessor; + + /** + * @brief %Accessor base class for Table<2, T> and TransposeTable. + * + * This class is compatible with the requirements for an %Accessor described + * in LinearIndexIterator: See the documentation of that class for a + * description of the split between iterators and accessors. + * + * @tparam TableType the type of the %Table, e.g., Table<2, T> or + * TransposeTable. + * + * @tparam Constness whether or not this object stores a constant pointer + * and can modify the provided table object. + * + * @tparam storage_order The storage scheme of the underlying table, e.g., + * Storage::row_major for Table<2, T>. + */ + template + class AccessorBase + { + public: + /** + * Type of the stored pointer to the table. + */ + using container_pointer_type = typename std:: + conditional::type; + + /** + * Value type of the underlying container. + */ + using value_type = typename TableType::value_type; + + /** + * Numerical type of the row and column indices of the table. + */ + using size_type = typename TableType::size_type; + + /** + * Default constructor. + */ + AccessorBase(); + + /** + * Constructor setting up the end iterator. + */ + AccessorBase(const container_pointer_type table); + + /** + * Copy constructor from a non-const Accessor. + */ + AccessorBase(const AccessorBase &); + + /** + * Constructor taking an array index. + */ + AccessorBase(const container_pointer_type table, + const std::ptrdiff_t linear_index); + + /** + * Get a constant reference to the value of the element represented by + * this accessor. + */ + const value_type & + value() const; + + /** + * Return the row of the current entry. + */ + size_type + row() const; + + /** + * Return the column of the current entry. + */ + size_type + column() const; + + protected: + /** + * Pointer to the table. + */ + container_pointer_type container; + + /** + * Current index. + */ + std::ptrdiff_t linear_index; + + /** + * Check that linear_index corresponds to an entry that is + * actually stored by the table (i.e., assert that + * linear_index is nonnegative and less than + * container->size()). + */ + void + assert_valid_linear_index() const; + + /** + * Make the const version a friend for copying. + */ + friend class AccessorBase; + + /** + * Make the underlying iterator class a friend. + */ + friend class LinearIndexIterator< + Iterator, + Accessor>; + }; + + /** + * @brief %Accessor class offering read-only access to elements of a + * table. This is the same as the base class. + */ + template + class Accessor + : public AccessorBase + { + public: + /** + * Use the base class value type. + */ + using value_type = + typename AccessorBase::value_type; + + /** + * Use the base class size type. + */ + using size_type = + typename AccessorBase::size_type; + + /** + * Inherit the base class constructors. + */ + using AccessorBase::AccessorBase; + }; + + /** + * @brief %Accessor class offering read and write access to the elements of + * a table. + */ + template + class Accessor + : public AccessorBase + { + public: + /** + * Use the base class value type. + */ + using value_type = + typename AccessorBase::value_type; + + /** + * Use the base class size type. + */ + using size_type = + typename AccessorBase::size_type; + + /** + * Inherit the base class constructors. + */ + using AccessorBase::AccessorBase; + + /** + * Assignment operator. This assigns a new value to the table entry at the + * current row and column coordinates. + */ + const Accessor & + operator=(const value_type &) const; + + /** + * Move assignment operator. This assigns a new value to the table entry at + * the current row and column coordinates. + */ + const Accessor & + operator=(value_type &&) const; + + /** + * Since we overload value() we have to explicitly use the base class + * version. + */ + using AccessorBase::value; + + /** + * Get a reference to the value of the element represented by + * this accessor. + */ + value_type & + value() const; + }; + + /** + * @brief %Iterator class for both matrix-like tables, i.e., Table<2, T> and + * TransposeTable. + * + * @tparam TableType the type of the %Table, e.g., Table<2, T> or + * TransposeTable. + * + * @tparam Constness whether or not this is a constant iterator. + * + * @tparam storage_order The storage scheme of the underlying table, e.g., + * Storage::row_major for Table<2, T>. + */ + template + class Iterator + : public LinearIndexIterator, + Accessor> + { + public: + /** + * Size type used by the underlying table. + */ + using size_type = typename TableType::size_type; + + /** + * Type of the stored pointer to the table. + */ + using container_pointer_type = typename std:: + conditional::type; + + /** + * Constructor from an accessor. + */ + Iterator(const Accessor &accessor); + + /** + * Constructor. Create the end iterator for a table. + */ + Iterator(const container_pointer_type object); + + /** + * Constructor for a particular table entry. + */ + Iterator(const container_pointer_type object, + const size_type row, + const size_type column); + + /** + * Copy constructor from a non-const iterator. + */ + Iterator(const Iterator &i); + + /** + * Constructor for an entry with a particular linear index. + */ + Iterator(const container_pointer_type container, + const std::ptrdiff_t linear_index); + }; +} // namespace MatrixTableIterators + + + /** * A class representing a two-dimensional table, i.e. a matrix of objects (not * necessarily only numbers). @@ -850,6 +1138,34 @@ public: */ using size_type = typename TableBase<2, T>::size_type; + /** + * Typedef for the values in the table. + */ + using value_type = typename AlignedVector::value_type; + + /** + * Typedef for the references in the table. + */ + using reference = typename AlignedVector::reference; + + /** + * Typedef for the constant references in the table. + */ + using const_reference = typename AlignedVector::const_reference; + + /** + * Typedef for a constant iterator that traverses the table in column-major + * order. + */ + using const_iterator = MatrixTableIterators:: + Iterator, true, MatrixTableIterators::Storage::row_major>; + + /** + * Typedef for an iterator that traverses the table in column-major order. + */ + using iterator = MatrixTableIterators:: + Iterator, false, MatrixTableIterators::Storage::row_major>; + /** * Default constructor. Set all dimensions to zero. */ @@ -982,6 +1298,30 @@ public: size_type n_cols() const; + /** + * Return an iterator pointing to the first entry. + */ + iterator + begin(); + + /** + * Return a constant iterator pointing to the first entry. + */ + const_iterator + begin() const; + + /** + * Return an iterator pointing to one past the last entry. + */ + iterator + end(); + + /** + * Return a constant iterator pointing to one past the last entry. + */ + const_iterator + end() const; + protected: /** * Return a read-write reference to the element (i,j). @@ -1012,6 +1352,22 @@ protected: */ typename AlignedVector::const_reference el(const size_type i, const size_type j) const; + + /** + * Make the AccessorBase class a friend so that it may directly index into + * the array. + */ + friend class MatrixTableIterators:: + AccessorBase, true, MatrixTableIterators::Storage::row_major>; + friend class MatrixTableIterators:: + AccessorBase, false, MatrixTableIterators::Storage::row_major>; + + /** + * Make the mutable accessor class a friend so that we can write to array + * entries with iterators. + */ + friend class MatrixTableIterators:: + Accessor, false, MatrixTableIterators::Storage::row_major>; }; @@ -1538,221 +1894,33 @@ public: }; - /** * A namespace for iterators for TransposeTable. TransposeTable is unique in * that it stores entries in column-major order. + * + * @warning The classes defined in this namespace have been deprecated in + * favor of the more general versions in MatrixTableIterators. Use those + * instead. */ namespace TransposeTableIterators { - // Forward declaration of the iterator class. template - class Iterator; + using AccessorBase DEAL_II_DEPRECATED = MatrixTableIterators::AccessorBase< + TransposeTable, + Constness, + MatrixTableIterators::Storage::column_major>; - /** - * Accessor class template. This will be partially specialized for both - * values of Constness. - */ template - class Accessor; + using Accessor DEAL_II_DEPRECATED = + MatrixTableIterators::Accessor, + Constness, + MatrixTableIterators::Storage::column_major>; - /** - * Accessor base class. Stores a const pointer when Constness == - * true and a non-const pointer when Constness == false. - */ template - class AccessorBase - { - public: - /** - * Type of the stored pointer to the TransposeTable. - */ - using container_pointer_type = - typename std::conditional *, - TransposeTable *>::type; - - /** - * Default constructor. - */ - AccessorBase(); - - /** - * Constructor setting up the end iterator. - */ - AccessorBase(const container_pointer_type table); - - /** - * Copy constructor from a non-const Accessor. - */ - AccessorBase(const AccessorBase &); - - /** - * Constructor taking an array index. Note that TransposeTable stores data - * in column-major order. - */ - AccessorBase(const container_pointer_type table, - const std::ptrdiff_t linear_index); - - /** - * Get a constant reference to the value of the element represented by - * this object. - */ - const T & - value() const; - - /** - * Numerical type of the row and column indices of the TransposeTable. - */ - using size_type = typename TransposeTable::size_type; - - /** - * Access the row of the current entry. - */ - size_type - row() const; - - /** - * Access the column of the current entry. - */ - size_type - column() const; - - protected: - /** - * Pointer to the TransposeTable. - */ - container_pointer_type container; - - /** - * Current index. - */ - std::ptrdiff_t linear_index; - - /** - * Check that linear_index corresponds to an entry that is - * actually stored by the Table (i.e., assert that - * linear_index is nonnegative and less than - * container->size()). - */ - void - assert_valid_linear_index() const; - - /** - * Make the const version a friend for copying. - */ - friend class AccessorBase; - - /** - * Make the underlying iterator class a friend. - */ - friend class LinearIndexIterator, - Accessor>; - }; - - /** - * Constant Accessor. Same as the base class. - */ - template - class Accessor : public AccessorBase - { - public: - /** - * Inherit the base class constructors. - */ - using AccessorBase::AccessorBase; - }; - - /** - * Accessor class for non-constant tables. Defines a few extra operations - * for modifying values. - */ - template - class Accessor : public AccessorBase - { - public: - /** - * Inherit the base class constructors. - */ - using AccessorBase::AccessorBase; - - /** - * Assignment operator. This assigns a new value to the table entry at the - * current row and column coordinates. - */ - const Accessor & - operator=(const T &) const; - - /** - * Move assignment operator. This assigns a new value to the table entry at - * the current row and column coordinates. - */ - const Accessor & - operator=(T &&) const; - - /** - * Since we overload value() we have to explicitly use the base class - * version. - */ - using AccessorBase::value; - - /** - * Get a reference to the value of the element represented by - * this object. - */ - T & - value() const; - }; - - /** - * Iterator class for TransposeTable. - */ - template - class Iterator - : public LinearIndexIterator, Accessor> - { - public: - /** - * Size type used by TransposeTable. - */ - using size_type = typename TransposeTable::size_type; - - /** - * Type of the stored pointer to the TransposeTable. - */ - using container_pointer_type = - typename std::conditional *, - TransposeTable *>::type; - - /** - * Constructor from an accessor. - */ - Iterator(const Accessor &accessor); - - /** - * Constructor. Create the end iterator for a table. - */ - Iterator(const container_pointer_type object); - - /** - * Constructor for a particular table entry. - */ - Iterator(const container_pointer_type object, - const size_type row, - const size_type column); - - /** - * Copy constructor from a non-const iterator. - */ - Iterator(const Iterator &i); - - /** - * Constructor for an entry with a particular linear index. - */ - Iterator(const container_pointer_type container, - const std::ptrdiff_t linear_index); - }; + using Iterator DEAL_II_DEPRECATED = + MatrixTableIterators::Iterator, + Constness, + MatrixTableIterators::Storage::column_major>; } // namespace TransposeTableIterators @@ -1797,12 +1965,18 @@ public: * Typedef for a constant iterator that traverses the table in column-major * order. */ - using const_iterator = TransposeTableIterators::Iterator; + using const_iterator = + MatrixTableIterators::Iterator, + true, + MatrixTableIterators::Storage::column_major>; /** * Typedef for an iterator that traverses the table in column-major order. */ - using iterator = TransposeTableIterators::Iterator; + using iterator = + MatrixTableIterators::Iterator, + false, + MatrixTableIterators::Storage::column_major>; /** * Default constructor. Set all dimensions to zero. @@ -1915,14 +2089,23 @@ protected: * Make the AccessorBase class a friend so that it may directly index into * the array. */ - friend class TransposeTableIterators::AccessorBase; - friend class TransposeTableIterators::AccessorBase; + friend class MatrixTableIterators::AccessorBase< + TransposeTable, + true, + MatrixTableIterators::Storage::column_major>; + friend class MatrixTableIterators::AccessorBase< + TransposeTable, + false, + MatrixTableIterators::Storage::column_major>; /** * Make the mutable accessor class a friend so that we can write to array * entries with iterators. */ - friend class TransposeTableIterators::Accessor; + friend class MatrixTableIterators::Accessor< + TransposeTable, + false, + MatrixTableIterators::Storage::column_major>; }; @@ -2617,19 +2800,98 @@ Table<2, T>::n_cols() const +template +inline typename Table<2, T>::iterator +Table<2, T>::begin() +{ + return typename Table<2, T>::iterator(this, 0, 0); +} + + + +template +inline typename Table<2, T>::const_iterator +Table<2, T>::begin() const +{ + return typename Table<2, T>::const_iterator(this, 0, 0); +} + + + +template +inline typename Table<2, T>::iterator +Table<2, T>::end() +{ + return typename Table<2, T>::iterator(this); +} + + + +template +inline typename Table<2, T>::const_iterator +Table<2, T>::end() const +{ + return typename Table<2, T>::const_iterator(this); +} + + + //--------------------------------------------------------------------------- -namespace TransposeTableIterators +namespace MatrixTableIterators { - template - inline AccessorBase::AccessorBase() + namespace internal + { + // Internal calculation routines for AccessorBase. These do not do any + // checking whatsoever. + template + inline std::ptrdiff_t + get_row_index(const std::ptrdiff_t linear_index, + const TableType *const container) + { + switch (storage_order) + { + case Storage::row_major: + return linear_index / container->n_cols(); + case Storage::column_major: + return linear_index % container->n_rows(); + default: + Assert(false, ExcInternalError()); + } + return {}; + } + + + + template + inline std::ptrdiff_t + get_column_index(const std::ptrdiff_t linear_index, + const TableType *const container) + { + switch (storage_order) + { + case Storage::row_major: + return linear_index % container->n_cols(); + case Storage::column_major: + return linear_index / container->n_rows(); + default: + Assert(false, ExcInternalError()); + } + return {}; + } + } // namespace internal + + + + template + inline AccessorBase::AccessorBase() : container(nullptr) , linear_index(std::numeric_limits::max()) {} - template - inline AccessorBase::AccessorBase( + template + inline AccessorBase::AccessorBase( const container_pointer_type table) : container(table) , linear_index(container->values.size()) @@ -2637,17 +2899,17 @@ namespace TransposeTableIterators - template - inline AccessorBase::AccessorBase( - const AccessorBase &a) + template + inline AccessorBase::AccessorBase( + const AccessorBase &a) : container(a.container) , linear_index(a.linear_index) {} - template - inline AccessorBase::AccessorBase( + template + inline AccessorBase::AccessorBase( const container_pointer_type table, const std::ptrdiff_t index) : container(table) @@ -2661,9 +2923,10 @@ namespace TransposeTableIterators - template - inline const T & - AccessorBase::value() const + template + inline const typename AccessorBase:: + value_type & + AccessorBase::value() const { assert_valid_linear_index(); return this->container->values[linear_index]; @@ -2671,40 +2934,52 @@ namespace TransposeTableIterators - template - inline typename AccessorBase::size_type - AccessorBase::AccessorBase::row() const + template + inline typename AccessorBase::size_type + AccessorBase::row() const { assert_valid_linear_index(); - return linear_index % container->n_rows(); + return static_cast( + internal::get_row_index(linear_index, + container)); } - template - inline typename AccessorBase::size_type - AccessorBase::AccessorBase::column() const + template + inline typename AccessorBase::size_type + AccessorBase::column() const { assert_valid_linear_index(); - return linear_index / container->n_rows(); + return static_cast( + internal::get_column_index(linear_index, + container)); } - template + template inline void - AccessorBase::AccessorBase::assert_valid_linear_index() const + AccessorBase::assert_valid_linear_index() + const { # ifdef DEBUG // avoid unused variable warnings by guarding everything + Assert(container != nullptr, + ExcMessage("This accessor has been default-constructed and does not " + "have a corresponding table.")); Assert(!container->empty(), ExcMessage("An empty table has no rows or columns.")); Assert(0 <= linear_index && std::size_t(linear_index) < container->values.size(), ExcMessage("The current iterator points outside of the table.")); - const auto column_n = linear_index / container->n_rows(); + const std::ptrdiff_t row_n = + internal::get_row_index(linear_index, + container); + const std::ptrdiff_t column_n = + internal::get_column_index(linear_index, + container); Assert(0 <= column_n && std::size_t(column_n) < container->n_cols(), ExcMessage("The current iterator points outside the table.")); - const auto row_n = linear_index % container->n_rows(); Assert(0 <= row_n && std::size_t(row_n) < container->n_rows(), ExcMessage("The current iterator points outside the table.")); # endif @@ -2712,9 +2987,11 @@ namespace TransposeTableIterators - template - inline const Accessor & - Accessor::operator=(const T &t) const + template + inline const Accessor & + Accessor::operator=( + const typename Accessor::value_type &t) + const { this->assert_valid_linear_index(); this->container->values[this->linear_index] = t; @@ -2723,9 +3000,10 @@ namespace TransposeTableIterators - template - inline const Accessor & - Accessor::operator=(T &&t) const + template + inline const Accessor & + Accessor::operator=( + typename Accessor::value_type &&t) const { this->assert_valid_linear_index(); this->container->values[this->linear_index] = t; @@ -2734,9 +3012,9 @@ namespace TransposeTableIterators - template - inline T & - Accessor::value() const + template + inline typename Accessor::value_type & + Accessor::value() const { this->assert_valid_linear_index(); return this->container->values[this->linear_index]; @@ -2744,44 +3022,56 @@ namespace TransposeTableIterators - template - Iterator::Iterator(const Accessor &a) - : LinearIndexIterator, Accessor>(a) + template + inline Iterator::Iterator( + const Accessor &a) + : LinearIndexIterator, + Accessor>(a) {} - template - Iterator::Iterator(const container_pointer_type table) - : LinearIndexIterator, Accessor>( - Accessor(table)) + template + inline Iterator::Iterator( + const container_pointer_type table) + : LinearIndexIterator, + Accessor>( + Accessor(table)) {} - template - Iterator::Iterator(const Iterator &i) - : LinearIndexIterator, Accessor>(*i) + template + inline Iterator::Iterator( + const Iterator &i) + : LinearIndexIterator, + Accessor>(*i) {} - template - Iterator::Iterator(const container_pointer_type table, - const size_type row_n, - const size_type col_n) - : Iterator(table, table->n_rows() * col_n + row_n) + template + inline Iterator::Iterator( + const container_pointer_type table, + const size_type row_n, + const size_type col_n) + : Iterator(table, + storage_order == Storage::row_major ? + table->n_cols() * row_n + col_n : + table->n_rows() * col_n + row_n) {} - template - Iterator::Iterator(const container_pointer_type table, - const std::ptrdiff_t linear_index) - : LinearIndexIterator, Accessor>( - Accessor(table, linear_index)) + template + inline Iterator::Iterator( + const container_pointer_type table, + const std::ptrdiff_t linear_index) + : LinearIndexIterator, + Accessor>( + Accessor(table, linear_index)) {} -} // namespace TransposeTableIterators +} // namespace MatrixTableIterators diff --git a/include/deal.II/lac/full_matrix.h b/include/deal.II/lac/full_matrix.h index 77600face8..c8a06cfe19 100644 --- a/include/deal.II/lac/full_matrix.h +++ b/include/deal.II/lac/full_matrix.h @@ -87,6 +87,25 @@ public: */ using value_type = number; + /** + * Use the base class mutable iterator type. + */ + using iterator = typename Table<2, number>::iterator; + + /** + * Use the base class constant iterator type. + */ + using const_iterator = typename Table<2, number>::const_iterator; + + /** + * Use the base class iterator functions. + */ + using Table<2, number>::begin; + + /** + * Use the base class iterator functions + */ + using Table<2, number>::end; /** * Declare a type that has holds real-valued numbers with the same precision @@ -99,129 +118,6 @@ public: */ using real_type = typename numbers::NumberTraits::real_type; - - class const_iterator; - - /** - * Accessor class for iterators - */ - class Accessor - { - public: - /** - * Constructor. Since we use accessors only for read access, a const - * matrix pointer is sufficient. - */ - Accessor(const FullMatrix *matrix, - const size_type row, - const size_type col); - - /** - * Row number of the element represented by this object. - */ - size_type - row() const; - - /** - * Column number of the element represented by this object. - */ - size_type - column() const; - - /** - * Value of this matrix entry. - */ - number - value() const; - - protected: - /** - * The matrix accessed. - */ - const FullMatrix *matrix; - - /** - * Current row number. - */ - size_type a_row; - - /** - * Current column number. - */ - size_type a_col; - - /* - * Make enclosing class a friend. - */ - friend class const_iterator; - }; - - /** - * Standard-conforming iterator. - */ - class const_iterator - { - public: - /** - * Constructor. - */ - const_iterator(const FullMatrix *matrix, - const size_type row, - const size_type col); - - /** - * Prefix increment. - */ - const_iterator & - operator++(); - - /** - * Postfix increment. - */ - const_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 const_iterator &) const; - /** - * Inverse of ==. - */ - bool - operator!=(const 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. - */ - bool - operator<(const const_iterator &) const; - - /** - * Comparison operator. Compares just the other way around than the - * operator above. - */ - bool - operator>(const const_iterator &) const; - - private: - /** - * Store an object of the accessor class. - */ - Accessor accessor; - }; /** * @name Constructors and initialization. See also the base class Table. */ @@ -627,25 +523,25 @@ public: //@{ /** - * Iterator starting at the first entry. + * Mutable iterator starting at the first entry of row r. */ - const_iterator - begin() const; + iterator + begin(const size_type r); /** - * Final iterator. + * One past the end mutable iterator of row r. */ - const_iterator - end() const; + iterator + end(const size_type r); /** - * Iterator starting at the first entry of row r. + * Constant iterator starting at the first entry of row r. */ const_iterator begin(const size_type r) const; /** - * Final iterator of row r. + * One past the end constant iterator of row r. */ const_iterator end(const size_type r) const; @@ -1202,8 +1098,6 @@ public: */ DeclException0(ExcMatrixNotPositiveDefinite); //@} - - friend class Accessor; }; /**@}*/ @@ -1376,147 +1270,25 @@ FullMatrix::Tvmult_add(Vector & w, //--------------------------------------------------------------------------- - - -template -inline FullMatrix::Accessor::Accessor(const FullMatrix *matrix, - const size_type r, - const size_type c) - : matrix(matrix) - , a_row(r) - , a_col(c) -{} - - -template -inline typename FullMatrix::size_type -FullMatrix::Accessor::row() const -{ - return a_row; -} - - -template -inline typename FullMatrix::size_type -FullMatrix::Accessor::column() const -{ - return a_col; -} - - -template -inline number -FullMatrix::Accessor::value() const -{ - AssertIsFinite(matrix->el(a_row, a_col)); - return matrix->el(a_row, a_col); -} - - -template -inline FullMatrix::const_iterator::const_iterator( - const FullMatrix *matrix, - const size_type r, - const size_type c) - : accessor(matrix, r, c) -{} - - -template -inline typename FullMatrix::const_iterator & -FullMatrix::const_iterator::operator++() -{ - Assert(accessor.a_row < accessor.matrix->m(), ExcIteratorPastEnd()); - - ++accessor.a_col; - if (accessor.a_col >= accessor.matrix->n()) - { - accessor.a_col = 0; - accessor.a_row++; - } - return *this; -} - - -template -inline typename FullMatrix::const_iterator -FullMatrix::const_iterator::operator++(int) -{ - const typename FullMatrix::const_iterator current = *this; - ++(*this); - - return current; -} - - -template -inline const typename FullMatrix::Accessor & - FullMatrix::const_iterator::operator*() const -{ - return accessor; -} - - -template -inline const typename FullMatrix::Accessor * - FullMatrix::const_iterator::operator->() const -{ - return &accessor; -} - - template -inline bool -FullMatrix::const_iterator:: -operator==(const const_iterator &other) const +inline typename FullMatrix::iterator +FullMatrix::begin(const size_type r) { - return (accessor.row() == other.accessor.row() && - accessor.column() == other.accessor.column()); -} - - -template -inline bool -FullMatrix::const_iterator:: -operator!=(const const_iterator &other) const -{ - return !(*this == other); -} - - -template -inline bool -FullMatrix::const_iterator::operator<(const const_iterator &other) const -{ - return (accessor.row() < other.accessor.row() || - (accessor.row() == other.accessor.row() && - accessor.column() < other.accessor.column())); + AssertIndexRange(r, m()); + return iterator(this, r, 0); } -template -inline bool -FullMatrix::const_iterator::operator>(const const_iterator &other) const -{ - return (other < *this); -} - template -inline typename FullMatrix::const_iterator -FullMatrix::begin() const +inline typename FullMatrix::iterator +FullMatrix::end(const size_type r) { - return const_iterator(this, 0, 0); + AssertIndexRange(r, m()); + return iterator(this, r + 1, 0); } -template -inline typename FullMatrix::const_iterator -FullMatrix::end() const -{ - return const_iterator(this, m(), 0); -} - template inline typename FullMatrix::const_iterator diff --git a/tests/full_matrix/full_matrix_iterator_01.cc b/tests/full_matrix/full_matrix_iterator_01.cc index 9c525b17b1..ee8ac21cdb 100644 --- a/tests/full_matrix/full_matrix_iterator_01.cc +++ b/tests/full_matrix/full_matrix_iterator_01.cc @@ -22,17 +22,21 @@ #include "../tests.h" +template void test() { FullMatrix A(3, 3); // test prefix operator - const FullMatrix::const_iterator k = A.begin(), j = ++A.begin(); + const IteratorType k = A.begin(), j = ++A.begin(); AssertThrow(k < j, ExcInternalError()); AssertThrow(j > k, ExcInternalError()); + AssertThrow(k <= j, ExcInternalError()); + AssertThrow(j >= k, ExcInternalError()); + AssertThrow(!(j < k), ExcInternalError()); AssertThrow(!(k > j), ExcInternalError()); @@ -43,15 +47,12 @@ test() AssertThrow(!(k != k), ExcInternalError()); // test postfix operator - FullMatrix::const_iterator l = A.begin(); - FullMatrix::const_iterator m = l++; + IteratorType l = A.begin(); + IteratorType m = l++; AssertThrow(m == k, ExcInternalError()); AssertThrow(l > m, ExcInternalError()); AssertThrow(l->column() == k->column() + 1, ExcInternalError()); - - - deallog << "OK" << std::endl; } @@ -63,7 +64,9 @@ main() try { - test(); + test::iterator>(); + test::const_iterator>(); + deallog << "OK" << std::endl; } catch (std::exception &exc) { -- 2.39.5