]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove DoFHandlerType from mg_tools.h
authorPeter Munch <peterrmuench@gmail.com>
Mon, 29 Jun 2020 06:57:55 +0000 (08:57 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Mon, 29 Jun 2020 06:57:55 +0000 (08:57 +0200)
include/deal.II/multigrid/mg_tools.h
source/multigrid/mg_tools.cc
source/multigrid/mg_tools.inst.in

index 142ea23fbe215f4d86f7c08238a264f03f5ba2d1..d6d821c26a1f6e238c0d3715452d6305ea31b1ac 100644 (file)
@@ -80,11 +80,11 @@ namespace MGTools
    * There is no need to consider hanging nodes here, since only one level is
    * considered.
    */
-  template <typename DoFHandlerType, typename SparsityPatternType>
+  template <int dim, int spacedim, typename SparsityPatternType>
   void
-  make_sparsity_pattern(const DoFHandlerType &dof_handler,
-                        SparsityPatternType & sparsity,
-                        const unsigned int    level);
+  make_sparsity_pattern(const DoFHandler<dim, spacedim> &dof_handler,
+                        SparsityPatternType &            sparsity,
+                        const unsigned int               level);
 
   /**
    * Make a sparsity pattern including fluxes of discontinuous Galerkin
@@ -151,9 +151,9 @@ namespace MGTools
    * degrees of freedom on a refinement edge to those not on the refinement edge
    * of a certain level.
    */
-  template <typename DoFHandlerType, typename SparsityPatternType>
+  template <int dim, int spacedim, typename SparsityPatternType>
   void
-  make_interface_sparsity_pattern(const DoFHandlerType &   dof_handler,
+  make_interface_sparsity_pattern(const DoFHandler<dim, spacedim> &dof_handler,
                                   const MGConstrainedDoFs &mg_constrained_dofs,
                                   SparsityPatternType &    sparsity,
                                   const unsigned int       level);
@@ -165,10 +165,10 @@ namespace MGTools
    * Result is a vector containing for each level a vector containing the
    * number of dofs for each block (access is <tt>result[level][block]</tt>).
    */
-  template <typename DoFHandlerType>
+  template <int dim, int spacedim>
   void
   count_dofs_per_block(
-    const DoFHandlerType &                             dof_handler,
+    const DoFHandler<dim, spacedim> &                  dof_handler,
     std::vector<std::vector<types::global_dof_index>> &dofs_per_block,
     std::vector<unsigned int>                          target_block = {});
 
index e33c86e778487cb332e2666653548494aeef7fb2..bcfd0c89b53542b6fa487ba8b80c6e5841600fc1 100644 (file)
@@ -560,11 +560,11 @@ namespace MGTools
 
 
 
-  template <typename DoFHandlerType, typename SparsityPatternType>
+  template <int dim, int spacedim, typename SparsityPatternType>
   void
-  make_sparsity_pattern(const DoFHandlerType &dof,
-                        SparsityPatternType & sparsity,
-                        const unsigned int    level)
+  make_sparsity_pattern(const DoFHandler<dim, spacedim> &dof,
+                        SparsityPatternType &            sparsity,
+                        const unsigned int               level)
   {
     const types::global_dof_index n_dofs = dof.n_dofs(level);
     (void)n_dofs;
@@ -575,9 +575,9 @@ namespace MGTools
            ExcDimensionMismatch(sparsity.n_cols(), n_dofs));
 
     const unsigned int dofs_per_cell = dof.get_fe().dofs_per_cell;
-    std::vector<types::global_dof_index>   dofs_on_this_cell(dofs_per_cell);
-    typename DoFHandlerType::cell_iterator cell = dof.begin(level),
-                                           endc = dof.end(level);
+    std::vector<types::global_dof_index> dofs_on_this_cell(dofs_per_cell);
+    typename DoFHandler<dim, spacedim>::cell_iterator cell = dof.begin(level),
+                                                      endc = dof.end(level);
     for (; cell != endc; ++cell)
       if (dof.get_triangulation().locally_owned_subdomain() ==
             numbers::invalid_subdomain_id ||
@@ -1018,9 +1018,9 @@ namespace MGTools
 
 
 
-  template <typename DoFHandlerType, typename SparsityPatternType>
+  template <int dim, int spacedim, typename SparsityPatternType>
   void
-  make_interface_sparsity_pattern(const DoFHandlerType &   dof,
+  make_interface_sparsity_pattern(const DoFHandler<dim, spacedim> &dof,
                                   const MGConstrainedDoFs &mg_constrained_dofs,
                                   SparsityPatternType &    sparsity,
                                   const unsigned int       level)
@@ -1034,9 +1034,9 @@ namespace MGTools
            ExcDimensionMismatch(sparsity.n_cols(), n_dofs));
 
     const unsigned int dofs_per_cell = dof.get_fe().dofs_per_cell;
-    std::vector<types::global_dof_index>   dofs_on_this_cell(dofs_per_cell);
-    typename DoFHandlerType::cell_iterator cell = dof.begin(level),
-                                           endc = dof.end(level);
+    std::vector<types::global_dof_index> dofs_on_this_cell(dofs_per_cell);
+    typename DoFHandler<dim, spacedim>::cell_iterator cell = dof.begin(level),
+                                                      endc = dof.end(level);
     for (; cell != endc; ++cell)
       if (cell->is_locally_owned_on_level())
         {
@@ -1153,18 +1153,16 @@ namespace MGTools
 
 
 
-  template <typename DoFHandlerType>
+  template <int dim, int spacedim>
   void
   count_dofs_per_block(
-    const DoFHandlerType &                             dof_handler,
+    const DoFHandler<dim, spacedim> &                  dof_handler,
     std::vector<std::vector<types::global_dof_index>> &dofs_per_block,
     std::vector<unsigned int>                          target_block)
   {
-    const FiniteElement<DoFHandlerType::dimension,
-                        DoFHandlerType::space_dimension> &fe =
-      dof_handler.get_fe();
-    const unsigned int n_blocks = fe.n_blocks();
-    const unsigned int n_levels =
+    const FiniteElement<dim, spacedim> &fe       = dof_handler.get_fe();
+    const unsigned int                  n_blocks = fe.n_blocks();
+    const unsigned int                  n_levels =
       dof_handler.get_triangulation().n_global_levels();
 
     AssertDimension(dofs_per_block.size(), n_levels);
@@ -1213,11 +1211,10 @@ namespace MGTools
         for (unsigned int i = 0; i < n_blocks; ++i)
           {
             void (*fun_ptr)(const unsigned int level,
-                            const DoFHandlerType &,
+                            const DoFHandler<dim, spacedim> &,
                             const BlockMask &,
                             std::vector<bool> &) =
-              &DoFTools::extract_level_dofs<DoFHandlerType::dimension,
-                                            DoFHandlerType::space_dimension>;
+              &DoFTools::extract_level_dofs<dim, spacedim>;
 
             std::vector<bool> tmp(n_blocks, false);
             tmp[i]          = true;
index 9d70c819f18fc06fcb05c356f33eb19a840c6b91..5aaff41760a4784dba1570a4fd540f14aae9c993 100644 (file)
@@ -15,8 +15,7 @@
 
 
 
-for (DoFHandler : DOFHANDLER_TEMPLATES; PATTERN : SPARSITY_PATTERNS;
-     deal_II_dimension : DIMENSIONS;
+for (PATTERN : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS;
      deal_II_space_dimension : SPACE_DIMENSIONS)
   {
     namespace MGTools
@@ -24,11 +23,12 @@ for (DoFHandler : DOFHANDLER_TEMPLATES; PATTERN : SPARSITY_PATTERNS;
 
 #if deal_II_dimension <= deal_II_space_dimension
       template void
-      make_sparsity_pattern<
-        DoFHandler<deal_II_dimension, deal_II_space_dimension>,
-        PATTERN>(const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
-                 PATTERN &,
-                 const unsigned int);
+      make_sparsity_pattern<deal_II_dimension,
+                            deal_II_space_dimension,
+                            PATTERN>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        PATTERN &,
+        const unsigned int);
 #endif
     \}
   }
@@ -42,12 +42,13 @@ for (PATTERN : SPARSITY_PATTERNS; deal_II_dimension : DIMENSIONS;
 
 #if deal_II_dimension <= deal_II_space_dimension
       template void
-      make_interface_sparsity_pattern<
-        DoFHandler<deal_II_dimension, deal_II_space_dimension>,
-        PATTERN>(const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
-                 const MGConstrainedDoFs &,
-                 PATTERN &,
-                 const unsigned int);
+      make_interface_sparsity_pattern<deal_II_dimension,
+                                      deal_II_space_dimension,
+                                      PATTERN>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const MGConstrainedDoFs &,
+        PATTERN &,
+        const unsigned int);
 #endif
 
 #if deal_II_dimension == deal_II_space_dimension

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.