]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
New compute_* functions, that compute and return the renumbering vectors, but do...
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 3 Feb 2003 10:19:47 +0000 (10:19 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 3 Feb 2003 10:19:47 +0000 (10:19 +0000)
git-svn-id: https://svn.dealii.org/trunk@7002 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/dof_renumbering.h
deal.II/deal.II/source/dofs/dof_renumbering.cc
deal.II/doc/news/2002/c-3-4.html

index 8847d9502ad0c3834233dd480fa339a5d4310130..bbc8d57e6fc018bbdacebf5e5aab44a2af234ca1 100644 (file)
@@ -216,6 +216,23 @@ class DoFRenumbering
                   const bool                  use_constraints    = false,
                   const std::vector<unsigned int> &starting_indices   = std::vector<unsigned int>());
 
+                                    /**
+                                     * Computes the renumbering
+                                     * vector needed by the
+                                     * @p{Cuthill_McKee}
+                                     * function. Does not perform the
+                                     * renumbering on the
+                                     * @p{DoFHandler} dofs but
+                                     * returns the renumbering
+                                     * vector.
+                                     */    
+    template <int dim>
+    static std::vector<unsigned int>
+    compute_Cuthill_McKee (DoFHandler<dim>                 &dof_handler,
+                          const bool                       reversed_numbering = false,
+                          const bool                       use_constraints    = false,
+                          const std::vector<unsigned int> &starting_indices   = std::vector<unsigned int>());
+
                                     /**
                                      * Renumber the degrees of freedom
                                      * according to the Cuthill-McKee method,
@@ -294,9 +311,24 @@ class DoFRenumbering
                                      */
     template <int dim>
     static void
-    component_wise (DoFHandler<dim>            &dof_handler,
+    component_wise (DoFHandler<dim>                 &dof_handler,
                    const std::vector<unsigned int> &component_order = std::vector<unsigned int>());
 
+                                    /**
+                                     * Computes the renumbering
+                                     * vector needed by the
+                                     * @p{component_wise}
+                                     * function. Does not perform the
+                                     * renumbering on the
+                                     * @p{DoFHandler} dofs but
+                                     * returns the renumbering
+                                     * vector.
+                                     */    
+    template <int dim>
+    static std::vector<unsigned int>
+    compute_component_wise (DoFHandler<dim>                 &dof_handler,
+                           const std::vector<unsigned int> &component_order = std::vector<unsigned int>());
+    
                                     /**
                                      * Sort the degrees of freedom by
                                      * component. It does the same
@@ -332,9 +364,23 @@ class DoFRenumbering
                                      */
     template <int dim>
     static void
-    cell_wise_dg (DoFHandler<dim>                     &dof_handler,
+    cell_wise_dg (DoFHandler<dim> &dof_handler,
                  const std::vector<typename DoFHandler<dim>::cell_iterator> &cell_order);
-    
+
+                                    /**
+                                     * Computes the renumbering
+                                     * vector needed by the
+                                     * @p{cell_wise_dg}
+                                     * function. Does not perform the
+                                     * renumbering on the
+                                     * @p{DoFHandler} dofs but
+                                     * returns the renumbering
+                                     * vector.
+                                     */
+    template <int dim>
+    static std::vector<unsigned int>
+    compute_cell_wise_dg (DoFHandler<dim> &dof_handler,
+                         const std::vector<typename DoFHandler<dim>::cell_iterator> &cell_order);
 
                                     /**
                                      * Cell-wise renumbering on one
@@ -381,6 +427,21 @@ class DoFRenumbering
     downstream_dg (DoFHandler<dim>  &dof_handler,
                   const Point<dim> &direction);
 
+                                    /**
+                                     * Computes the renumbering
+                                     * vector needed by the
+                                     * @p{downstream_dg}
+                                     * function. Does not perform the
+                                     * renumbering on the
+                                     * @p{DoFHandler} dofs but
+                                     * returns the renumbering
+                                     * vector.
+                                     */
+    template <int dim>
+    static std::vector<unsigned int>
+    compute_downstream_dg (DoFHandler<dim>  &dof_handler,
+                          const Point<dim> &direction);
+
                                     /**
                                      * Cell-wise downstream numbering
                                      * with respect to a constant
@@ -412,6 +473,21 @@ class DoFRenumbering
     sort_selected_dofs_back (DoFHandler<dim>         &dof_handler,
                             const std::vector<bool> &selected_dofs);
 
+                                    /**
+                                     * Computes the renumbering
+                                     * vector needed by the
+                                     * @p{sort_selected_dofs_back}
+                                     * function. Does not perform the
+                                     * renumbering on the
+                                     * @p{DoFHandler} dofs but
+                                     * returns the renumbering
+                                     * vector.
+                                     */
+    template <int dim>
+    static std::vector<unsigned int>
+    compute_sort_selected_dofs_back (DoFHandler<dim>         &dof_handler,
+                                    const std::vector<bool> &selected_dofs);
+
                                     /**
                                      * Renumber the degrees of
                                      * freedom in a random way.
@@ -419,7 +495,20 @@ class DoFRenumbering
     template <int dim>
     static void
     random (DoFHandler<dim> &dof_handler);
-    
+
+                                    /**
+                                     * Computes the renumbering
+                                     * vector needed by the
+                                     * @p{random} function. Does not
+                                     * perform the renumbering on the
+                                     * @p{DoFHandler} dofs but
+                                     * returns the renumbering
+                                     * vector.
+                                     */   
+    template <int dim>
+    static std::vector<unsigned int>
+    compute_random (DoFHandler<dim> &dof_handler);
+
                                     /**
                                      * Exception
                                      */
index 5f6b159ec7adcbb30416b934e9bafb8a6a226de9..7239aad6236f598859a3a514c4e91f4ac3a786b5 100644 (file)
@@ -50,11 +50,32 @@ extern "C" long int lrand48 (void);
 #endif
 
 
+
 template <int dim>
-void DoFRenumbering::Cuthill_McKee (DoFHandler<dim>                 &dof_handler,
-                                   const bool                       reversed_numbering,
-                                   const bool                       use_constraints,
-                                   const std::vector<unsigned int> &starting_indices)
+void
+DoFRenumbering::Cuthill_McKee (DoFHandler<dim>                 &dof_handler,
+                              const bool                       reversed_numbering,
+                              const bool                       use_constraints,
+                              const std::vector<unsigned int> &starting_indices)
+{
+  std::vector<unsigned int> renumbering
+    =compute_Cuthill_McKee(dof_handler, reversed_numbering,
+                          use_constraints, starting_indices);
+
+                                  // actually perform renumbering;
+                                  // this is dimension specific and
+                                  // thus needs an own function
+  dof_handler.renumber_dofs (renumbering);
+}
+
+
+
+template <int dim>
+std::vector<unsigned int>
+DoFRenumbering::compute_Cuthill_McKee (DoFHandler<dim>                 &dof_handler,
+                                      const bool                       reversed_numbering,
+                                      const bool                       use_constraints,
+                                      const std::vector<unsigned int> &starting_indices)
 {
                                   // make the connection graph
   SparsityPattern sparsity (dof_handler.n_dofs(),
@@ -253,10 +274,7 @@ void DoFRenumbering::Cuthill_McKee (DoFHandler<dim>                 &dof_handler
         i!=new_number.end(); ++i)
       *i = n_dofs-*i-1;
 
-                                  // actually perform renumbering;
-                                  // this is dimension specific and
-                                  // thus needs an own function
-  dof_handler.renumber_dofs (new_number);
+  return new_number;
 }
 
 
@@ -432,14 +450,34 @@ template <int dim>
 void
 DoFRenumbering::component_wise (DoFHandler<dim>                 &dof_handler,
                                 const std::vector<unsigned int> &component_order_arg)
+{
+  std::vector<unsigned int> renumbering
+    =compute_component_wise(dof_handler, component_order_arg);
+
+  if (renumbering.size()!=0)
+    dof_handler.renumber_dofs (renumbering);
+}
+  
+
+
+template <int dim>
+std::vector<unsigned int>
+DoFRenumbering::compute_component_wise (DoFHandler<dim>                 &dof_handler,
+                                       const std::vector<unsigned int> &component_order_arg)
 {
   const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell;
   const FiniteElement<dim> &fe     = dof_handler.get_fe();
 
+  std::vector<unsigned int> new_indices (dof_handler.n_dofs(),
+                                        DoFHandler<dim>::invalid_dof_index);
+
                                   // do nothing if the FE has only
                                   // one component
   if (fe.n_components() == 1)
-    return;
+    {
+      new_indices.resize(0);
+      return new_indices;
+    }
   
   std::vector<unsigned int> component_order (component_order_arg);
   if (component_order.size() == 0)
@@ -556,8 +594,6 @@ DoFRenumbering::component_wise (DoFHandler<dim>                 &dof_handler,
                                    // components in the order the user
                                    // desired to see
   unsigned int next_free_index = 0;
-  std::vector<unsigned int> new_indices (dof_handler.n_dofs(),
-                                        DoFHandler<dim>::invalid_dof_index);
   for (unsigned int c=0; c<fe.n_components(); ++c)
     {
       const unsigned int component = component_order[c];
@@ -575,7 +611,7 @@ DoFRenumbering::component_wise (DoFHandler<dim>                 &dof_handler,
   Assert (next_free_index == dof_handler.n_dofs(),
          ExcInternalError());
 
-  dof_handler.renumber_dofs (new_indices);
+  return new_indices;
 }
 
 
@@ -735,6 +771,19 @@ template <int dim>
 void
 DoFRenumbering::sort_selected_dofs_back (DoFHandler<dim>         &dof_handler,
                                         const std::vector<bool> &selected_dofs)
+{
+  std::vector<unsigned int> renumbering
+    =compute_sort_selected_dofs_back(dof_handler, selected_dofs);
+
+  dof_handler.renumber_dofs(renumbering);
+}
+
+
+
+template <int dim>
+std::vector<unsigned int>
+DoFRenumbering::compute_sort_selected_dofs_back (DoFHandler<dim>         &dof_handler,
+                                                const std::vector<bool> &selected_dofs)
 {
   const unsigned int n_dofs = dof_handler.n_dofs();
   Assert (selected_dofs.size() == n_dofs,
@@ -763,15 +812,30 @@ DoFRenumbering::sort_selected_dofs_back (DoFHandler<dim>         &dof_handler,
   Assert (next_unselected == n_selected_dofs, ExcInternalError());
   Assert (next_selected == n_dofs, ExcInternalError());
 
-                                  // now perform the renumbering
-  dof_handler.renumber_dofs (new_dof_indices);
+  return new_dof_indices;
 }
 
 
+
 template <int dim>
 void
-DoFRenumbering::cell_wise_dg (DoFHandler<dim>& dof,
-                             const typename std::vector<typename DoFHandler<dim>::cell_iterator>& cells)
+DoFRenumbering::cell_wise_dg (
+  DoFHandler<dim>& dof,
+  const typename std::vector<typename DoFHandler<dim>::cell_iterator>& cells)
+{
+  std::vector<unsigned int> renumbering
+    =compute_cell_wise_dg(dof, cells);
+  
+  dof.renumber_dofs(renumbering);
+}
+
+
+
+template <int dim>
+std::vector<unsigned int>
+DoFRenumbering::compute_cell_wise_dg (
+  DoFHandler<dim>& dof,
+  const typename std::vector<typename DoFHandler<dim>::cell_iterator>& cells)
 {
   Assert(cells.size() == dof.get_tria().n_active_cells(),
         ExcDimensionMismatch(cells.size(),
@@ -816,16 +880,15 @@ DoFRenumbering::cell_wise_dg (DoFHandler<dim>& dof,
   for (unsigned int i=0;i<new_order.size(); ++i)
     reverse[new_order[i]] = i;
 
-  dof.renumber_dofs(reverse);
+  return reverse;
 }
 
 
 #ifdef ENABLE_MULTIGRID
 template <int dim>
-void
-DoFRenumbering::cell_wise_dg (MGDoFHandler<dim>& dof,
-                             const unsigned int level,
-                             const typename std::vector<typename MGDoFHandler<dim>::cell_iterator>& cells)
+void DoFRenumbering::cell_wise_dg (MGDoFHandler<dim>& dof,
+                                  const unsigned int level,
+                                  const typename std::vector<typename MGDoFHandler<dim>::cell_iterator>& cells)
 {
   Assert(cells.size() == dof.get_tria().n_cells(level),
         ExcDimensionMismatch(cells.size(),
@@ -905,11 +968,24 @@ struct CompCells
       }
 };
 
-    
+
 template <int dim>
 void
 DoFRenumbering::downstream_dg (DoFHandler<dim>& dof,
                               const Point<dim>& direction)
+{
+  std::vector<unsigned int> renumbering
+    =compute_downstream_dg(dof, direction);
+
+  dof.renumber_dofs(renumbering);
+}
+
+
+
+template <int dim>
+std::vector<unsigned int>
+DoFRenumbering::compute_downstream_dg (DoFHandler<dim>& dof,
+                                      const Point<dim>& direction)
 {
   std::vector<typename DoFHandler<dim>::cell_iterator>
     ordered_cells(dof.get_tria().n_active_cells());
@@ -921,16 +997,15 @@ DoFRenumbering::downstream_dg (DoFHandler<dim>& dof,
   copy (begin, end, ordered_cells.begin());
   sort (ordered_cells.begin(), ordered_cells.end(), comparator);
 
-  cell_wise_dg(dof, ordered_cells);
+  return compute_cell_wise_dg(dof, ordered_cells);
 }
 
 
 #ifdef ENABLE_MULTIGRID
 template <int dim>
-void
-DoFRenumbering::downstream_dg (MGDoFHandler<dim>& dof,
-                              const unsigned int level,
-                              const Point<dim>& direction)
+void DoFRenumbering::downstream_dg (MGDoFHandler<dim>& dof,
+                                   const unsigned int level,
+                                   const Point<dim>& direction)
 {
   std::vector<typename MGDoFHandler<dim>::cell_iterator>
     ordered_cells(dof.get_tria().n_cells(level));
@@ -947,9 +1022,21 @@ DoFRenumbering::downstream_dg (MGDoFHandler<dim>& dof,
 #endif
 
 
+
 template <int dim>
 void
 DoFRenumbering::random (DoFHandler<dim> &dof_handler)
+{
+  std::vector<unsigned int> renumbering
+    =compute_random(dof_handler);
+
+  dof_handler.renumber_dofs(renumbering);
+}
+
+
+template <int dim>
+std::vector<unsigned int>
+DoFRenumbering::compute_random (DoFHandler<dim> &dof_handler)
 {
   const unsigned int n_dofs = dof_handler.n_dofs();
   
@@ -958,7 +1045,7 @@ DoFRenumbering::random (DoFHandler<dim> &dof_handler)
     new_indices[i] = i;
   
   std::random_shuffle (new_indices.begin(), new_indices.end());
-  dof_handler.renumber_dofs(new_indices);  
+  return new_indices;
 }
 
 
@@ -972,11 +1059,25 @@ void DoFRenumbering::Cuthill_McKee<deal_II_dimension>
  const bool,
  const std::vector<unsigned int>&);
 
+template
+std::vector<unsigned int>
+DoFRenumbering::compute_Cuthill_McKee<deal_II_dimension>
+(DoFHandler<deal_II_dimension>&,
+ const bool,
+ const bool,
+ const std::vector<unsigned int>&);
+
 template
 void DoFRenumbering::component_wise<deal_II_dimension>
 (DoFHandler<deal_II_dimension>&,
  const std::vector<unsigned int>&);
 
+template
+std::vector<unsigned int>
+DoFRenumbering::compute_component_wise<deal_II_dimension>
+(DoFHandler<deal_II_dimension>&,
+ const std::vector<unsigned int>&);
+
 template
 void DoFRenumbering::component_wise<deal_II_dimension>
 (MGDoFHandler<deal_II_dimension>&,
@@ -985,12 +1086,26 @@ void DoFRenumbering::component_wise<deal_II_dimension>
 
 
 template
-void DoFRenumbering::cell_wise_dg<deal_II_dimension>
+void
+DoFRenumbering::cell_wise_dg<deal_II_dimension>
 (DoFHandler<deal_II_dimension>&,
  const std::vector<DoFHandler<deal_II_dimension>::cell_iterator>&);
 
 template
-void DoFRenumbering::downstream_dg<deal_II_dimension>
+std::vector<unsigned int>
+DoFRenumbering::compute_cell_wise_dg<deal_II_dimension>
+(DoFHandler<deal_II_dimension>&,
+ const std::vector<DoFHandler<deal_II_dimension>::cell_iterator>&);
+
+template
+void
+DoFRenumbering::downstream_dg<deal_II_dimension>
+(DoFHandler<deal_II_dimension>&,
+ const Point<deal_II_dimension>&);
+
+template
+std::vector<unsigned int>
+DoFRenumbering::compute_downstream_dg<deal_II_dimension>
 (DoFHandler<deal_II_dimension>&,
  const Point<deal_II_dimension>&);
 
@@ -999,10 +1114,21 @@ void DoFRenumbering::sort_selected_dofs_back<deal_II_dimension>
 (DoFHandler<deal_II_dimension> &,
  const std::vector<bool> &);
 
+template
+std::vector<unsigned int>
+DoFRenumbering::compute_sort_selected_dofs_back<deal_II_dimension>
+(DoFHandler<deal_II_dimension> &,
+ const std::vector<bool> &);
+
 template
 void DoFRenumbering::random<deal_II_dimension>
 (DoFHandler<deal_II_dimension> &);
 
+template
+std::vector<unsigned int>
+DoFRenumbering::compute_random<deal_II_dimension>
+(DoFHandler<deal_II_dimension> &);
+
 #ifdef ENABLE_MULTIGRID
 template
 void DoFRenumbering::Cuthill_McKee<deal_II_dimension>
@@ -1010,6 +1136,7 @@ void DoFRenumbering::Cuthill_McKee<deal_II_dimension>
  const unsigned int,
  const bool,
  const std::vector<unsigned int>&);
+
 template
 void DoFRenumbering::cell_wise_dg<deal_II_dimension>
 (MGDoFHandler<deal_II_dimension>&,
index 0798c328ad4f14279c30b340d0f7ca863420d252..74737544d63ef9f0074cc4e0587f3a5328d6c95d 100644 (file)
@@ -635,6 +635,18 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
 <h3>deal.II</h3>
 
 <ol>
+  <li> <p>
+       New: For each of the renumbering functions in the <code
+       class="class">DoFRenumbering</code> class there is now an
+       additional <code class="member">compute_*</code>
+       function. These new functions compute and return the
+       corresponding renumbering vectors but do not perform the actual
+       renumbering on the <code class="class">DoFHandler</code>
+       object. The behaviour of the old functions is not changed.
+       <br>
+       (RH 2003/02/03)
+       </p>
+
   <li> <p> 
        New: The <code class="class">ThreadMutex</code> classes now have a
        member class <code class="class">ScopedLock</code> that implements the

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.