]> https://gitweb.dealii.org/ - dealii.git/commitdiff
(Dynamic)SparsityPattern: inherit from SparsityPatternBase. 14396/head
authorDavid Wells <drwells@email.unc.edu>
Thu, 13 Oct 2022 17:58:34 +0000 (13:58 -0400)
committerDavid Wells <drwells@email.unc.edu>
Sun, 6 Nov 2022 12:40:13 +0000 (07:40 -0500)
include/deal.II/lac/dynamic_sparsity_pattern.h
include/deal.II/lac/sparsity_pattern.h
include/deal.II/lac/sparsity_pattern_base.h [new file with mode: 0644]
source/lac/dynamic_sparsity_pattern.cc
source/lac/sparsity_pattern.cc

index f260aaa24575e90590e2384c646fc6849becf7c6..d3634254312cb6309853589b5ec8da142aeab138 100644 (file)
@@ -24,6 +24,7 @@
 #include <deal.II/base/utilities.h>
 
 #include <deal.II/lac/exceptions.h>
+#include <deal.II/lac/sparsity_pattern_base.h>
 
 #include <algorithm>
 #include <iostream>
@@ -318,7 +319,7 @@ namespace DynamicSparsityPatternIterators
  * sp.copy_from (dynamic_pattern);
  * @endcode
  */
-class DynamicSparsityPattern : public Subscriptor
+class DynamicSparsityPattern : public SparsityPatternBase
 {
 public:
   /**
@@ -439,6 +440,15 @@ public:
               ForwardIterator end,
               const bool      indices_are_unique_and_sorted = false);
 
+  virtual void
+  add_row_entries(const size_type &                 row,
+                  const ArrayView<const size_type> &columns,
+                  const bool indices_are_sorted = false) override;
+
+  virtual void
+  add_entries(const ArrayView<const size_type> &rows,
+              const ArrayView<const size_type> &columns) override;
+
   /**
    * Check if a value at a certain position may be non-zero.
    */
@@ -506,19 +516,6 @@ public:
   void
   print_gnuplot(std::ostream &out) const;
 
-  /**
-   * Return the number of rows, which equals the dimension of the image space.
-   */
-  size_type
-  n_rows() const;
-
-  /**
-   * Return the number of columns, which equals the dimension of the range
-   * space.
-   */
-  size_type
-  n_cols() const;
-
   /**
    * Number of entries in a specific row. This function can only be called if
    * the given row is a member of the index set of rows that we want to store.
@@ -675,16 +672,6 @@ private:
    */
   bool have_entries;
 
-  /**
-   * Number of rows that this sparsity structure shall represent.
-   */
-  size_type rows;
-
-  /**
-   * Number of columns that this sparsity structure shall represent.
-   */
-  size_type cols;
-
   /**
    * A set that contains the valid rows.
    */
@@ -1014,27 +1001,11 @@ DynamicSparsityPattern::Line::add(const size_type j)
 
 
 
-inline DynamicSparsityPattern::size_type
-DynamicSparsityPattern::n_rows() const
-{
-  return rows;
-}
-
-
-
-inline types::global_dof_index
-DynamicSparsityPattern::n_cols() const
-{
-  return cols;
-}
-
-
-
 inline void
 DynamicSparsityPattern::add(const size_type i, const size_type j)
 {
-  AssertIndexRange(i, rows);
-  AssertIndexRange(j, cols);
+  AssertIndexRange(i, n_rows());
+  AssertIndexRange(j, n_cols());
 
   if (rowset.size() > 0 && !rowset.is_element(i))
     return;
index adc2459e6b4a64ea62b4acfd1de441b4f008a004..78d907d43143044b5b37f417df35c90c36fc2ddd 100644 (file)
@@ -24,6 +24,8 @@
 #include <deal.II/base/linear_index_iterator.h>
 #include <deal.II/base/subscriptor.h>
 
+#include <deal.II/lac/sparsity_pattern_base.h>
+
 // boost::serialization::make_array used to be in array.hpp, but was
 // moved to a different file in BOOST 1.64
 #include <boost/version.hpp>
@@ -338,7 +340,7 @@ namespace SparsityPatternIterators
  * SparsityPatternIterators::Accessor, like <tt>it->column()</tt> and
  * <tt>it->row()</tt>.
  */
-class SparsityPattern : public Subscriptor
+class SparsityPattern : public SparsityPatternBase
 {
 public:
   /**
@@ -728,6 +730,16 @@ public:
               ForwardIterator end,
               const bool      indices_are_sorted = false);
 
+  virtual void
+  add_row_entries(const size_type &                 row,
+                  const ArrayView<const size_type> &columns,
+                  const bool indices_are_sorted = false) override;
+
+  virtual void
+  add_entries(const ArrayView<const size_type> &rows,
+              const ArrayView<const size_type> &columns) override;
+
+
   /**
    * Make the sparsity pattern symmetric by adding the sparsity pattern of the
    * transpose object.
@@ -761,20 +773,6 @@ public:
   std::size_t
   n_nonzero_elements() const;
 
-  /**
-   * Return number of rows of this matrix, which equals the dimension of the
-   * image space.
-   */
-  size_type
-  n_rows() const;
-
-  /**
-   * Return number of columns of this matrix, which equals the dimension of
-   * the range space.
-   */
-  size_type
-  n_cols() const;
-
   /**
    * Return whether the object is empty. It is empty if no memory is
    * allocated, which is the same as that both dimensions are zero.
@@ -1096,16 +1094,6 @@ private:
    */
   size_type max_dim;
 
-  /**
-   * Number of rows that this sparsity structure shall represent.
-   */
-  size_type rows;
-
-  /**
-   * Number of columns that this sparsity structure shall represent.
-   */
-  size_type cols;
-
   /**
    * Size of the actually allocated array #colnums. Here, the same applies as
    * for the #rowstart array, i.e. it may be larger than the actually used
@@ -1331,22 +1319,6 @@ SparsityPattern::n_nonzero_elements() const
 
 
 
-inline SparsityPattern::size_type
-SparsityPattern::n_rows() const
-{
-  return rows;
-}
-
-
-
-inline SparsityPattern::size_type
-SparsityPattern::n_cols() const
-{
-  return cols;
-}
-
-
-
 inline bool
 SparsityPattern::is_compressed() const
 {
diff --git a/include/deal.II/lac/sparsity_pattern_base.h b/include/deal.II/lac/sparsity_pattern_base.h
new file mode 100644 (file)
index 0000000..8fcf366
--- /dev/null
@@ -0,0 +1,190 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2022 - 2022 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_sparsity_pattern_base_h
+#define dealii_sparsity_pattern_base_h
+
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/base/array_view.h>
+#include <deal.II/base/exceptions.h>
+#include <deal.II/base/subscriptor.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+/**
+ * @addtogroup Sparsity
+ * @{
+ */
+
+/**
+ * Base class for all sparsity patterns, defining a common interface by which
+ * new values can be added.
+ */
+class SparsityPatternBase : public Subscriptor
+{
+public:
+  /**
+   * Declare type for container size.
+   */
+  using size_type = types::global_dof_index;
+
+  /**
+   * Constructor. Sets up an empty (zero-by-zero) sparsity pattern.
+   */
+  SparsityPatternBase();
+
+  /**
+   * Constructor. Sets up a @p rows by @p cols sparsity pattern.
+   */
+  SparsityPatternBase(const size_type rows, const size_type cols);
+
+  /**
+   * Copy constructor.
+   */
+  SparsityPatternBase(const SparsityPatternBase &sparsity_pattern) = default;
+
+  /**
+   * Move constructor.
+   */
+  SparsityPatternBase(SparsityPatternBase &&sparsity_pattern) noexcept =
+    default;
+
+  /**
+   * Assignment operator.
+   */
+  SparsityPatternBase &
+  operator=(const SparsityPatternBase &sparsity_pattern) = default;
+
+  /**
+   * Move assignment operator.
+   */
+  SparsityPatternBase &
+  operator=(SparsityPatternBase &&sparsity_pattern) noexcept = default;
+
+  /**
+   * Return number of rows of this matrix, which equals the dimension of the
+   * image space.
+   */
+  size_type
+  n_rows() const;
+
+  /**
+   * Return number of columns of this matrix, which equals the dimension of
+   * the range space.
+   */
+  size_type
+  n_cols() const;
+
+  /**
+   * Optimized function for adding new entries to a given row.
+   */
+  virtual void
+  add_row_entries(const size_type &                 row,
+                  const ArrayView<const size_type> &columns,
+                  const bool indices_are_sorted = false) = 0;
+
+  /**
+   * General function for adding new entries. By default this function calls
+   * add_row_entries() for each new entry: inheriting classes may want to add a
+   * more optimized implementation.
+   */
+  virtual void
+  add_entries(const ArrayView<const size_type> &rows,
+              const ArrayView<const size_type> &columns);
+
+protected:
+  /**
+   * Internal function for updating the stored size of the sparsity pattern.
+   */
+  virtual void
+  resize(const size_type rows, const size_type cols);
+
+  /**
+   * Number of rows that this sparsity pattern shall represent.
+   */
+  size_type rows;
+
+  /**
+   * Number of columns that this sparsity pattern shall represent.
+   */
+  size_type cols;
+};
+
+/**
+ * @}
+ */
+
+
+/* ---------------------------- Inline functions ---------------------------- */
+
+#ifndef DOXYGEN
+
+inline SparsityPatternBase::SparsityPatternBase()
+  : rows(0)
+  , cols(0)
+{}
+
+
+
+inline SparsityPatternBase::SparsityPatternBase(const size_type rows,
+                                                const size_type cols)
+  : rows(rows)
+  , cols(cols)
+{}
+
+
+
+inline SparsityPatternBase::size_type
+SparsityPatternBase::n_rows() const
+{
+  return rows;
+}
+
+
+
+inline SparsityPatternBase::size_type
+SparsityPatternBase::n_cols() const
+{
+  return cols;
+}
+
+
+
+inline void
+SparsityPatternBase::add_entries(const ArrayView<const size_type> &rows,
+                                 const ArrayView<const size_type> &columns)
+{
+  AssertDimension(rows.size(), columns.size());
+  for (std::size_t i = 0; i < rows.size(); ++i)
+    add_row_entries(rows[i],
+                    make_array_view(columns.begin() + i,
+                                    columns.begin() + i + 1));
+}
+
+
+
+inline void
+SparsityPatternBase::resize(const size_type rows, const size_type cols)
+{
+  this->rows = rows;
+  this->cols = cols;
+}
+#endif
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index 3ed8948833615e2a52d304183260fa2a3b40e940..0dabe209fe51eb73399ce91bbba9d192ca8ba658 100644 (file)
@@ -220,19 +220,16 @@ DynamicSparsityPattern::Line::memory_consumption() const
 
 
 DynamicSparsityPattern::DynamicSparsityPattern()
-  : have_entries(false)
-  , rows(0)
-  , cols(0)
+  : SparsityPatternBase()
+  , have_entries(false)
   , rowset(0)
 {}
 
 
 
 DynamicSparsityPattern::DynamicSparsityPattern(const DynamicSparsityPattern &s)
-  : Subscriptor()
+  : SparsityPatternBase()
   , have_entries(false)
-  , rows(0)
-  , cols(0)
   , rowset(0)
 {
   (void)s;
@@ -248,9 +245,8 @@ DynamicSparsityPattern::DynamicSparsityPattern(const DynamicSparsityPattern &s)
 DynamicSparsityPattern::DynamicSparsityPattern(const size_type m,
                                                const size_type n,
                                                const IndexSet &rowset_)
-  : have_entries(false)
-  , rows(0)
-  , cols(0)
+  : SparsityPatternBase()
+  , have_entries(false)
   , rowset(0)
 {
   reinit(m, n, rowset_);
@@ -258,9 +254,8 @@ DynamicSparsityPattern::DynamicSparsityPattern(const size_type m,
 
 
 DynamicSparsityPattern::DynamicSparsityPattern(const IndexSet &rowset_)
-  : have_entries(false)
-  , rows(0)
-  , cols(0)
+  : SparsityPatternBase()
+  , have_entries(false)
   , rowset(0)
 {
   reinit(rowset_.size(), rowset_.size(), rowset_);
@@ -268,9 +263,8 @@ DynamicSparsityPattern::DynamicSparsityPattern(const IndexSet &rowset_)
 
 
 DynamicSparsityPattern::DynamicSparsityPattern(const size_type n)
-  : have_entries(false)
-  , rows(0)
-  , cols(0)
+  : SparsityPatternBase()
+  , have_entries(false)
   , rowset(0)
 {
   reinit(n, n);
@@ -282,13 +276,13 @@ DynamicSparsityPattern &
 DynamicSparsityPattern::operator=(const DynamicSparsityPattern &s)
 {
   (void)s;
-  Assert(s.rows == 0 && s.cols == 0,
+  Assert(s.n_rows() == 0 && s.n_cols() == 0,
          ExcMessage(
            "This operator can only be called if the provided argument "
            "is the sparsity pattern for an empty matrix. This operator can "
            "not be used to copy a non-empty sparsity pattern."));
 
-  Assert(rows == 0 && cols == 0,
+  Assert(n_rows() == 0 && n_cols() == 0,
          ExcMessage("This operator can only be called if the current object is "
                     "empty."));
 
@@ -302,9 +296,8 @@ DynamicSparsityPattern::reinit(const size_type m,
                                const size_type n,
                                const IndexSet &rowset_)
 {
+  resize(m, n);
   have_entries = false;
-  rows         = m;
-  cols         = n;
   rowset       = rowset_;
 
   Assert(rowset.size() == 0 || rowset.size() == m,
@@ -316,7 +309,8 @@ DynamicSparsityPattern::reinit(const size_type m,
            "of indices in this IndexSet may be less than the number "
            "of rows, but the *size* of the IndexSet must be equal.)"));
 
-  std::vector<Line> new_lines(rowset.size() == 0 ? rows : rowset.n_elements());
+  std::vector<Line> new_lines(rowset.size() == 0 ? n_rows() :
+                                                   rowset.n_elements());
   lines.swap(new_lines);
 }
 
@@ -353,6 +347,28 @@ DynamicSparsityPattern::max_entries_per_row() const
 
 
 
+void
+DynamicSparsityPattern::add_row_entries(
+  const size_type &                 row,
+  const ArrayView<const size_type> &columns,
+  const bool                        indices_are_sorted)
+{
+  add_entries(row, columns.begin(), columns.end(), indices_are_sorted);
+}
+
+
+
+void
+DynamicSparsityPattern::add_entries(const ArrayView<const size_type> &rows,
+                                    const ArrayView<const size_type> &columns)
+{
+  AssertDimension(rows.size(), columns.size());
+  for (std::size_t i = 0; i < rows.size(); ++i)
+    add(rows[i], columns[i]);
+}
+
+
+
 bool
 DynamicSparsityPattern::exists(const size_type i, const size_type j) const
 {
index d6958cbbfdb8b0bd20c6dbeaa07cc6aa4c534988..10f0add3e55e8f2a756df7df663c445a844542ab 100644 (file)
@@ -36,10 +36,9 @@ constexpr SparsityPattern::size_type SparsityPattern::invalid_entry;
 
 
 SparsityPattern::SparsityPattern()
-  : store_diagonal_first_in_row(false)
+  : SparsityPatternBase()
+  , store_diagonal_first_in_row(false)
   , max_dim(0)
-  , rows(0)
-  , cols(0)
   , max_vec_len(0)
   , max_row_length(0)
   , compressed(false)
@@ -107,10 +106,10 @@ SparsityPattern::SparsityPattern(const SparsityPattern &original,
                                  const size_type        extra_off_diagonals)
   : SparsityPattern()
 {
-  Assert(original.rows == original.cols, ExcNotQuadratic());
+  Assert(original.n_rows() == original.n_cols(), ExcNotQuadratic());
   Assert(original.is_compressed(), ExcNotCompressed());
 
-  reinit(original.rows, original.cols, max_per_row);
+  reinit(original.n_rows(), original.n_cols(), max_per_row);
 
   // now copy the entries from the other object
   for (size_type row = 0; row < original.rows; ++row)
@@ -206,9 +205,7 @@ SparsityPattern::reinit(const size_type                      m,
                         const ArrayView<const unsigned int> &row_lengths)
 {
   AssertDimension(row_lengths.size(), m);
-
-  rows = m;
-  cols = n;
+  resize(m, n);
 
   // delete empty matrices
   if ((m == 0) || (n == 0))
@@ -216,7 +213,7 @@ SparsityPattern::reinit(const size_type                      m,
       rowstart.reset();
       colnums.reset();
 
-      max_vec_len = max_dim = rows = cols = 0;
+      max_vec_len = max_dim = 0;
       // if dimension is zero: ignore max_per_row
       max_row_length = 0;
       compressed     = false;
@@ -301,7 +298,7 @@ SparsityPattern::reinit(const size_type                      m,
   // if diagonal elements are special: let the first entry in each row be the
   // diagonal value
   if (store_diagonal_first_in_row)
-    for (size_type i = 0; i < rows; ++i)
+    for (size_type i = 0; i < n_rows(); ++i)
       colnums[rowstart[i]] = i;
 
   compressed = false;
@@ -361,7 +358,7 @@ SparsityPattern::compress()
   std::vector<size_type> tmp_entries(max_row_length);
 
   // Traverse all rows
-  for (size_type line = 0; line < rows; ++line)
+  for (size_type line = 0; line < n_rows(); ++line)
     {
       // copy used entries, break if first unused entry is reached
       row_length = 0;
@@ -616,11 +613,11 @@ SparsityPattern::empty() const
   // and freeing memory was not present in the original implementation and I
   // don't know at how many places I missed something in adding it, so I try
   // to be cautious. wb)
-  if ((rowstart == nullptr) || (rows == 0) || (cols == 0))
+  if ((rowstart == nullptr) || (n_rows() == 0) || (n_cols() == 0))
     {
       Assert(rowstart == nullptr, ExcInternalError());
-      Assert(rows == 0, ExcInternalError());
-      Assert(cols == 0, ExcInternalError());
+      Assert(n_rows() == 0, ExcInternalError());
+      Assert(n_cols() == 0, ExcInternalError());
       Assert(colnums == nullptr, ExcInternalError());
       Assert(max_vec_len == 0, ExcInternalError());
 
@@ -635,8 +632,8 @@ bool
 SparsityPattern::exists(const size_type i, const size_type j) const
 {
   Assert((rowstart != nullptr) && (colnums != nullptr), ExcEmptyObject());
-  AssertIndexRange(i, rows);
-  AssertIndexRange(j, cols);
+  AssertIndexRange(i, n_rows());
+  AssertIndexRange(j, n_cols());
 
   for (size_type k = rowstart[i]; k < rowstart[i + 1]; ++k)
     {
@@ -678,7 +675,7 @@ SparsityPattern::bandwidth() const
 {
   Assert((rowstart != nullptr) && (colnums != nullptr), ExcEmptyObject());
   size_type b = 0;
-  for (size_type i = 0; i < rows; ++i)
+  for (size_type i = 0; i < n_rows(); ++i)
     for (size_type j = rowstart[i]; j < rowstart[i + 1]; ++j)
       if (colnums[j] != invalid_entry)
         {
@@ -717,8 +714,8 @@ SparsityPattern::size_type
 SparsityPattern::operator()(const size_type i, const size_type j) const
 {
   Assert((rowstart != nullptr) && (colnums != nullptr), ExcEmptyObject());
-  AssertIndexRange(i, rows);
-  AssertIndexRange(j, cols);
+  AssertIndexRange(i, n_rows());
+  AssertIndexRange(j, n_cols());
   Assert(compressed, ExcNotCompressed());
 
   // let's see whether there is something in this line
@@ -755,8 +752,8 @@ void
 SparsityPattern::add(const size_type i, const size_type j)
 {
   Assert((rowstart != nullptr) && (colnums != nullptr), ExcEmptyObject());
-  AssertIndexRange(i, rows);
-  AssertIndexRange(j, cols);
+  AssertIndexRange(i, n_rows());
+  AssertIndexRange(j, n_cols());
   Assert(compressed == false, ExcMatrixIsCompressed());
 
   for (std::size_t k = rowstart[i]; k < rowstart[i + 1]; ++k)
@@ -829,6 +826,27 @@ SparsityPattern::add_entries(const size_type row,
 
 
 
+void
+SparsityPattern::add_row_entries(const size_type &                 row,
+                                 const ArrayView<const size_type> &columns,
+                                 const bool indices_are_sorted)
+{
+  add_entries(row, columns.begin(), columns.end(), indices_are_sorted);
+}
+
+
+
+void
+SparsityPattern::add_entries(const ArrayView<const size_type> &rows,
+                             const ArrayView<const size_type> &columns)
+{
+  AssertDimension(rows.size(), columns.size());
+  for (std::size_t i = 0; i < rows.size(); ++i)
+    add(rows[i], columns[i]);
+}
+
+
+
 void
 SparsityPattern::symmetrize()
 {
@@ -836,7 +854,7 @@ SparsityPattern::symmetrize()
   Assert(compressed == false, ExcMatrixIsCompressed());
   // Note that we only require a quadratic matrix here, no special treatment
   // of diagonals
-  Assert(rows == cols, ExcNotQuadratic());
+  Assert(n_rows() == n_cols(), ExcNotQuadratic());
 
   // loop over all elements presently in the sparsity pattern and add the
   // transpose element. note:
@@ -846,7 +864,7 @@ SparsityPattern::symmetrize()
   //
   // 2. that the @p{add} function can be called on elements that already exist
   // without any harm
-  for (size_type row = 0; row < rows; ++row)
+  for (size_type row = 0; row < n_rows(); ++row)
     for (size_type k = rowstart[row]; k < rowstart[row + 1]; ++k)
       {
         // check whether we are at the end of the entries of this row. if so,
@@ -867,8 +885,8 @@ SparsityPattern::size_type
 SparsityPattern::row_position(const size_type i, const size_type j) const
 {
   Assert((rowstart != nullptr) && (colnums != nullptr), ExcEmptyObject());
-  AssertIndexRange(i, rows);
-  AssertIndexRange(j, cols);
+  AssertIndexRange(i, n_rows());
+  AssertIndexRange(j, n_cols());
 
   for (size_type k = rowstart[i]; k < rowstart[i + 1]; ++k)
     {
@@ -888,7 +906,7 @@ SparsityPattern::print(std::ostream &out) const
 
   AssertThrow(out.fail() == false, ExcIO());
 
-  for (size_type i = 0; i < rows; ++i)
+  for (size_type i = 0; i < n_rows(); ++i)
     {
       out << '[' << i;
       for (size_type j = rowstart[i]; j < rowstart[i + 1]; ++j)
@@ -909,7 +927,7 @@ SparsityPattern::print_gnuplot(std::ostream &out) const
 
   AssertThrow(out.fail() == false, ExcIO());
 
-  for (size_type i = 0; i < rows; ++i)
+  for (size_type i = 0; i < n_rows(); ++i)
     for (size_type j = rowstart[i]; j < rowstart[i + 1]; ++j)
       if (colnums[j] != invalid_entry)
         // while matrix entries are usually written (i,j), with i vertical and
@@ -961,8 +979,8 @@ SparsityPattern::block_write(std::ostream &out) const
   AssertThrow(out.fail() == false, ExcIO());
 
   // first the simple objects, bracketed in [...]
-  out << '[' << max_dim << ' ' << rows << ' ' << cols << ' ' << max_vec_len
-      << ' ' << max_row_length << ' ' << compressed << ' '
+  out << '[' << max_dim << ' ' << n_rows() << ' ' << n_cols() << ' '
+      << max_vec_len << ' ' << max_row_length << ' ' << compressed << ' '
       << store_diagonal_first_in_row << "][";
   // then write out real data
   out.write(reinterpret_cast<const char *>(rowstart.get()),

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.