]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Patch from Markus Buerg: Implement DoFRenumbering::component_wise for hp::DoFHandler.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 2 Oct 2009 00:26:52 +0000 (00:26 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 2 Oct 2009 00:26:52 +0000 (00:26 +0000)
git-svn-id: https://svn.dealii.org/trunk@19657 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_renumbering.h
deal.II/deal.II/source/dofs/dof_renumbering.cc

index 4ba0ebe0453a73487276cf4b0b624163e16ee776..9c91a358b4812a0140abce2d676ce6959548b176 100644 (file)
@@ -18,6 +18,7 @@
 #include <base/exceptions.h>
 #include <base/point.h>
 #include <dofs/dof_handler.h>
+#include <hp/dof_handler.h>
 #include <multigrid/mg_dof_handler.h>
 
 #include <vector>
@@ -693,6 +694,18 @@ namespace DoFRenumbering
                  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
@@ -724,7 +737,6 @@ namespace DoFRenumbering
   void
   component_wise (MGDoFHandler<dim>&               dof_handler,
                  const std::vector<unsigned int>& target_component = std::vector<unsigned int>());
-
     
                                   /**
                                    * Computes the renumbering
index a79a813d002ce90b8ec55695fda4d1068d0a46f8..ab4b2683424fe46a576ffc992fa65ca7ea41d3f7 100644 (file)
@@ -23,6 +23,7 @@
 #include <grid/tria.h>
 #include <dofs/dof_handler.h>
 #include <hp/dof_handler.h>
+#include <hp/fe_collection.h>
 #include <dofs/dof_tools.h>
 #include <fe/fe.h>
 #include <dofs/dof_renumbering.h>
@@ -112,7 +113,7 @@ namespace DoFRenumbering
        T::operator==;
     };
   }
-  
+
   namespace boost
   {
     namespace types
@@ -120,7 +121,7 @@ namespace DoFRenumbering
       using namespace ::boost;
       using namespace std;
 
-      typedef adjacency_list<vecS, vecS, undirectedS, 
+      typedef adjacency_list<vecS, vecS, undirectedS,
                             property<vertex_color_t, default_color_type,
                                      property<vertex_degree_t,int> > > Graph;
       typedef graph_traits<Graph>::vertex_descriptor Vertex;
@@ -128,7 +129,7 @@ namespace DoFRenumbering
 
       typedef std::pair<size_type, size_type> Pair;
     }
-  
+
 
     namespace internal
     {
@@ -143,7 +144,7 @@ namespace DoFRenumbering
                                   // (faster than directly submitting
                                   // indices)
          ConstraintMatrix constraints;
-         if (use_constraints) 
+         if (use_constraints)
            DoFTools::make_hanging_node_constraints (dof_handler, constraints);
          constraints.close ();
          CompressedSimpleSparsityPattern csp (dof_handler.n_dofs(),
@@ -165,7 +166,7 @@ namespace DoFRenumbering
       }
     }
 
-  
+
     template <class DH>
     void
     Cuthill_McKee (DH&              dof_handler,
@@ -197,7 +198,7 @@ namespace DoFRenumbering
        graph_degree;
 
       internal::create_graph (dof_handler, use_constraints, graph, graph_degree);
-    
+
       types::property_map<types::Graph, types::vertex_index_t>::type
        index_map = get(::boost::vertex_index, graph);
 
@@ -212,7 +213,7 @@ namespace DoFRenumbering
        ::boost::cuthill_mckee_ordering(graph, inv_perm.begin(),
                                        get(::boost::vertex_color, graph),
                                        make_degree_map(graph));
-    
+
       for (types::size_type c = 0; c != inv_perm.size(); ++c)
        new_dof_indices[index_map[inv_perm[c]]] = c;
 
@@ -254,7 +255,7 @@ namespace DoFRenumbering
        graph_degree;
 
       internal::create_graph (dof_handler, use_constraints, graph, graph_degree);
-    
+
       types::property_map<types::Graph, types::vertex_index_t>::type
        index_map = get(::boost::vertex_index, graph);
 
@@ -265,7 +266,7 @@ namespace DoFRenumbering
        ::boost::king_ordering(graph, inv_perm.rbegin());
       else
        ::boost::king_ordering(graph, inv_perm.begin());
-    
+
       for (types::size_type c = 0; c != inv_perm.size(); ++c)
        new_dof_indices[index_map[inv_perm[c]]] = c;
 
@@ -344,8 +345,8 @@ namespace DoFRenumbering
                  add_edge (dofs_on_this_cell[j], dofs_on_this_cell[i], G);
                }
        }
-  
-  
+
+
       typedef std::vector<int> Vector;
 
 
@@ -393,7 +394,7 @@ namespace DoFRenumbering
     }
 
   }  // namespace boost
-  
+
 
 
   template <class DH>
@@ -439,10 +440,10 @@ namespace DoFRenumbering
                                     // empty, and calling condense with
                                     // it is a no-op
     ConstraintMatrix constraints;
-    if (use_constraints) 
+    if (use_constraints)
       DoFTools::make_hanging_node_constraints (dof_handler, constraints);
     constraints.close ();
-    
+
     SparsityPattern sparsity;
     if (DH::dimension < 2)
       {
@@ -462,7 +463,7 @@ namespace DoFRenumbering
 
                                     // constraints are not needed anymore
     constraints.clear ();
-  
+
     Assert(new_indices.size() == sparsity.n_rows(),
           ExcDimensionMismatch(new_indices.size(),
                                sparsity.n_rows()));
@@ -482,12 +483,12 @@ namespace DoFRenumbering
                      const std::vector<unsigned int> &starting_indices)
   {
 //TODO: we should be doing the same here as in the other compute_CMK function to preserve some memory
-  
+
                                     // make the connection graph
     SparsityPattern sparsity (dof_handler.n_dofs(level),
                              dof_handler.max_couplings_between_dofs());
     MGTools::make_sparsity_pattern (dof_handler, sparsity, level);
-    
+
     std::vector<unsigned int> new_indices(sparsity.n_rows());
     SparsityTools::reorder_Cuthill_McKee (sparsity, new_indices,
                                          starting_indices);
@@ -504,7 +505,6 @@ namespace DoFRenumbering
 
   template <int dim>
   void
-
   component_wise (DoFHandler<dim>                 &dof_handler,
                  const std::vector<unsigned int> &component_order_arg)
   {
@@ -520,24 +520,56 @@ namespace DoFRenumbering
     const typename DoFHandler<dim>::cell_iterator
       end = dof_handler.end();
 
-    unsigned int result =
+    const unsigned int result =
       compute_component_wise<dim, ITERATOR,
       typename DoFHandler<dim>::cell_iterator>(renumbering,
                                               start, end,
                                               component_order_arg);
 
     if (result == 0) return;
-  
+
     Assert (result == dof_handler.n_dofs(),
            ExcRenumberingIncomplete());
-  
+
     dof_handler.renumber_dofs (renumbering);
   }
 
 
 
   template <int dim>
-  void 
+  void
+  component_wise (hp::DoFHandler<dim>                 &dof_handler,
+                 const std::vector<unsigned int> &component_order_arg)
+  {
+    std::vector<unsigned int> renumbering (dof_handler.n_dofs(),
+                                          hp::DoFHandler<dim>::invalid_dof_index);
+
+    typedef
+      internal::WrapDoFIterator<typename hp::DoFHandler<dim>::active_cell_iterator> ITERATOR;
+
+    typename hp::DoFHandler<dim>::active_cell_iterator
+      istart = dof_handler.begin_active();
+    ITERATOR start = istart;
+    const typename hp::DoFHandler<dim>::cell_iterator
+      end = dof_handler.end();
+
+    const unsigned int result =
+      compute_component_wise<dim, ITERATOR,
+      typename hp::DoFHandler<dim>::cell_iterator>(renumbering,
+                                                  start, end,
+                                                  component_order_arg);
+
+    if (result == 0) return;
+
+    Assert (result == dof_handler.n_dofs(),
+           ExcRenumberingIncomplete());
+
+    dof_handler.renumber_dofs (renumbering);
+  }
+
+
+  template <int dim>
+  void
   component_wise (MGDoFHandler<dim> &dof_handler,
                  const unsigned int level,
                  const std::vector<unsigned int> &component_order_arg)
@@ -547,7 +579,7 @@ namespace DoFRenumbering
 
     typedef
       internal::WrapMGDoFIterator<typename MGDoFHandler<dim>::cell_iterator> ITERATOR;
-  
+
     typename MGDoFHandler<dim>::cell_iterator
       istart =dof_handler.begin(level);
     ITERATOR start = istart;
@@ -560,7 +592,7 @@ namespace DoFRenumbering
        renumbering, start, end, component_order_arg);
 
     if (result == 0) return;
-  
+
     Assert (result == dof_handler.n_dofs(level),
            ExcRenumberingIncomplete());
 
@@ -569,9 +601,8 @@ namespace DoFRenumbering
   }
 
 
-
   template <int dim>
-  void 
+  void
   component_wise (MGDoFHandler<dim> &dof_handler,
                  const std::vector<unsigned int> &component_order_arg)
   {
@@ -587,7 +618,7 @@ namespace DoFRenumbering
       = &component_wise<dim>;
     Threads::Task<>
       task = Threads::new_task (non_mg_part, dof_handler, component_order_arg);
-    
+
     for (unsigned int level=0; level<dof_handler.get_tria().n_levels(); ++level)
       component_wise (dof_handler, level, component_order_arg);
 
@@ -604,22 +635,23 @@ namespace DoFRenumbering
                          const ENDITERATOR& end,
                          const std::vector<unsigned int> &component_order_arg)
   {
-//TODO: Modify for hp
-    const FiniteElement<dim>& fe = start->get_fe();
-    const unsigned int dofs_per_cell = fe.dofs_per_cell;
-
-//    Assert (new_indices.size() ==  start->dof_handler().n_dofs(),
-//       ExcDimensionMismatch(new_indices.size(),
-//                            start->dof_handler.n_dofs()));
+    const hp::FECollection<dim> fe_collection (start->get_fe ());
 
                                     // do nothing if the FE has only
                                     // one component
-    if (fe.n_components() == 1)
+    if (fe_collection.n_components() == 1)
       {
        new_indices.resize(0);
        return 0;
       }
 
+    const FiniteElement<dim>& fe = fe_collection[0];
+    const unsigned int dofs_per_cell = fe.dofs_per_cell;
+
+//    Assert (new_indices.size() ==  start->dof_handler().n_dofs(),
+//       ExcDimensionMismatch(new_indices.size(),
+//                            start->dof_handler.n_dofs()));
+
                                     // Copy last argument into a
                                     // writable vector.
     std::vector<unsigned int> component_order (component_order_arg);
@@ -633,11 +665,11 @@ namespace DoFRenumbering
 
     Assert (component_order.size() == fe.n_components(),
            ExcDimensionMismatch(component_order.size(), fe.n_components()));
-  
+
     for (unsigned int i=0; i< component_order.size(); ++i)
       Assert(component_order[i] < fe.n_components(),
             ExcIndexRange(component_order[i], 0, fe.n_components()));
-  
+
                                     // vector to hold the dof indices on
                                     // the cell we visit at a time
     std::vector<unsigned int> local_dof_indices(dofs_per_cell);
@@ -677,7 +709,7 @@ namespace DoFRenumbering
                                           // component
          component_list[i] = component_order[c];
        };
-  
+
                                     // set up a map where for each
                                     // component the respective degrees
                                     // of freedom are collected.
@@ -700,7 +732,7 @@ namespace DoFRenumbering
        for (unsigned int i=0; i<dofs_per_cell; ++i)
          component_to_dof_map[component_list[i]].push_back (local_dof_indices[i]);
       };
-  
+
                                     // now we've got all indices sorted
                                     // into buckets labelled with their
                                     // target component number. we've
@@ -741,11 +773,11 @@ namespace DoFRenumbering
     for (unsigned int c=0; c<fe.n_components(); ++c)
       {
        const unsigned int component = c;
-      
+
        const typename std::vector<unsigned int>::const_iterator
          begin_of_component = component_to_dof_map[component].begin(),
          end_of_component   = component_to_dof_map[component].end();
-            
+
        for (typename std::vector<unsigned int>::const_iterator
               dof_index = begin_of_component;
             dof_index != end_of_component; ++dof_index)
@@ -785,11 +817,11 @@ namespace DoFRenumbering
                                     // their selection state
     Assert (new_indices.size() == n_dofs,
            ExcDimensionMismatch(new_indices.size(), n_dofs));
-  
+
     const unsigned int   n_selected_dofs = std::count (selected_dofs.begin(),
                                                       selected_dofs.end(),
                                                       false);
-  
+
     unsigned int next_unselected = 0;
     unsigned int next_selected   = n_selected_dofs;
     for (unsigned int i=0; i<n_dofs; ++i)
@@ -817,11 +849,11 @@ namespace DoFRenumbering
     std::vector<unsigned int> renumbering(dof.n_dofs());
     std::vector<unsigned int> reverse(dof.n_dofs());
     compute_cell_wise_dg(renumbering, reverse, dof, cells);
-  
+
     dof.renumber_dofs(renumbering);
   }
 
-  
+
 
   template <class DH>
   void
@@ -832,7 +864,7 @@ namespace DoFRenumbering
     std::vector<unsigned int> renumbering(dof.n_dofs());
     std::vector<unsigned int> reverse(dof.n_dofs());
     compute_cell_wise(renumbering, reverse, dof, cells);
-  
+
     dof.renumber_dofs(renumbering);
   }
 
@@ -858,16 +890,16 @@ namespace DoFRenumbering
                                     // Later, irs inverse is computed
                                     // into new_indices, which is the
                                     // return argument.
-  
+
     Assert(new_indices.size() == n_global_dofs,
           ExcDimensionMismatch(new_indices.size(), n_global_dofs));
     Assert(reverse.size() == n_global_dofs,
           ExcDimensionMismatch(reverse.size(), n_global_dofs));
-  
+
     std::vector<unsigned int> cell_dofs;
 
     unsigned int global_index = 0;
-  
+
     typename std::vector<typename DH::cell_iterator>::const_iterator cell;
 
     for(cell = cells.begin(); cell != cells.end(); ++cell)
@@ -879,7 +911,7 @@ namespace DoFRenumbering
                                         // numbers.
        unsigned int n_cell_dofs = (*cell)->get_fe().n_dofs_per_cell();
        cell_dofs.resize(n_cell_dofs);
-      
+
        (*cell)->get_dof_indices(cell_dofs);
 
                                         // Sort here to make sure that
@@ -913,14 +945,14 @@ namespace DoFRenumbering
                                dof.get_tria().n_active_cells()));
 
     unsigned int n_global_dofs = dof.n_dofs();
-  
+
                                     // Actually, we compute the
                                     // inverse of the reordering
                                     // vector, called reverse here.
                                     // Later, irs inverse is computed
                                     // into new_indices, which is the
                                     // return argument.
-  
+
     Assert(new_indices.size() == n_global_dofs,
           ExcDimensionMismatch(new_indices.size(), n_global_dofs));
     Assert(reverse.size() == n_global_dofs,
@@ -933,7 +965,7 @@ namespace DoFRenumbering
     std::vector<unsigned int> cell_dofs;
 
     unsigned int global_index = 0;
-  
+
     typename std::vector<typename DH::cell_iterator>::const_iterator cell;
 
     for(cell = cells.begin(); cell != cells.end(); ++cell)
@@ -944,7 +976,7 @@ namespace DoFRenumbering
                                         // numbers.
        unsigned int n_cell_dofs = (*cell)->get_fe().n_dofs_per_cell();
        cell_dofs.resize(n_cell_dofs);
-      
+
        (*cell)->get_dof_indices(cell_dofs);
 
                                         // Sort here to make sure that
@@ -968,7 +1000,7 @@ namespace DoFRenumbering
       new_indices[reverse[i]] = i;
   }
 
-  
+
 
   template <int dim>
   void cell_wise_dg (
@@ -978,12 +1010,12 @@ namespace DoFRenumbering
   {
     std::vector<unsigned int> renumbering(dof.n_dofs(level));
     std::vector<unsigned int> reverse(dof.n_dofs(level));
-  
+
     compute_cell_wise_dg(renumbering, reverse, dof, level, cells);
     dof.renumber_dofs(level, renumbering);
   }
 
-  
+
 
   template <int dim>
   void cell_wise (
@@ -993,12 +1025,12 @@ namespace DoFRenumbering
   {
     std::vector<unsigned int> renumbering(dof.n_dofs(level));
     std::vector<unsigned int> reverse(dof.n_dofs(level));
-  
+
     compute_cell_wise(renumbering, reverse, dof, level, cells);
     dof.renumber_dofs(level, renumbering);
   }
 
-  
+
 
   template <int dim>
   void compute_cell_wise_dg (
@@ -1023,19 +1055,19 @@ namespace DoFRenumbering
              Assert(dof.get_fe().n_dofs_per_vertex()==0,
                     ExcNotDGFEM());
       }
-  
+
     Assert (new_order.size() == dof.n_dofs(level),
            ExcDimensionMismatch(new_order.size(), dof.n_dofs(level)));
     Assert (reverse.size() == dof.n_dofs(level),
            ExcDimensionMismatch(reverse.size(), dof.n_dofs(level)));
-  
+
     unsigned int n_global_dofs = dof.n_dofs(level);
     unsigned int n_cell_dofs = dof.get_fe().n_dofs_per_cell();
-  
+
     std::vector<unsigned int> cell_dofs(n_cell_dofs);
-  
+
     unsigned int global_index = 0;
-  
+
     typename std::vector<typename MGDoFHandler<dim>::cell_iterator>::const_iterator cell;
 
     for(cell = cells.begin(); cell != cells.end(); ++cell)
@@ -1053,7 +1085,7 @@ namespace DoFRenumbering
     Assert(global_index == n_global_dofs, ExcRenumberingIncomplete());
 
     for (unsigned int i=0;i<new_order.size(); ++i)
-      new_order[reverse[i]] = i;  
+      new_order[reverse[i]] = i;
   }
 
 
@@ -1073,15 +1105,15 @@ namespace DoFRenumbering
            ExcDimensionMismatch(new_order.size(), dof.n_dofs(level)));
     Assert (reverse.size() == dof.n_dofs(level),
            ExcDimensionMismatch(reverse.size(), dof.n_dofs(level)));
-  
+
     unsigned int n_global_dofs = dof.n_dofs(level);
     unsigned int n_cell_dofs = dof.get_fe().n_dofs_per_cell();
-  
+
     std::vector<bool> already_sorted(n_global_dofs, false);
     std::vector<unsigned int> cell_dofs(n_cell_dofs);
-  
+
     unsigned int global_index = 0;
-  
+
     typename std::vector<typename MGDoFHandler<dim>::cell_iterator>::const_iterator cell;
 
     for(cell = cells.begin(); cell != cells.end(); ++cell)
@@ -1090,7 +1122,7 @@ namespace DoFRenumbering
 
        (*cell)->get_mg_dof_indices(cell_dofs);
        std::sort(cell_dofs.begin(), cell_dofs.end());
-      
+
        for (unsigned int i=0;i<n_cell_dofs;++i)
          {
            if (!already_sorted[cell_dofs[i]])
@@ -1103,7 +1135,7 @@ namespace DoFRenumbering
     Assert(global_index == n_global_dofs, ExcRenumberingIncomplete());
 
     for (unsigned int i=0;i<new_order.size(); ++i)
-      new_order[reverse[i]] = i;  
+      new_order[reverse[i]] = i;
   }
 
 
@@ -1143,10 +1175,10 @@ namespace DoFRenumbering
     std::vector<typename DH::cell_iterator>
       ordered_cells(dof.get_tria().n_active_cells());
     const CompareDownstream<typename DH::cell_iterator, dim> comparator(direction);
-  
+
     typename DH::active_cell_iterator begin = dof.begin_active();
     typename DH::active_cell_iterator end = dof.end();
-  
+
     copy (begin, end, ordered_cells.begin());
     std::sort (ordered_cells.begin(), ordered_cells.end(), comparator);
 
@@ -1154,7 +1186,7 @@ namespace DoFRenumbering
     compute_cell_wise_dg(new_indices, reverse, dof, ordered_cells);
   }
 
-  
+
 
   template <class DH, int dim>
   void
@@ -1167,13 +1199,13 @@ namespace DoFRenumbering
     std::vector<typename DH::cell_iterator>
       ordered_cells(dof.get_tria().n_active_cells());
     const CompareDownstream<typename DH::cell_iterator, dim> comparator(direction);
-  
+
     typename DH::active_cell_iterator begin = dof.begin_active();
     typename DH::active_cell_iterator end = dof.end();
-  
+
     copy (begin, end, ordered_cells.begin());
     std::sort (ordered_cells.begin(), ordered_cells.end(), comparator);
-  
+
     compute_cell_wise_dg(new_indices, reverse, dof, ordered_cells);
   }
 
@@ -1189,18 +1221,18 @@ namespace DoFRenumbering
     std::vector<typename DH::cell_iterator>
       ordered_cells(dof.get_tria().n_active_cells());
     const CompareDownstream<typename DH::cell_iterator, dim> comparator(direction);
-  
+
     typename DH::active_cell_iterator begin = dof.begin_active();
     typename DH::active_cell_iterator end = dof.end();
-  
+
     copy (begin, end, ordered_cells.begin());
     std::sort (ordered_cells.begin(), ordered_cells.end(), comparator);
-  
+
     compute_cell_wise(new_indices, reverse, dof, ordered_cells);
   }
 
 
-  
+
   template <int dim>
   void downstream_dg (MGDoFHandler<dim>& dof,
                      const unsigned int level,
@@ -1209,11 +1241,11 @@ namespace DoFRenumbering
     std::vector<unsigned int> renumbering(dof.n_dofs(level));
     std::vector<unsigned int> reverse(dof.n_dofs(level));
     compute_downstream_dg(renumbering, reverse, dof, level, direction);
-  
+
     dof.renumber_dofs(level, renumbering);
   }
 
-  
+
 
   template <int dim>
   void downstream (MGDoFHandler<dim>& dof,
@@ -1223,11 +1255,11 @@ namespace DoFRenumbering
     std::vector<unsigned int> renumbering(dof.n_dofs(level));
     std::vector<unsigned int> reverse(dof.n_dofs(level));
     compute_downstream(renumbering, reverse, dof, level, direction);
-  
+
     dof.renumber_dofs(level, renumbering);
   }
 
-  
+
 
   template <int dim>
   void
@@ -1241,17 +1273,17 @@ namespace DoFRenumbering
     std::vector<typename MGDoFHandler<dim>::cell_iterator>
       ordered_cells(dof.get_tria().n_cells(level));
     const CompareDownstream<typename MGDoFHandler<dim>::cell_iterator, dim> comparator(direction);
-  
+
     typename MGDoFHandler<dim>::cell_iterator begin = dof.begin(level);
     typename MGDoFHandler<dim>::cell_iterator end = dof.end(level);
-  
+
     std::copy (begin, end, ordered_cells.begin());
     std::sort (ordered_cells.begin(), ordered_cells.end(), comparator);
-  
+
     compute_cell_wise_dg(new_indices, reverse, dof, level, ordered_cells);
   }
 
-  
+
 
   template <int dim>
   void
@@ -1265,13 +1297,13 @@ namespace DoFRenumbering
     std::vector<typename MGDoFHandler<dim>::cell_iterator>
       ordered_cells(dof.get_tria().n_cells(level));
     const CompareDownstream<typename MGDoFHandler<dim>::cell_iterator, dim> comparator(direction);
-  
+
     typename MGDoFHandler<dim>::cell_iterator begin = dof.begin(level);
     typename MGDoFHandler<dim>::cell_iterator end = dof.end(level);
-  
+
     std::copy (begin, end, ordered_cells.begin());
     std::sort (ordered_cells.begin(), ordered_cells.end(), comparator);
-  
+
     compute_cell_wise(new_indices, reverse, dof, level, ordered_cells);
   }
 
@@ -1293,7 +1325,7 @@ namespace DoFRenumbering
                                          * Revert sorting order.
                                          */
        bool counter;
-    
+
                                         /**
                                          * Constructor.
                                          */
@@ -1308,7 +1340,7 @@ namespace DoFRenumbering
        bool operator () (const DHCellIterator& c1,
                          const DHCellIterator& c2) const
          {
-       
+
            const Point<dim> v1 = c1->center() - center;
            const Point<dim> v2 = c2->center() - center;
            const double s1 = std::atan2(v1(0), v1(1));
@@ -1318,7 +1350,7 @@ namespace DoFRenumbering
     };
   }
 
-  
+
 
   template <class DH, int dim>
   void
@@ -1346,10 +1378,10 @@ namespace DoFRenumbering
     std::vector<typename DH::cell_iterator>
       ordered_cells(dof.get_tria().n_active_cells());
     internal::ClockCells<dim> comparator(center, counter);
-  
+
     typename DH::active_cell_iterator begin = dof.begin_active();
     typename DH::active_cell_iterator end = dof.end();
-  
+
     std::copy (begin, end, ordered_cells.begin());
     std::sort (ordered_cells.begin(), ordered_cells.end(), comparator);
 
@@ -1357,7 +1389,7 @@ namespace DoFRenumbering
     compute_cell_wise_dg(new_indices, reverse, dof, ordered_cells);
   }
 
-  
+
 
   template <int dim>
   void clockwise_dg (MGDoFHandler<dim>& dof,
@@ -1368,10 +1400,10 @@ namespace DoFRenumbering
     std::vector<typename MGDoFHandler<dim>::cell_iterator>
       ordered_cells(dof.get_tria().n_cells(level));
     internal::ClockCells<dim> comparator(center, counter);
-  
+
     typename MGDoFHandler<dim>::cell_iterator begin = dof.begin(level);
     typename MGDoFHandler<dim>::cell_iterator end = dof.end(level);
-  
+
     std::copy (begin, end, ordered_cells.begin());
     std::sort (ordered_cells.begin(), ordered_cells.end(), comparator);
 
@@ -1391,7 +1423,7 @@ namespace DoFRenumbering
     dof_handler.renumber_dofs(renumbering);
   }
 
-  
+
 
   template <class DH>
   void
@@ -1402,10 +1434,10 @@ namespace DoFRenumbering
     const unsigned int n_dofs = dof_handler.n_dofs();
     Assert(new_indices.size() == n_dofs,
           ExcDimensionMismatch(new_indices.size(), n_dofs));
-  
+
     for (unsigned i=0; i<n_dofs; ++i)
       new_indices[i] = i;
-  
+
     std::random_shuffle (new_indices.begin(), new_indices.end());
   }
 
@@ -1422,7 +1454,7 @@ namespace DoFRenumbering
     dof_handler.renumber_dofs(renumbering);
   }
 
-  
+
 
   template <class DH>
   void
@@ -1442,7 +1474,7 @@ namespace DoFRenumbering
     const unsigned int n_subdomains
       = *std::max_element (subdomain_association.begin(),
                           subdomain_association.end()) + 1;
-  
+
                                     // then renumber the subdomains by first
                                     // looking at those belonging to subdomain
                                     // 0, then those of subdomain 1, etc. note
@@ -1505,7 +1537,7 @@ namespace DoFRenumbering
     compute_minimum_degree (std::vector<unsigned int> &,
                            const DoFHandler<deal_II_dimension> &, bool, bool);
   }
-  
+
 
 // non-boost functions:
 
@@ -1548,6 +1580,11 @@ namespace DoFRenumbering
   (DoFHandler<deal_II_dimension>&,
    const std::vector<unsigned int>&);
 
+  template
+  void component_wise<deal_II_dimension>
+  (hp::DoFHandler<deal_II_dimension>&,
+   const std::vector<unsigned int>&);
+
 //  template
 //  void
 //  compute_component_wise<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.