]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Redo SparsityPatternBase::add_entries() to use pairs. 14441/head
authorDavid Wells <drwells@email.unc.edu>
Mon, 21 Nov 2022 20:47:59 +0000 (15:47 -0500)
committerDavid Wells <drwells@email.unc.edu>
Mon, 21 Nov 2022 21:30:18 +0000 (16:30 -0500)
15 files changed:
include/deal.II/lac/affine_constraints.h
include/deal.II/lac/affine_constraints.templates.h
include/deal.II/lac/block_sparsity_pattern.h
include/deal.II/lac/dynamic_sparsity_pattern.h
include/deal.II/lac/sparsity_pattern.h
include/deal.II/lac/sparsity_pattern_base.h
include/deal.II/lac/trilinos_sparsity_pattern.h
source/dofs/dof_tools_sparsity.cc
source/lac/CMakeLists.txt
source/lac/dynamic_sparsity_pattern.cc
source/lac/sparsity_pattern.cc
source/lac/sparsity_pattern_base.cc [new file with mode: 0644]
source/lac/trilinos_sparsity_pattern.cc
tests/sparsity/sparsity_pattern_base_01.cc [new file with mode: 0644]
tests/sparsity/sparsity_pattern_base_01.output [new file with mode: 0644]

index cb640ff5cfac48e1f4d4d61ee0f09a30bff0512c..816a61f6f22eba55129c1effb2388b4ce49267ac 100644 (file)
@@ -358,6 +358,11 @@ namespace internal
        */
       bool in_use;
 
+      /**
+       * Temporary array for pairs of indices
+       */
+      std::vector<std::pair<size_type, size_type>> new_entries;
+
       /**
        * Temporary array for row indices
        */
index babfd8ef061028328833b4ec81cca257d1319383..584ceb39704b496fcbe4fac35a3520b501d785de 100644 (file)
@@ -3703,12 +3703,10 @@ namespace internal
       // that have been left out above
       if (global_rows.n_constraints() > 0)
         {
-          scratch_data.rows.resize(0);
-          scratch_data.rows.reserve(global_rows.n_constraints() *
-                                    local_dof_indices.size());
-          scratch_data.columns.resize(0);
-          scratch_data.columns.reserve(global_rows.n_constraints() *
-                                       local_dof_indices.size());
+          std::vector<std::pair<size_type, size_type>> &cell_entries =
+            scratch_data.new_entries;
+          cell_entries.resize(0);
+          cell_entries.reserve(local_dof_indices.size());
           for (size_type i = 0; i < global_rows.n_constraints(); ++i)
             {
               const size_type local_row  = global_rows.constraint_origin(i);
@@ -3718,26 +3716,18 @@ namespace internal
                   for (size_type j = 0; j < local_dof_indices.size(); ++j)
                     {
                       if (dof_mask(local_row, j) == true)
-                        {
-                          scratch_data.rows.push_back(global_row);
-                          scratch_data.columns.push_back(local_dof_indices[j]);
-                        }
+                        cell_entries.emplace_back(global_row,
+                                                  local_dof_indices[j]);
                       if (dof_mask(j, local_row) == true)
-                        {
-                          scratch_data.rows.push_back(local_dof_indices[j]);
-                          scratch_data.columns.push_back(global_row);
-                        }
+                        cell_entries.emplace_back(local_dof_indices[j],
+                                                  global_row);
                     }
                 }
               else
-                {
-                  // don't keep constrained entries - just add the diagonal.
-                  scratch_data.rows.push_back(global_row);
-                  scratch_data.columns.push_back(global_row);
-                }
+                // don't keep constrained entries - just add the diagonal.
+                cell_entries.emplace_back(global_row, global_row);
             }
-          sparsity_pattern.add_entries(make_array_view(scratch_data.rows),
-                                       make_array_view(scratch_data.columns));
+          sparsity_pattern.add_entries(make_array_view(cell_entries));
         }
     }
 
@@ -4331,32 +4321,29 @@ AffineConstraints<number>::add_entries_local_to_global(
                                          true);
 
       // need to add the whole row and column structure in case we keep
-      // constrained entries. Unfortunately, we can't use the nice matrix
-      // structure we use elsewhere, so use the more general update function
-      actual_dof_indices.resize(0);
-      actual_dof_indices.reserve(n_local_dofs * n_local_dofs);
-      std::vector<size_type> &rows = scratch_data->rows;
-      rows.resize(0);
-      rows.reserve(n_local_dofs * n_local_dofs);
+      // constrained entries.
+      std::vector<std::pair<size_type, size_type>> &cell_entries =
+        scratch_data->new_entries;
+      cell_entries.resize(0);
+      cell_entries.reserve(n_local_dofs);
       for (size_type i = 0; i < n_local_dofs; ++i)
         if (is_constrained(local_dof_indices[i]))
           {
             if (keep_constrained_entries == true)
               for (size_type j = 0; j < n_local_dofs; ++j)
                 {
-                  rows.push_back(local_dof_indices[i]);
-                  actual_dof_indices.push_back(local_dof_indices[j]);
-                  rows.push_back(local_dof_indices[j]);
-                  actual_dof_indices.push_back(local_dof_indices[i]);
+                  cell_entries.emplace_back(local_dof_indices[i],
+                                            local_dof_indices[j]);
+                  cell_entries.emplace_back(local_dof_indices[j],
+                                            local_dof_indices[i]);
                 }
             else
               {
-                rows.push_back(local_dof_indices[i]);
-                actual_dof_indices.push_back(local_dof_indices[i]);
+                cell_entries.emplace_back(local_dof_indices[i],
+                                          local_dof_indices[i]);
               }
           }
-      sparsity_pattern.add_entries(make_array_view(rows),
-                                   make_array_view(actual_dof_indices));
+      sparsity_pattern.add_entries(make_array_view(cell_entries));
 
       return;
     }
@@ -4413,13 +4400,11 @@ AffineConstraints<number>::add_entries_local_to_global(
   const size_type n_local_cols = col_indices.size();
 
   typename internal::AffineConstraints::ScratchDataAccessor<number>
-                          scratch_data(this->scratch_data);
-  std::vector<size_type> &rows = scratch_data->rows;
-  std::vector<size_type> &cols = scratch_data->columns;
-  rows.resize(0);
-  rows.reserve(row_indices.size() * col_indices.size());
-  cols.resize(0);
-  cols.reserve(row_indices.size() * col_indices.size());
+    scratch_data(this->scratch_data);
+  std::vector<std::pair<size_type, size_type>> &cell_entries =
+    scratch_data->new_entries;
+  cell_entries.resize(0);
+  cell_entries.reserve(row_indices.size() * col_indices.size());
 
   // if constrained entries should be kept, need to add rows and columns of
   // those to the sparsity pattern
@@ -4428,19 +4413,12 @@ AffineConstraints<number>::add_entries_local_to_global(
       for (const size_type row_index : row_indices)
         if (is_constrained(row_index))
           for (const size_type col_index : col_indices)
-            {
-              rows.push_back(row_index);
-              cols.push_back(col_index);
-            }
+            cell_entries.emplace_back(row_index, col_index);
       for (const size_type col_index : col_indices)
         if (col_constraints.is_constrained(col_index))
           for (const size_type row_index : row_indices)
-            {
-              rows.push_back(row_index);
-              cols.push_back(col_index);
-            }
-      sparsity_pattern.add_entries(make_array_view(rows),
-                                   make_array_view(cols));
+            cell_entries.emplace_back(row_index, col_index);
+      sparsity_pattern.add_entries(make_array_view(cell_entries));
     }
 
   // if the dof mask is not active, all we have to do is to add some indices
@@ -4451,6 +4429,8 @@ AffineConstraints<number>::add_entries_local_to_global(
     dof_mask.n_rows() == n_local_rows && dof_mask.n_cols() == n_local_cols;
   if (dof_mask_is_active == false)
     {
+      std::vector<size_type> &rows = scratch_data->rows;
+      std::vector<size_type> &cols = scratch_data->columns;
       rows.resize(n_local_rows);
       cols.resize(n_local_cols);
       // TODO these fills may not be necessary: previously we assumed all zeros
index aac3dd42c4b8f16dd2c61df71dcb6bf20b354e49..156dc64b7679763f6614c1d0fa7d2b64dcc40fe0 100644 (file)
@@ -253,9 +253,7 @@ public:
                   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;
+  using SparsityPatternBase::add_entries;
 
   /**
    * Return number of rows of this matrix, which equals the dimension of the
@@ -996,19 +994,6 @@ BlockSparsityPatternBase<SparsityPatternType>::add_row_entries(
 
 
 
-template <typename SparsityPatternType>
-inline void
-BlockSparsityPatternBase<SparsityPatternType>::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]);
-}
-
-
-
 template <typename SparsityPatternType>
 inline bool
 BlockSparsityPatternBase<SparsityPatternType>::exists(const size_type i,
index d3634254312cb6309853589b5ec8da142aeab138..7b032fbe9491e9d62e8dc89f0fa3f04b4be05031 100644 (file)
@@ -445,9 +445,7 @@ public:
                   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;
+  using SparsityPatternBase::add_entries;
 
   /**
    * Check if a value at a certain position may be non-zero.
index 78d907d43143044b5b37f417df35c90c36fc2ddd..7e7985e9bf32a0e113540149430018ae9b423829 100644 (file)
@@ -735,10 +735,7 @@ public:
                   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;
-
+  using SparsityPatternBase::add_entries;
 
   /**
    * Make the sparsity pattern symmetric by adding the sparsity pattern of the
index 8fcf366567d58dc4a76e11e70e418289db47b304..d8fa38325bf7db099f0c0be5bccfb54cd1649471 100644 (file)
@@ -23,6 +23,8 @@
 #include <deal.II/base/exceptions.h>
 #include <deal.II/base/subscriptor.h>
 
+#include <utility>
+
 DEAL_II_NAMESPACE_OPEN
 
 /**
@@ -98,13 +100,14 @@ public:
                   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.
+   * General function for adding new entries from an unsorted list of pairs.
+   *
+   * This function is useful when multiple entries need to be added which do
+   * not correspond to a particular row, e.g., when assembling a flux sparsity
+   * pattern.
    */
   virtual void
-  add_entries(const ArrayView<const size_type> &rows,
-              const ArrayView<const size_type> &columns);
+  add_entries(const ArrayView<const std::pair<size_type, size_type>> &entries);
 
 protected:
   /**
@@ -164,19 +167,6 @@ SparsityPatternBase::n_cols() const
 
 
 
-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)
 {
index 4c6e575106f19f4dba9c4f56a6f5727c27dc0c48..0f918eea89b898e7d1ee18f03bd420ee5ebf6741 100644 (file)
@@ -778,9 +778,7 @@ namespace TrilinosWrappers
                     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;
+    using SparsityPatternBase::add_entries;
 
     /** @} */
     /**
index 68a2be4ac8223bf4e0e516c6e4c518112e3b6202..310ee84492f33675337838cb191b6e15ed2bcbd0 100644 (file)
@@ -774,13 +774,9 @@ namespace DoFTools
                const unsigned int)> &face_has_flux_coupling)
       {
         std::vector<types::global_dof_index> rows;
-        std::vector<types::global_dof_index> cols;
-
-        auto append_pair = [&](const types::global_dof_index i,
-                               const types::global_dof_index j) {
-          rows.push_back(i);
-          cols.push_back(j);
-        };
+        std::vector<std::pair<SparsityPatternBase::size_type,
+                              SparsityPatternBase::size_type>>
+          cell_entries;
 
         if (dof.has_hp_capabilities() == false)
           {
@@ -849,14 +845,13 @@ namespace DoFTools
                                   if (flux_dof_mask(i, j) == always ||
                                       (flux_dof_mask(i, j) == nonzero &&
                                        i_non_zero_i && j_non_zero_i))
-                                    append_pair(dofs_on_this_cell[i],
-                                                dofs_on_this_cell[j]);
+                                    cell_entries.emplace_back(
+                                      dofs_on_this_cell[i],
+                                      dofs_on_this_cell[j]);
                                 }
                             }
-                          sparsity.add_entries(make_array_view(rows),
-                                               make_array_view(cols));
-                          rows.clear();
-                          cols.clear();
+                          sparsity.add_entries(make_array_view(cell_entries));
+                          cell_entries.clear();
                         }
                       else
                         {
@@ -948,14 +943,16 @@ namespace DoFTools
 
                                           if (flux_dof_mask(i, j) == always)
                                             {
-                                              append_pair(
+                                              cell_entries.emplace_back(
                                                 dofs_on_this_cell[i],
                                                 dofs_on_other_cell[j]);
-                                              append_pair(dofs_on_other_cell[i],
-                                                          dofs_on_this_cell[j]);
-                                              append_pair(dofs_on_this_cell[i],
-                                                          dofs_on_this_cell[j]);
-                                              append_pair(
+                                              cell_entries.emplace_back(
+                                                dofs_on_other_cell[i],
+                                                dofs_on_this_cell[j]);
+                                              cell_entries.emplace_back(
+                                                dofs_on_this_cell[i],
+                                                dofs_on_this_cell[j]);
+                                              cell_entries.emplace_back(
                                                 dofs_on_other_cell[i],
                                                 dofs_on_other_cell[j]);
                                             }
@@ -963,33 +960,35 @@ namespace DoFTools
                                                    nonzero)
                                             {
                                               if (i_non_zero_i && j_non_zero_e)
-                                                append_pair(
+                                                cell_entries.emplace_back(
                                                   dofs_on_this_cell[i],
                                                   dofs_on_other_cell[j]);
                                               if (i_non_zero_e && j_non_zero_i)
-                                                append_pair(
+                                                cell_entries.emplace_back(
                                                   dofs_on_other_cell[i],
                                                   dofs_on_this_cell[j]);
                                               if (i_non_zero_i && j_non_zero_i)
-                                                append_pair(
+                                                cell_entries.emplace_back(
                                                   dofs_on_this_cell[i],
                                                   dofs_on_this_cell[j]);
                                               if (i_non_zero_e && j_non_zero_e)
-                                                append_pair(
+                                                cell_entries.emplace_back(
                                                   dofs_on_other_cell[i],
                                                   dofs_on_other_cell[j]);
                                             }
 
                                           if (flux_dof_mask(j, i) == always)
                                             {
-                                              append_pair(
+                                              cell_entries.emplace_back(
                                                 dofs_on_this_cell[j],
                                                 dofs_on_other_cell[i]);
-                                              append_pair(dofs_on_other_cell[j],
-                                                          dofs_on_this_cell[i]);
-                                              append_pair(dofs_on_this_cell[j],
-                                                          dofs_on_this_cell[i]);
-                                              append_pair(
+                                              cell_entries.emplace_back(
+                                                dofs_on_other_cell[j],
+                                                dofs_on_this_cell[i]);
+                                              cell_entries.emplace_back(
+                                                dofs_on_this_cell[j],
+                                                dofs_on_this_cell[i]);
+                                              cell_entries.emplace_back(
                                                 dofs_on_other_cell[j],
                                                 dofs_on_other_cell[i]);
                                             }
@@ -997,29 +996,28 @@ namespace DoFTools
                                                    nonzero)
                                             {
                                               if (j_non_zero_i && i_non_zero_e)
-                                                append_pair(
+                                                cell_entries.emplace_back(
                                                   dofs_on_this_cell[j],
                                                   dofs_on_other_cell[i]);
                                               if (j_non_zero_e && i_non_zero_i)
-                                                append_pair(
+                                                cell_entries.emplace_back(
                                                   dofs_on_other_cell[j],
                                                   dofs_on_this_cell[i]);
                                               if (j_non_zero_i && i_non_zero_i)
-                                                append_pair(
+                                                cell_entries.emplace_back(
                                                   dofs_on_this_cell[j],
                                                   dofs_on_this_cell[i]);
                                               if (j_non_zero_e && i_non_zero_e)
-                                                append_pair(
+                                                cell_entries.emplace_back(
                                                   dofs_on_other_cell[j],
                                                   dofs_on_other_cell[i]);
                                             }
                                         }
                                     }
                                 }
-                              sparsity.add_entries(make_array_view(rows),
-                                                   make_array_view(cols));
-                              rows.clear();
-                              cols.clear();
+                              sparsity.add_entries(
+                                make_array_view(cell_entries));
+                              cell_entries.clear();
                             }
                           else
                             {
@@ -1041,63 +1039,78 @@ namespace DoFTools
                                         support_on_face(j, neighbor_face_n);
                                       if (flux_dof_mask(i, j) == always)
                                         {
-                                          append_pair(dofs_on_this_cell[i],
-                                                      dofs_on_other_cell[j]);
-                                          append_pair(dofs_on_other_cell[i],
-                                                      dofs_on_this_cell[j]);
-                                          append_pair(dofs_on_this_cell[i],
-                                                      dofs_on_this_cell[j]);
-                                          append_pair(dofs_on_other_cell[i],
-                                                      dofs_on_other_cell[j]);
+                                          cell_entries.emplace_back(
+                                            dofs_on_this_cell[i],
+                                            dofs_on_other_cell[j]);
+                                          cell_entries.emplace_back(
+                                            dofs_on_other_cell[i],
+                                            dofs_on_this_cell[j]);
+                                          cell_entries.emplace_back(
+                                            dofs_on_this_cell[i],
+                                            dofs_on_this_cell[j]);
+                                          cell_entries.emplace_back(
+                                            dofs_on_other_cell[i],
+                                            dofs_on_other_cell[j]);
                                         }
                                       if (flux_dof_mask(i, j) == nonzero)
                                         {
                                           if (i_non_zero_i && j_non_zero_e)
-                                            append_pair(dofs_on_this_cell[i],
-                                                        dofs_on_other_cell[j]);
+                                            cell_entries.emplace_back(
+                                              dofs_on_this_cell[i],
+                                              dofs_on_other_cell[j]);
                                           if (i_non_zero_e && j_non_zero_i)
-                                            append_pair(dofs_on_other_cell[i],
-                                                        dofs_on_this_cell[j]);
+                                            cell_entries.emplace_back(
+                                              dofs_on_other_cell[i],
+                                              dofs_on_this_cell[j]);
                                           if (i_non_zero_i && j_non_zero_i)
-                                            append_pair(dofs_on_this_cell[i],
-                                                        dofs_on_this_cell[j]);
+                                            cell_entries.emplace_back(
+                                              dofs_on_this_cell[i],
+                                              dofs_on_this_cell[j]);
                                           if (i_non_zero_e && j_non_zero_e)
-                                            append_pair(dofs_on_other_cell[i],
-                                                        dofs_on_other_cell[j]);
+                                            cell_entries.emplace_back(
+                                              dofs_on_other_cell[i],
+                                              dofs_on_other_cell[j]);
                                         }
 
                                       if (flux_dof_mask(j, i) == always)
                                         {
-                                          append_pair(dofs_on_this_cell[j],
-                                                      dofs_on_other_cell[i]);
-                                          append_pair(dofs_on_other_cell[j],
-                                                      dofs_on_this_cell[i]);
-                                          append_pair(dofs_on_this_cell[j],
-                                                      dofs_on_this_cell[i]);
-                                          append_pair(dofs_on_other_cell[j],
-                                                      dofs_on_other_cell[i]);
+                                          cell_entries.emplace_back(
+                                            dofs_on_this_cell[j],
+                                            dofs_on_other_cell[i]);
+                                          cell_entries.emplace_back(
+                                            dofs_on_other_cell[j],
+                                            dofs_on_this_cell[i]);
+                                          cell_entries.emplace_back(
+                                            dofs_on_this_cell[j],
+                                            dofs_on_this_cell[i]);
+                                          cell_entries.emplace_back(
+                                            dofs_on_other_cell[j],
+                                            dofs_on_other_cell[i]);
                                         }
                                       if (flux_dof_mask(j, i) == nonzero)
                                         {
                                           if (j_non_zero_i && i_non_zero_e)
-                                            append_pair(dofs_on_this_cell[j],
-                                                        dofs_on_other_cell[i]);
+                                            cell_entries.emplace_back(
+                                              dofs_on_this_cell[j],
+                                              dofs_on_other_cell[i]);
                                           if (j_non_zero_e && i_non_zero_i)
-                                            append_pair(dofs_on_other_cell[j],
-                                                        dofs_on_this_cell[i]);
+                                            cell_entries.emplace_back(
+                                              dofs_on_other_cell[j],
+                                              dofs_on_this_cell[i]);
                                           if (j_non_zero_i && i_non_zero_i)
-                                            append_pair(dofs_on_this_cell[j],
-                                                        dofs_on_this_cell[i]);
+                                            cell_entries.emplace_back(
+                                              dofs_on_this_cell[j],
+                                              dofs_on_this_cell[i]);
                                           if (j_non_zero_e && i_non_zero_e)
-                                            append_pair(dofs_on_other_cell[j],
-                                                        dofs_on_other_cell[i]);
+                                            cell_entries.emplace_back(
+                                              dofs_on_other_cell[j],
+                                              dofs_on_other_cell[i]);
                                         }
                                     }
                                 }
-                              sparsity.add_entries(make_array_view(rows),
-                                                   make_array_view(cols));
-                              rows.clear();
-                              cols.clear();
+                              sparsity.add_entries(
+                                make_array_view(cell_entries));
+                              cell_entries.clear();
                             }
                         }
                     }
@@ -1287,7 +1300,7 @@ namespace DoFTools
                                           if ((flux_mask(ii, jj) == always) ||
                                               (flux_mask(ii, jj) == nonzero))
                                             {
-                                              append_pair(
+                                              cell_entries.emplace_back(
                                                 dofs_on_this_cell[i],
                                                 dofs_on_other_cell[j]);
                                             }
@@ -1295,8 +1308,9 @@ namespace DoFTools
                                           if ((flux_mask(jj, ii) == always) ||
                                               (flux_mask(jj, ii) == nonzero))
                                             {
-                                              append_pair(dofs_on_other_cell[j],
-                                                          dofs_on_this_cell[i]);
+                                              cell_entries.emplace_back(
+                                                dofs_on_other_cell[j],
+                                                dofs_on_this_cell[i]);
                                             }
                                         }
                                     }
@@ -1343,25 +1357,25 @@ namespace DoFTools
                                       if ((flux_mask(ii, jj) == always) ||
                                           (flux_mask(ii, jj) == nonzero))
                                         {
-                                          append_pair(dofs_on_this_cell[i],
-                                                      dofs_on_other_cell[j]);
+                                          cell_entries.emplace_back(
+                                            dofs_on_this_cell[i],
+                                            dofs_on_other_cell[j]);
                                         }
 
                                       if ((flux_mask(jj, ii) == always) ||
                                           (flux_mask(jj, ii) == nonzero))
                                         {
-                                          append_pair(dofs_on_other_cell[j],
-                                                      dofs_on_this_cell[i]);
+                                          cell_entries.emplace_back(
+                                            dofs_on_other_cell[j],
+                                            dofs_on_this_cell[i]);
                                         }
                                     }
                                 }
                             }
                         }
                     }
-                  sparsity.add_entries(make_array_view(rows),
-                                       make_array_view(cols));
-                  rows.clear();
-                  cols.clear();
+                  sparsity.add_entries(make_array_view(cell_entries));
+                  cell_entries.clear();
                 }
           }
       }
index 57373d13d18c6aec7918f1920ed8156de858f30f..3d907f8be184a3f92b7c082ce07906b2938ee4f5 100644 (file)
@@ -42,6 +42,7 @@ SET(_unity_include_src
   sparse_matrix_ez.cc
   sparse_mic.cc
   sparse_vanka.cc
+  sparsity_pattern_base.cc
   sparsity_pattern.cc
   sparsity_tools.cc
   tensor_product_matrix.cc
index 0dabe209fe51eb73399ce91bbba9d192ca8ba658..e43e2dbc61b7da2b9caef7398a3a21b8e27896e8 100644 (file)
@@ -358,17 +358,6 @@ DynamicSparsityPattern::add_row_entries(
 
 
 
-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 8094e9098a51fd1ff64f329729359218f6a57daa..8b022eaa29ddfd7783d7637e6356d21c94795697 100644 (file)
@@ -838,17 +838,6 @@ SparsityPattern::add_row_entries(const size_type &                 row,
 
 
 
-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()
 {
diff --git a/source/lac/sparsity_pattern_base.cc b/source/lac/sparsity_pattern_base.cc
new file mode 100644 (file)
index 0000000..78d46aa
--- /dev/null
@@ -0,0 +1,64 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/lac/sparsity_pattern_base.h>
+
+DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
+#include <boost/container/small_vector.hpp>
+DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
+
+#include <algorithm>
+#include <utility>
+
+DEAL_II_NAMESPACE_OPEN
+
+void
+SparsityPatternBase::add_entries(
+  const ArrayView<const std::pair<size_type, size_type>> &inputs)
+{
+  // We always want to sort so that we can use the optimized add_row_entries()
+  // function
+  boost::container::small_vector<std::pair<size_type, size_type>, 128> entries(
+    inputs.begin(), inputs.end());
+  std::sort(entries.begin(), entries.end());
+  boost::container::small_vector<size_type, 128> columns;
+
+  auto entry = entries.begin();
+  while (entry != entries.end())
+    {
+      const auto row     = entry->first;
+      auto       row_end = entry;
+      while (row_end != entries.end() && row_end->first == row)
+        ++row_end;
+
+      columns.resize(0);
+      columns.reserve(row_end - entry);
+      // Simultaneously transform and uniquify
+      columns.push_back(entry->second);
+      ++entry;
+      while (entry != row_end)
+        {
+          if (columns.back() != entry->second)
+            columns.push_back(entry->second);
+          ++entry;
+        }
+
+      add_row_entries(row,
+                      make_array_view(columns.begin(), columns.end()),
+                      true);
+    }
+}
+
+DEAL_II_NAMESPACE_CLOSE
index fe31d1da8cd79e2387d1206e4cf731206dbde278..b84c7bc78da41baa80fda3aaa3396e75c8d9cae5 100644 (file)
@@ -990,17 +990,6 @@ namespace TrilinosWrappers
 
 
 
-  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]);
-  }
-
-
-
   const Epetra_Map &
   SparsityPattern::domain_partitioner() const
   {
diff --git a/tests/sparsity/sparsity_pattern_base_01.cc b/tests/sparsity/sparsity_pattern_base_01.cc
new file mode 100644 (file)
index 0000000..2cd55e6
--- /dev/null
@@ -0,0 +1,66 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2001 - 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.
+//
+// ---------------------------------------------------------------------
+
+
+// Test SparsityPatternBase::add_entries()
+
+
+#include <deal.II/lac/sparsity_pattern_base.h>
+
+#include "../tests.h"
+
+class TestPattern : public SparsityPatternBase
+{
+  using SparsityPatternBase::SparsityPatternBase;
+
+  using SparsityPatternBase::size_type;
+
+  virtual void
+  add_row_entries(const size_type &                 row,
+                  const ArrayView<const size_type> &columns,
+                  const bool indices_are_sorted = false) override
+  {
+    deallog << "row = " << row;
+
+    AssertThrow(std::is_sorted(columns.begin(), columns.end()),
+                ExcInternalError());
+
+    for (const auto &col : columns)
+      deallog << "    " << col;
+
+    deallog << std::endl;
+  }
+};
+
+
+
+int
+main()
+{
+  initlog();
+
+  std::vector<
+    std::pair<SparsityPatternBase::size_type, SparsityPatternBase::size_type>>
+    entries;
+
+  for (unsigned int i = 0; i < 100; ++i)
+    entries.emplace_back(rand() % 10, rand() % 100);
+
+  for (auto &entry : entries)
+    deallog << entry.first << ", " << entry.second << std::endl;
+
+  TestPattern test(100, 100);
+  test.add_entries(make_array_view(entries));
+}
diff --git a/tests/sparsity/sparsity_pattern_base_01.output b/tests/sparsity/sparsity_pattern_base_01.output
new file mode 100644 (file)
index 0000000..f70a8ae
--- /dev/null
@@ -0,0 +1,111 @@
+
+DEAL::6, 83
+DEAL::5, 77
+DEAL::5, 93
+DEAL::2, 86
+DEAL::1, 49
+DEAL::7, 62
+DEAL::9, 90
+DEAL::6, 63
+DEAL::6, 40
+DEAL::6, 72
+DEAL::8, 11
+DEAL::9, 67
+DEAL::0, 82
+DEAL::3, 62
+DEAL::5, 67
+DEAL::2, 29
+DEAL::8, 22
+DEAL::7, 69
+DEAL::6, 93
+DEAL::2, 11
+DEAL::3, 29
+DEAL::9, 21
+DEAL::7, 84
+DEAL::4, 98
+DEAL::0, 15
+DEAL::6, 13
+DEAL::0, 91
+DEAL::3, 56
+DEAL::0, 62
+DEAL::1, 96
+DEAL::5, 5
+DEAL::7, 84
+DEAL::5, 36
+DEAL::9, 46
+DEAL::7, 13
+DEAL::5, 24
+DEAL::5, 82
+DEAL::7, 14
+DEAL::4, 34
+DEAL::0, 43
+DEAL::8, 87
+DEAL::8, 76
+DEAL::4, 88
+DEAL::1, 3
+DEAL::9, 54
+DEAL::0, 32
+DEAL::8, 76
+DEAL::2, 39
+DEAL::6, 26
+DEAL::9, 94
+DEAL::0, 95
+DEAL::8, 34
+DEAL::1, 67
+DEAL::2, 97
+DEAL::2, 17
+DEAL::6, 52
+DEAL::0, 1
+DEAL::1, 86
+DEAL::9, 65
+DEAL::9, 44
+DEAL::9, 40
+DEAL::7, 31
+DEAL::1, 97
+DEAL::5, 81
+DEAL::7, 9
+DEAL::6, 67
+DEAL::3, 97
+DEAL::5, 86
+DEAL::3, 6
+DEAL::4, 19
+DEAL::1, 28
+DEAL::9, 32
+DEAL::9, 3
+DEAL::8, 70
+DEAL::5, 8
+DEAL::9, 40
+DEAL::3, 96
+DEAL::5, 18
+DEAL::1, 46
+DEAL::5, 21
+DEAL::8, 79
+DEAL::8, 64
+DEAL::0, 41
+DEAL::0, 93
+DEAL::4, 34
+DEAL::4, 24
+DEAL::6, 87
+DEAL::1, 43
+DEAL::5, 27
+DEAL::6, 59
+DEAL::1, 32
+DEAL::8, 37
+DEAL::7, 75
+DEAL::1, 74
+DEAL::5, 58
+DEAL::7, 29
+DEAL::3, 35
+DEAL::8, 18
+DEAL::1, 43
+DEAL::9, 28
+DEAL::row = 0    1    15    32    41    43    62    82    91    93    95
+DEAL::row = 1    3    28    32    43    46    49    67    74    86    96    97
+DEAL::row = 2    11    17    29    39    86    97
+DEAL::row = 3    6    29    35    56    62    96    97
+DEAL::row = 4    19    24    34    88    98
+DEAL::row = 5    5    8    18    21    24    27    36    58    67    77    81    82    86    93
+DEAL::row = 6    13    26    40    52    59    63    67    72    83    87    93
+DEAL::row = 7    9    13    14    29    31    62    69    75    84
+DEAL::row = 8    11    18    22    34    37    64    70    76    79    87
+DEAL::row = 9    3    21    28    32    40    44    46    54    65    67    90    94

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.