]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Unify the implementation of two functions.
authorWolfgang Bangerth <bangerth@colostate.edu>
Thu, 6 Jul 2017 21:27:26 +0000 (15:27 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Fri, 7 Jul 2017 20:05:23 +0000 (14:05 -0600)
include/deal.II/dofs/dof_renumbering.h
source/dofs/dof_renumbering.cc
source/dofs/dof_renumbering.inst.in

index 42635eba43745737fb1ef42a0caaf8af9a92f4a2..094cc4e0bedcca43df91ecce35a792747544a6d5 100644 (file)
@@ -614,22 +614,13 @@ namespace DoFRenumbering
    * For finite elements with only one component, or a single non-primitive
    * base element, this function is the identity operation.
    */
-  template <int dim, int spacedim>
+  template <typename DoFHandlerType>
   void
-  component_wise (DoFHandler<dim,spacedim>        &dof_handler,
+  component_wise (DoFHandlerType                  &dof_handler,
                   const std::vector<unsigned int> &target_component
                   = std::vector<unsigned int>());
 
 
-  /**
-   * Sort the degrees of freedom by component. It does the same thing as the
-   * above function.
-   */
-  template <int dim>
-  void
-  component_wise (hp::DoFHandler<dim>             &dof_handler,
-                  const std::vector<unsigned int> &target_component = std::vector<unsigned int> ());
-
   /**
    * Sort the degrees of freedom by component. It does the same thing as the
    * above function, only that it does this for one single level of a
index bbaf17d46a4dc82bab29221e30a4c9692dacedeb..fe2ac84df706327f566fc3cee07f7eaeea155618 100644 (file)
@@ -496,23 +496,23 @@ namespace DoFRenumbering
 
 
 
-  template <int dim, int spacedim>
+  template <typename DoFHandlerType>
   void
-  component_wise (DoFHandler<dim,spacedim>        &dof_handler,
+  component_wise (DoFHandlerType                  &dof_handler,
                   const std::vector<unsigned int> &component_order_arg)
   {
     std::vector<types::global_dof_index> renumbering (dof_handler.n_locally_owned_dofs(),
                                                       numbers::invalid_dof_index);
 
-    typename DoFHandler<dim,spacedim>::active_cell_iterator
+    typename DoFHandlerType::active_cell_iterator
     start = dof_handler.begin_active();
-    const typename DoFHandler<dim,spacedim>::level_cell_iterator
+    const typename DoFHandlerType::level_cell_iterator
     end = dof_handler.end();
 
     const types::global_dof_index result =
-      compute_component_wise<dim,spacedim,
-      typename DoFHandler<dim,spacedim>::active_cell_iterator,
-      typename DoFHandler<dim,spacedim>::level_cell_iterator>
+      compute_component_wise<DoFHandlerType::dimension,DoFHandlerType::space_dimension,
+      typename DoFHandlerType::active_cell_iterator,
+      typename DoFHandlerType::level_cell_iterator>
       (renumbering, start, end, component_order_arg, false);
     if (result == 0)
       return;
@@ -532,40 +532,6 @@ namespace DoFRenumbering
             ExcRenumberingIncomplete());
 
     dof_handler.renumber_dofs (renumbering);
-
-    // for (unsigned int level=0;level<dof_handler.get_triangulation().n_levels();++level)
-    //   if (dof_handler.n_dofs(level) != numbers::invalid_dof_index)
-    //  component_wise(dof_handler, level, component_order_arg);
-  }
-
-
-
-  template <int dim>
-  void
-  component_wise (hp::DoFHandler<dim>             &dof_handler,
-                  const std::vector<unsigned int> &component_order_arg)
-  {
-//TODO: Merge with previous function
-    std::vector<types::global_dof_index> renumbering (dof_handler.n_dofs(),
-                                                      numbers::invalid_dof_index);
-
-    typename hp::DoFHandler<dim>::active_cell_iterator
-    start = dof_handler.begin_active();
-    const typename hp::DoFHandler<dim>::level_cell_iterator
-    end = dof_handler.end();
-
-    const types::global_dof_index result =
-      compute_component_wise<hp::DoFHandler<dim>::dimension,hp::DoFHandler<dim>::space_dimension,
-      typename hp::DoFHandler<dim>::active_cell_iterator,
-      typename hp::DoFHandler<dim>::level_cell_iterator>
-      (renumbering, start, end, component_order_arg, false);
-
-    if (result == 0) return;
-
-    Assert (result == dof_handler.n_dofs(),
-            ExcRenumberingIncomplete());
-
-    dof_handler.renumber_dofs (renumbering);
   }
 
 
index f3546de3eb3c8eccbd863cd4b7cd876c8d76e6f0..8ba908e2e5e5f69292acf93d818f8b67dd7f8a89 100644 (file)
@@ -43,10 +43,15 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
      const std::vector<types::global_dof_index>&);
 
     template
-    void component_wise<deal_II_dimension,deal_II_space_dimension>
+    void component_wise<DoFHandler<deal_II_dimension,deal_II_space_dimension> >
     (DoFHandler<deal_II_dimension,deal_II_space_dimension>&,
      const std::vector<unsigned int>&);
 
+    template
+    void component_wise<hp::DoFHandler<deal_II_dimension,deal_II_space_dimension> >
+    (hp::DoFHandler<deal_II_dimension,deal_II_space_dimension>&,
+     const std::vector<unsigned int>&);
+
     template
     void block_wise<deal_II_dimension,deal_II_space_dimension>
     (DoFHandler<deal_II_dimension,deal_II_space_dimension>&);
@@ -149,11 +154,6 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
      const bool,
      const std::vector<types::global_dof_index>&);
 
-    template
-    void component_wise
-    (hp::DoFHandler<deal_II_dimension>&,
-     const std::vector<unsigned int>&);
-
     template
     void component_wise
     (DoFHandler<deal_II_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.