]> https://gitweb.dealii.org/ - dealii.git/commitdiff
multigrid tools: Remove template on sparsity pattern type. 14433/head
authorDavid Wells <drwells@email.unc.edu>
Sat, 5 Nov 2022 14:12:37 +0000 (10:12 -0400)
committerDavid Wells <drwells@email.unc.edu>
Tue, 22 Nov 2022 16:18:31 +0000 (11:18 -0500)
include/deal.II/multigrid/mg_tools.h
source/multigrid/mg_tools.cc
source/multigrid/mg_tools.inst.in

index a78c53b32c403d1bbab476b1d8fd652d99c2763c..2ea75eacfb5d76f0ac6fa46cdae4d048f1e81373 100644 (file)
@@ -76,8 +76,8 @@ namespace MGTools
   /**
    * Write the sparsity structure of the matrix belonging to the specified @p
    * level. The sparsity pattern is not compressed, so before creating the
-   * actual matrix you have to compress the matrix yourself, using
-   * <tt>SparsityPatternType::compress()</tt>.
+   * actual matrix you have to compress the matrix yourself, either using
+   * SparsityPattern::compress() or by copying to a SparsityPattern.
    *
    * The optional AffineConstraints argument allows to define constraints of
    * the level matrices like Dirichlet boundary conditions. Note that there is
@@ -85,14 +85,11 @@ namespace MGTools
    * one level is considered. See DoFTools::make_sparsity_pattern() for more
    * details about the arguments.
    */
-  template <int dim,
-            int spacedim,
-            typename SparsityPatternType,
-            typename number = double>
+  template <int dim, int spacedim, typename number = double>
   void
   make_sparsity_pattern(
     const DoFHandler<dim, spacedim> &dof_handler,
-    SparsityPatternType &            sparsity,
+    SparsityPatternBase &            sparsity,
     const unsigned int               level,
     const AffineConstraints<number> &constraints = AffineConstraints<number>(),
     const bool                       keep_constrained_dofs = true);
@@ -105,14 +102,11 @@ namespace MGTools
    * and
    * @ref DoFTools
    */
-  template <int dim,
-            int spacedim,
-            typename SparsityPatternType,
-            typename number = double>
+  template <int dim, int spacedim, typename number = double>
   void
   make_flux_sparsity_pattern(
     const DoFHandler<dim, spacedim> &dof_handler,
-    SparsityPatternType &            sparsity,
+    SparsityPatternBase &            sparsity,
     const unsigned int               level,
     const AffineConstraints<number> &constraints = AffineConstraints<number>(),
     const bool                       keep_constrained_dofs = true);
@@ -124,10 +118,10 @@ namespace MGTools
    *
    * make_flux_sparsity_pattern()
    */
-  template <int dim, typename SparsityPatternType, int spacedim>
+  template <int dim, int spacedim>
   void
   make_flux_sparsity_pattern_edge(const DoFHandler<dim, spacedim> &dof_handler,
-                                  SparsityPatternType &            sparsity,
+                                  SparsityPatternBase &            sparsity,
                                   const unsigned int               level);
   /**
    * This function does the same as the other with the same name, but it gets
@@ -138,10 +132,10 @@ namespace MGTools
    * There is one matrix for couplings in a cell and one for the couplings
    * occurring in fluxes.
    */
-  template <int dim, typename SparsityPatternType, int spacedim>
+  template <int dim, int spacedim>
   void
   make_flux_sparsity_pattern(const DoFHandler<dim, spacedim> &   dof,
-                             SparsityPatternType &               sparsity,
+                             SparsityPatternBase &               sparsity,
                              const unsigned int                  level,
                              const Table<2, DoFTools::Coupling> &int_mask,
                              const Table<2, DoFTools::Coupling> &flux_mask);
@@ -154,11 +148,11 @@ namespace MGTools
    *
    * make_flux_sparsity_pattern()
    */
-  template <int dim, typename SparsityPatternType, int spacedim>
+  template <int dim, int spacedim>
   void
   make_flux_sparsity_pattern_edge(
     const DoFHandler<dim, spacedim> &   dof_handler,
-    SparsityPatternType &               sparsity,
+    SparsityPatternBase &               sparsity,
     const unsigned int                  level,
     const Table<2, DoFTools::Coupling> &flux_mask);
 
@@ -169,11 +163,11 @@ namespace MGTools
    * degrees of freedom on a refinement edge to those not on the refinement edge
    * of a certain level.
    */
-  template <int dim, int spacedim, typename SparsityPatternType>
+  template <int dim, int spacedim>
   void
   make_interface_sparsity_pattern(const DoFHandler<dim, spacedim> &dof_handler,
                                   const MGConstrainedDoFs &mg_constrained_dofs,
-                                  SparsityPatternType &    sparsity,
+                                  SparsityPatternBase &    sparsity,
                                   const unsigned int       level);
 
 
index abc4b66614ab19f0c623e30d091eab5627d86505..568fb029ab8f8c3e94968b160be422c3da0230a7 100644 (file)
@@ -573,13 +573,10 @@ namespace MGTools
 
 
 
-  template <int dim,
-            int spacedim,
-            typename SparsityPatternType,
-            typename number>
+  template <int dim, int spacedim, typename number>
   void
   make_sparsity_pattern(const DoFHandler<dim, spacedim> &dof,
-                        SparsityPatternType &            sparsity,
+                        SparsityPatternBase &            sparsity,
                         const unsigned int               level,
                         const AffineConstraints<number> &constraints,
                         const bool                       keep_constrained_dofs)
@@ -606,13 +603,10 @@ namespace MGTools
 
 
 
-  template <int dim,
-            int spacedim,
-            typename SparsityPatternType,
-            typename number>
+  template <int dim, int spacedim, typename number>
   void
   make_flux_sparsity_pattern(const DoFHandler<dim, spacedim> &dof,
-                             SparsityPatternType &            sparsity,
+                             SparsityPatternBase &            sparsity,
                              const unsigned int               level,
                              const AffineConstraints<number> &constraints,
                              const bool keep_constrained_dofs)
@@ -683,10 +677,10 @@ namespace MGTools
 
 
 
-  template <int dim, typename SparsityPatternType, int spacedim>
+  template <int dim, int spacedim>
   void
   make_flux_sparsity_pattern_edge(const DoFHandler<dim, spacedim> &dof,
-                                  SparsityPatternType &            sparsity,
+                                  SparsityPatternBase &            sparsity,
                                   const unsigned int               level)
   {
     Assert((level >= 1) && (level < dof.get_triangulation().n_global_levels()),
@@ -733,16 +727,11 @@ namespace MGTools
                   cell->neighbor_or_periodic_neighbor(face);
                 neighbor->get_mg_dof_indices(dofs_on_other_cell);
 
+                std::sort(dofs_on_this_cell.begin(), dofs_on_this_cell.end());
                 for (unsigned int i = 0; i < dofs_per_cell; ++i)
-                  {
-                    for (unsigned int j = 0; j < dofs_per_cell; ++j)
-                      {
-                        sparsity.add(dofs_on_other_cell[i],
-                                     dofs_on_this_cell[j]);
-                        sparsity.add(dofs_on_other_cell[j],
-                                     dofs_on_this_cell[i]);
-                      }
-                  }
+                  sparsity.add_row_entries(dofs_on_other_cell[i],
+                                           make_array_view(dofs_on_this_cell),
+                                           true);
               }
           }
       }
@@ -750,10 +739,10 @@ namespace MGTools
 
 
 
-  template <int dim, typename SparsityPatternType, int spacedim>
+  template <int dim, int spacedim>
   void
   make_flux_sparsity_pattern(const DoFHandler<dim, spacedim> &   dof,
-                             SparsityPatternType &               sparsity,
+                             SparsityPatternBase &               sparsity,
                              const unsigned int                  level,
                              const Table<2, DoFTools::Coupling> &int_mask,
                              const Table<2, DoFTools::Coupling> &flux_mask)
@@ -780,7 +769,10 @@ namespace MGTools
     const unsigned int                   total_dofs = fe.n_dofs_per_cell();
     std::vector<types::global_dof_index> dofs_on_this_cell(total_dofs);
     std::vector<types::global_dof_index> dofs_on_other_cell(total_dofs);
-    Table<2, bool>                       support_on_face(total_dofs,
+    std::vector<std::pair<types::global_dof_index, types::global_dof_index>>
+      cell_entries;
+
+    Table<2, bool> support_on_face(total_dofs,
                                    GeometryInfo<dim>::faces_per_cell);
 
     const Table<2, DoFTools::Coupling>
@@ -807,7 +799,8 @@ namespace MGTools
         for (unsigned int i = 0; i < total_dofs; ++i)
           for (unsigned int j = 0; j < total_dofs; ++j)
             if (int_dof_mask[i][j] != DoFTools::none)
-              sparsity.add(dofs_on_this_cell[i], dofs_on_this_cell[j]);
+              cell_entries.emplace_back(dofs_on_this_cell[i],
+                                        dofs_on_this_cell[j]);
 
         // Loop over all interior neighbors
         for (const unsigned int face : GeometryInfo<dim>::face_indices())
@@ -827,12 +820,12 @@ namespace MGTools
                         const bool j_non_zero_i = support_on_face(j, face);
 
                         if (flux_dof_mask(i, j) == DoFTools::always)
-                          sparsity.add(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 (flux_dof_mask(i, j) == DoFTools::nonzero &&
                             i_non_zero_i && j_non_zero_i)
-                          sparsity.add(dofs_on_this_cell[i],
-                                       dofs_on_this_cell[j]);
+                          cell_entries.emplace_back(dofs_on_this_cell[i],
+                                                    dofs_on_this_cell[j]);
                       }
                   }
               }
@@ -861,71 +854,73 @@ namespace MGTools
                           support_on_face(j, neighbor_face);
                         if (flux_dof_mask(i, j) == DoFTools::always)
                           {
-                            sparsity.add(dofs_on_this_cell[i],
-                                         dofs_on_other_cell[j]);
-                            sparsity.add(dofs_on_other_cell[i],
-                                         dofs_on_this_cell[j]);
-                            sparsity.add(dofs_on_this_cell[i],
-                                         dofs_on_this_cell[j]);
-                            sparsity.add(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) == DoFTools::nonzero)
                           {
                             if (i_non_zero_i && j_non_zero_e)
-                              sparsity.add(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)
-                              sparsity.add(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)
-                              sparsity.add(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)
-                              sparsity.add(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) == DoFTools::always)
                           {
-                            sparsity.add(dofs_on_this_cell[j],
-                                         dofs_on_other_cell[i]);
-                            sparsity.add(dofs_on_other_cell[j],
-                                         dofs_on_this_cell[i]);
-                            sparsity.add(dofs_on_this_cell[j],
-                                         dofs_on_this_cell[i]);
-                            sparsity.add(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) == DoFTools::nonzero)
                           {
                             if (j_non_zero_i && i_non_zero_e)
-                              sparsity.add(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)
-                              sparsity.add(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)
-                              sparsity.add(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)
-                              sparsity.add(dofs_on_other_cell[j],
-                                           dofs_on_other_cell[i]);
+                              cell_entries.emplace_back(dofs_on_other_cell[j],
+                                                        dofs_on_other_cell[i]);
                           }
                       }
                   }
                 face_touched[neighbor->face(neighbor_face)->index()] = true;
               }
           }
+        sparsity.add_entries(make_array_view(cell_entries));
+        cell_entries.clear();
       }
   }
 
 
 
-  template <int dim, typename SparsityPatternType, int spacedim>
+  template <int dim, int spacedim>
   void
   make_flux_sparsity_pattern_edge(const DoFHandler<dim, spacedim> &   dof,
-                                  SparsityPatternType &               sparsity,
+                                  SparsityPatternBase &               sparsity,
                                   const unsigned int                  level,
                                   const Table<2, DoFTools::Coupling> &flux_mask)
   {
@@ -955,7 +950,10 @@ namespace MGTools
     const unsigned int dofs_per_cell = dof.get_fe().n_dofs_per_cell();
     std::vector<types::global_dof_index> dofs_on_this_cell(dofs_per_cell);
     std::vector<types::global_dof_index> dofs_on_other_cell(dofs_per_cell);
-    Table<2, bool>                       support_on_face(dofs_per_cell,
+    std::vector<std::pair<types::global_dof_index, types::global_dof_index>>
+      cell_entries;
+
+    Table<2, bool> support_on_face(dofs_per_cell,
                                    GeometryInfo<dim>::faces_per_cell);
 
     const Table<2, DoFTools::Coupling> flux_dof_mask =
@@ -997,25 +995,27 @@ namespace MGTools
                       {
                         if (flux_dof_mask(i, j) != DoFTools::none)
                           {
-                            sparsity.add(dofs_on_other_cell[i],
-                                         dofs_on_this_cell[j]);
-                            sparsity.add(dofs_on_other_cell[j],
-                                         dofs_on_this_cell[i]);
+                            cell_entries.emplace_back(dofs_on_other_cell[i],
+                                                      dofs_on_this_cell[j]);
+                            cell_entries.emplace_back(dofs_on_other_cell[j],
+                                                      dofs_on_this_cell[i]);
                           }
                       }
                   }
               }
           }
+        sparsity.add_entries(make_array_view(cell_entries));
+        cell_entries.clear();
       }
   }
 
 
 
-  template <int dim, int spacedim, typename SparsityPatternType>
+  template <int dim, int spacedim>
   void
   make_interface_sparsity_pattern(const DoFHandler<dim, spacedim> &dof,
                                   const MGConstrainedDoFs &mg_constrained_dofs,
-                                  SparsityPatternType &    sparsity,
+                                  SparsityPatternBase &    sparsity,
                                   const unsigned int       level)
   {
     const types::global_dof_index n_dofs = dof.n_dofs(level);
@@ -1028,15 +1028,25 @@ namespace MGTools
 
     const unsigned int dofs_per_cell = dof.get_fe().n_dofs_per_cell();
     std::vector<types::global_dof_index> dofs_on_this_cell(dofs_per_cell);
+    std::vector<types::global_dof_index> cols;
+    cols.reserve(dofs_per_cell);
+
     for (const auto &cell : dof.cell_iterators_on_level(level))
       if (cell->is_locally_owned_on_level())
         {
           cell->get_mg_dof_indices(dofs_on_this_cell);
+          std::sort(dofs_on_this_cell.begin(), dofs_on_this_cell.end());
           for (unsigned int i = 0; i < dofs_per_cell; ++i)
-            for (unsigned int j = 0; j < dofs_per_cell; ++j)
-              if (mg_constrained_dofs.is_interface_matrix_entry(
-                    level, dofs_on_this_cell[i], dofs_on_this_cell[j]))
-                sparsity.add(dofs_on_this_cell[i], dofs_on_this_cell[j]);
+            {
+              for (unsigned int j = 0; j < dofs_per_cell; ++j)
+                if (mg_constrained_dofs.is_interface_matrix_entry(
+                      level, dofs_on_this_cell[i], dofs_on_this_cell[j]))
+                  cols.push_back(dofs_on_this_cell[j]);
+              sparsity.add_row_entries(dofs_on_this_cell[i],
+                                       make_array_view(cols),
+                                       true);
+              cols.clear();
+            }
         }
   }
 
index b679e52ecbb1cc94015b7daff87370cc4e25f1e5..3ef092a5c5a2ac7578743570feaa22d72700c681 100644 (file)
 
 
 
-for (PATTERN : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS;
-     deal_II_space_dimension : SPACE_DIMENSIONS)
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
   {
     namespace MGTools
     \{
 
 #if deal_II_dimension <= deal_II_space_dimension
       template void
-      make_sparsity_pattern<deal_II_dimension,
-                            deal_II_space_dimension,
-                            PATTERN,
-                            double>(
+      make_sparsity_pattern<deal_II_dimension, deal_II_space_dimension, double>(
         const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
-        PATTERN &,
+        SparsityPatternBase &,
         const unsigned int,
         const AffineConstraints<double> &,
         const bool);
 
       template void
-      make_sparsity_pattern<deal_II_dimension,
-                            deal_II_space_dimension,
-                            PATTERN,
-                            float>(
+      make_sparsity_pattern<deal_II_dimension, deal_II_space_dimension, float>(
         const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
-        PATTERN &,
+        SparsityPatternBase &,
         const unsigned int,
         const AffineConstraints<float> &,
         const bool);
@@ -49,8 +42,7 @@ for (PATTERN : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS;
   }
 
 
-for (PATTERN : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS;
-     deal_II_space_dimension : SPACE_DIMENSIONS)
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
   {
     namespace MGTools
     \{
@@ -58,11 +50,10 @@ for (PATTERN : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS;
 #if deal_II_dimension <= deal_II_space_dimension
       template void
       make_interface_sparsity_pattern<deal_II_dimension,
-                                      deal_II_space_dimension,
-                                      PATTERN>(
+                                      deal_II_space_dimension>(
         const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
         const MGConstrainedDoFs &,
-        PATTERN &,
+        SparsityPatternBase &,
         const unsigned int);
 #endif
 
@@ -70,9 +61,8 @@ for (PATTERN : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS;
       template void
       make_flux_sparsity_pattern<deal_II_dimension,
                                  deal_II_space_dimension,
-                                 PATTERN,
                                  double>(const DoFHandler<deal_II_dimension> &,
-                                         PATTERN &,
+                                         SparsityPatternBase &,
                                          const unsigned int,
                                          const AffineConstraints<double> &,
                                          const bool);
@@ -80,9 +70,8 @@ for (PATTERN : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS;
       template void
       make_flux_sparsity_pattern<deal_II_dimension,
                                  deal_II_space_dimension,
-                                 PATTERN,
                                  float>(const DoFHandler<deal_II_dimension> &,
-                                        PATTERN &,
+                                        SparsityPatternBase &,
                                         const unsigned int,
                                         const AffineConstraints<float> &,
                                         const bool);
@@ -90,7 +79,7 @@ for (PATTERN : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS;
       template void
       make_flux_sparsity_pattern_edge<deal_II_dimension>(
         const DoFHandler<deal_II_dimension> &,
-        PATTERN &,
+        SparsityPatternBase &,
         const unsigned int);
 
 #  if deal_II_dimension > 1
@@ -98,7 +87,7 @@ for (PATTERN : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS;
       template void
       make_flux_sparsity_pattern<deal_II_dimension>(
         const DoFHandler<deal_II_dimension> &,
-        PATTERN &,
+        SparsityPatternBase &,
         const unsigned int,
         const Table<2, DoFTools::Coupling> &,
         const Table<2, DoFTools::Coupling> &);
@@ -106,7 +95,7 @@ for (PATTERN : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS;
       template void
       make_flux_sparsity_pattern_edge<deal_II_dimension>(
         const DoFHandler<deal_II_dimension> &,
-        PATTERN &,
+        SparsityPatternBase &,
         const unsigned int,
         const Table<2, DoFTools::Coupling> &);
 #  endif

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.