]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
add cellwise numbering for continuous methods
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 30 Apr 2007 03:13:31 +0000 (03:13 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 30 Apr 2007 03:13:31 +0000 (03:13 +0000)
git-svn-id: https://svn.dealii.org/trunk@14646 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 00f9c250e0d72fdd3c737c9a6f6682394234555b..eb7c547e06c50229e71ca43f717fc993299f3024 100644 (file)
@@ -171,20 +171,25 @@ DEAL_II_NAMESPACE_OPEN
  * block.
  *
  *
- * <h3>Cell-wise numbering for Discontinuous Galerkin FEM</h3>
+ * <h3>Cell-wise numbering</h3>
  *
- * One advantage of DGFEM is the fact, that it yields invertible
- * blocks on the diagonal of the global matrix. these blocks are in
- * fact the cell matrices and matrix elements outside these blocks are
- * due to fluxes.
+ * Given an ordered vector of cells, the function cell_wise()
+ * sorts the degrees of freedom such that degrees on earlier cells of
+ * this vector will occur before degrees on later cells.
  *
- * Still, it may be necessary to apply a downstream numbering of the
- * degrees of freedom. This renumbering may only exchange whole blocks
- * and must not destroy the block structure.
+ * This rule produces a well-defined ordering for discontinuous Galerkin
+ * methods (FE_DGP, FE_DGQ). For continuous methods, we use the
+ * additional rule that each degree of freedom is ordered according to
+ * the first cell in the ordered vector it belongs to.
  *
- * Given an ordered vector of cells, the function cell_wise_dg()
- * accomplishes this. Inside the cells, the previous ordering will be
- * preserved, so it may be useful to apply component_wise() first.
+ * Applications of this scheme are downstream() and
+ * clock_wise_dg(). The first orders the cells according to a
+ * downstream direction and then applies cell_wise().
+ *
+ * @note For DG elements, the internal numbering in each cell remains
+ * unaffected. This cannot be guaranteed for continuous elements
+ * anymore, since degrees of freedom shared with an earlier cell will
+ * be accounted for by the other cell.
  *
  *
  * <h3>Random renumbering</h3>
@@ -203,7 +208,7 @@ DEAL_II_NAMESPACE_OPEN
  * information on this.
  *
  * @ingroup dofs
- * @author Wolfgang Bangerth, Guido Kanschat, 1998, 1999, 2000, 2004
+ * @author Wolfgang Bangerth, Guido Kanschat, 1998, 1999, 2000, 2004, 2007
  */
 class DoFRenumbering 
 {
@@ -403,6 +408,16 @@ class DoFRenumbering
                                      */
     template <class DH>
     static void
+    cell_wise (DH &dof_handler,
+              const std::vector<typename DH::cell_iterator> &cell_order);
+
+                                    /**
+                                     * @deprecated Use cell_wise() instead.
+                                     *
+                                     * Cell-wise numbering only for DG.
+                                     */
+    template <class DH>
+    static void
     cell_wise_dg (DH &dof_handler,
                  const std::vector<typename DH::cell_iterator> &cell_order);
 
@@ -423,6 +438,35 @@ class DoFRenumbering
                          const std::vector<typename DH::cell_iterator> &cell_order);
 
                                     /**
+                                     * Computes the renumbering
+                                     * vector needed by the
+                                     * cell_wise() function. Does
+                                     * not perform the renumbering on
+                                     * the DoFHandler dofs but
+                                     * returns the renumbering
+                                     * vector.
+                                     */
+    template <class DH>
+    static void
+    compute_cell_wise (std::vector<unsigned int>& renumbering,
+                      std::vector<unsigned int>& inverse_renumbering,
+                      const DH &dof_handler,
+                      const std::vector<typename DH::cell_iterator> &cell_order);
+
+                                    /**
+                                     * Cell-wise renumbering on one
+                                     * level. See the other function
+                                     * with the same name.
+                                     */
+    template <int dim>
+    static void
+    cell_wise (MGDoFHandler<dim>   &dof_handler,
+              const unsigned int   level,
+              const std::vector<typename MGDoFHandler<dim>::cell_iterator> &cell_order);
+    
+                                    /**
+                                     * @deprecated Use cell_wise() instead.
+                                     *
                                      * Cell-wise renumbering on one
                                      * level for DG elements. See the
                                      * other function with the same
@@ -452,6 +496,24 @@ class DoFRenumbering
                          const std::vector<typename MGDoFHandler<dim>::cell_iterator>& cell_order);
 
 
+                                    /**
+                                     * Computes the renumbering
+                                     * vector needed by the
+                                     * cell_wise() level renumbering function. Does
+                                     * not perform the renumbering on
+                                     * the MGDoFHandler dofs but
+                                     * returns the renumbering
+                                     * vector.
+                                     */
+    template <int dim>
+    static void
+    compute_cell_wise (std::vector<unsigned int>& renumbering,
+                      std::vector<unsigned int>& inverse_renumbering,
+                      const MGDoFHandler<dim>&   dof_handler,
+                      const unsigned int         level,
+                      const std::vector<typename MGDoFHandler<dim>::cell_iterator>& cell_order);
+
+
                                     /**
                                      * Cell-wise downstream numbering
                                      * with respect to a constant
@@ -481,6 +543,14 @@ class DoFRenumbering
                                      */
     template <class DH, int dim>
     static void
+    downstream (DH&               dof_handler,
+               const Point<dim>& direction);
+
+                                    /**
+                                     * @deprecated Use downstream() instead.
+                                     */
+    template <class DH, int dim>
+    static void
     downstream_dg (DH&               dof_handler,
                   const Point<dim>& direction);
 
@@ -493,6 +563,15 @@ class DoFRenumbering
                                      */
     template <int dim>
     static void
+    downstream (MGDoFHandler<dim>  &dof_handler,
+               const unsigned int level,
+               const Point<dim> &direction);
+                                    /**
+                                     * @deprecated Use downstream()
+                                     * instead.
+                                     */
+    template <int dim>
+    static void
     downstream_dg (MGDoFHandler<dim>  &dof_handler,
                   const unsigned int level,
                   const Point<dim> &direction);
@@ -530,6 +609,17 @@ class DoFRenumbering
                                      */
     template <class DH, int dim>
     static void
+    compute_downstream (std::vector<unsigned int>& new_dof_indices,
+                       std::vector<unsigned int>& reverse,
+                       const DH&                  dof_handler,
+                       const Point<dim>&          direction);
+    
+                                    /**
+                                     * @deprecated Use
+                                     * compute_downstream() instead
+                                     */
+    template <class DH, int dim>
+    static void
     compute_downstream_dg (std::vector<unsigned int>& new_dof_indices,
                           std::vector<unsigned int>& reverse,
                           const DH&                  dof_handler,
@@ -546,6 +636,18 @@ class DoFRenumbering
                                      */
     template <int dim>
     static void
+    compute_downstream (std::vector<unsigned int>& new_dof_indices,
+                       std::vector<unsigned int>& reverse,
+                       const MGDoFHandler<dim>&   dof_handler,
+                       const unsigned int         level,
+                       const Point<dim>&          direction);
+    
+                                    /**
+                                     * @deprecated Use
+                                     * compute_downstream() instead
+                                     */
+    template <int dim>
+    static void
     compute_downstream_dg (std::vector<unsigned int>& new_dof_indices,
                           std::vector<unsigned int>& reverse,
                           const MGDoFHandler<dim>&   dof_handler,
index be8047068cb24b4060b862c30adc6daad92e6e66..3a9f5ec507c1ffa693669fff6f19a7eab3d6679a 100644 (file)
@@ -889,7 +889,21 @@ DoFRenumbering::cell_wise_dg (
 }
 
 
+template <class DH>
+void
+DoFRenumbering::cell_wise (
+  DH& dof,
+  const std::vector<typename DH::cell_iterator>& cells)
+{
+  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);
+}
 
+
+//TODO: Discuss if we cannot replace this function by the next
 template <class DH>
 void
 DoFRenumbering::compute_cell_wise_dg (
@@ -948,6 +962,71 @@ DoFRenumbering::compute_cell_wise_dg (
 }
 
 
+template <class DH>
+void
+DoFRenumbering::compute_cell_wise (
+  std::vector<unsigned int>& new_indices,
+  std::vector<unsigned int>& reverse,
+  const DH& dof,
+  const typename std::vector<typename DH::cell_iterator>& cells)
+{
+  Assert(cells.size() == dof.get_tria().n_active_cells(),
+        ExcDimensionMismatch(cells.size(),
+                             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,
+        ExcDimensionMismatch(reverse.size(), n_global_dofs));
+
+                                  // For continuous elements, we must
+                                  // make sure, that each dof is
+                                  // reordered only once.
+  std::vector<bool> already_sorted(n_global_dofs, false);
+  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)
+    {
+                                      // Determine the number of dofs
+                                      // on this cell and reinit the
+                                      // vector storing these
+                                      // 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
+                                      // degrees of freedom inside a
+                                      // single cell are in the same
+                                      // order after renumbering.
+      std::sort(cell_dofs.begin(), cell_dofs.end());
+
+      for (unsigned int i=0;i<n_cell_dofs;++i)
+       {
+         if (!already_sorted[cell_dofs[i]])
+           {
+             already_sorted[cell_dofs[i]] = true;
+             reverse[global_index++] = cell_dofs[i];
+           }
+       }
+    }
+  Assert(global_index == n_global_dofs, ExcRenumberingIncomplete());
+
+  for (unsigned int i=0;i<reverse.size(); ++i)
+    new_indices[reverse[i]] = i;
+}
+
+
 template <int dim>
 void DoFRenumbering::cell_wise_dg (
   MGDoFHandler<dim>& dof,
@@ -958,7 +1037,21 @@ void DoFRenumbering::cell_wise_dg (
   std::vector<unsigned int> reverse(dof.n_dofs(level));
   
   compute_cell_wise_dg(renumbering, reverse, dof, level, cells);
-  dof.renumber_dofs(level, reverse);
+  dof.renumber_dofs(level, renumbering);
+}
+
+
+template <int dim>
+void DoFRenumbering::cell_wise (
+  MGDoFHandler<dim>& dof,
+  const unsigned int level,
+  const typename std::vector<typename MGDoFHandler<dim>::cell_iterator>& cells)
+{
+  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);
 }
 
 
@@ -1009,13 +1102,63 @@ void DoFRenumbering::compute_cell_wise_dg (
 
       for (unsigned int i=0;i<n_cell_dofs;++i)
        {
-         new_order[global_index++] = cell_dofs[i];
+         reverse[global_index++] = cell_dofs[i];
+       }
+    }
+  Assert(global_index == n_global_dofs, ExcRenumberingIncomplete());
+
+  for (unsigned int i=0;i<new_order.size(); ++i)
+    new_order[reverse[i]] = i;  
+}
+
+
+
+template <int dim>
+void DoFRenumbering::compute_cell_wise (
+  std::vector<unsigned int>& new_order,
+  std::vector<unsigned int>& reverse,
+  const 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(),
+                             dof.get_tria().n_cells(level)));
+  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<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)
+    {
+      Assert ((*cell)->level() == (int) level, ExcInternalError());
+
+      (*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]])
+           {
+             already_sorted[cell_dofs[i]] = true;
+             reverse[global_index++] = cell_dofs[i];
+           }
        }
     }
   Assert(global_index == n_global_dofs, ExcRenumberingIncomplete());
 
   for (unsigned int i=0;i<new_order.size(); ++i)
-    reverse[new_order[i]] = i;  
+    new_order[reverse[i]] = i;  
 }
 
 
@@ -1032,6 +1175,19 @@ DoFRenumbering::downstream_dg (DH& dof, const Point<dim>& direction)
 
 
 
+template <class DH, int dim>
+void
+DoFRenumbering::downstream (DH& dof, const Point<dim>& direction)
+{
+  std::vector<unsigned int> renumbering(dof.n_dofs());
+  std::vector<unsigned int> reverse(dof.n_dofs());
+  compute_downstream(renumbering, reverse, dof, direction);
+
+  dof.renumber_dofs(renumbering);
+}
+
+
+
 template <class DH, int dim>
 void
 DoFRenumbering::compute_downstream_dg (
@@ -1076,10 +1232,62 @@ DoFRenumbering::compute_downstream_dg (
 }
 
 
+template <class DH, int dim>
+void
+DoFRenumbering::compute_downstream (
+  std::vector<unsigned int>& new_indices,
+  std::vector<unsigned int>& reverse,
+  const DH& dof,
+  const Point<dim>& direction)
+{
+  std::vector<typename DH::cell_iterator>
+    ordered_cells(dof.get_tria().n_active_cells());
+  const CompareDownstream<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());
+  sort (ordered_cells.begin(), ordered_cells.end(), comparator);
+  
+  compute_cell_wise(new_indices, reverse, dof, ordered_cells);
+}
+
+
 template <int dim>
 void DoFRenumbering::downstream_dg (MGDoFHandler<dim>& dof,
                                    const unsigned int level,
                                    const Point<dim>& direction)
+{
+  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 DoFRenumbering::downstream (MGDoFHandler<dim>& dof,
+                                const unsigned int level,
+                                const Point<dim>& direction)
+{
+  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
+DoFRenumbering::compute_downstream_dg (
+  std::vector<unsigned int>& new_indices,
+  std::vector<unsigned int>& reverse,
+  const 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));
@@ -1088,16 +1296,16 @@ void DoFRenumbering::downstream_dg (MGDoFHandler<dim>& dof,
   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);
-
-  cell_wise_dg (dof, level, ordered_cells);
+  copy (begin, end, ordered_cells.begin());
+  sort (ordered_cells.begin(), ordered_cells.end(), comparator);
+  
+  compute_cell_wise_dg(new_indices, reverse, dof, level, ordered_cells);
 }
 
 
 template <int dim>
 void
-DoFRenumbering::compute_downstream_dg (
+DoFRenumbering::compute_downstream (
   std::vector<unsigned int>& new_indices,
   std::vector<unsigned int>& reverse,
   const MGDoFHandler<dim>& dof,
@@ -1114,7 +1322,7 @@ DoFRenumbering::compute_downstream_dg (
   copy (begin, end, ordered_cells.begin());
   sort (ordered_cells.begin(), ordered_cells.end(), comparator);
   
-  compute_cell_wise_dg(new_indices, reverse, dof, level, ordered_cells);
+  compute_cell_wise(new_indices, reverse, dof, level, ordered_cells);
 }
 
 
@@ -1403,43 +1611,67 @@ DoFRenumbering::compute_downstream_dg<DoFHandler<deal_II_dimension> >
  const DoFHandler<deal_II_dimension>&,
  const Point<deal_II_dimension>&);
 
+template void
+DoFRenumbering::cell_wise<DoFHandler<deal_II_dimension> >
+(DoFHandler<deal_II_dimension>&,
+ const std::vector<DoFHandler<deal_II_dimension>::cell_iterator>&);
+
+template void
+DoFRenumbering::compute_cell_wise<DoFHandler<deal_II_dimension> >
+(std::vector<unsigned int>&, std::vector<unsigned int>&,
+ const DoFHandler<deal_II_dimension>&,
+ const std::vector<DoFHandler<deal_II_dimension>::cell_iterator>&);
+
+template void
+DoFRenumbering::downstream<DoFHandler<deal_II_dimension> >
+(DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&);
+
+template void
+DoFRenumbering::compute_downstream<DoFHandler<deal_II_dimension> >
+(std::vector<unsigned int>&,std::vector<unsigned int>&,
+ const DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&);
+
 template
 void
 DoFRenumbering::clockwise_dg<DoFHandler<deal_II_dimension> >
-(DoFHandler<deal_II_dimension>&,
- const Point<deal_II_dimension>&, bool);
+(DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&, bool);
 
 template
 void
 DoFRenumbering::compute_clockwise_dg<DoFHandler<deal_II_dimension> >
-(std::vector<unsigned int>&,
- const DoFHandler<deal_II_dimension>&,
- const Point<deal_II_dimension>&,
- const bool);
+(std::vector<unsigned int>&, const DoFHandler<deal_II_dimension>&,
+ const Point<deal_II_dimension>&, const bool);
 
-// DG renumbering for hp::DoFHandler
+// Renumbering for hp::DoFHandler
 
-template
-void
+template void
 DoFRenumbering::cell_wise_dg<hp::DoFHandler<deal_II_dimension> >
 (hp::DoFHandler<deal_II_dimension>&,
  const std::vector<hp::DoFHandler<deal_II_dimension>::cell_iterator>&);
 
-template
-void
+template void
 DoFRenumbering::compute_cell_wise_dg<hp::DoFHandler<deal_II_dimension> >
 (std::vector<unsigned int>&, std::vector<unsigned int>&,
  const hp::DoFHandler<deal_II_dimension>&,
  const std::vector<hp::DoFHandler<deal_II_dimension>::cell_iterator>&);
 
-template
-void
+template void
+DoFRenumbering::cell_wise<hp::DoFHandler<deal_II_dimension> >
+(hp::DoFHandler<deal_II_dimension>&,
+ const std::vector<hp::DoFHandler<deal_II_dimension>::cell_iterator>&);
+
+template void
+DoFRenumbering::compute_cell_wise<hp::DoFHandler<deal_II_dimension> >
+(std::vector<unsigned int>&, std::vector<unsigned int>&,
+ const hp::DoFHandler<deal_II_dimension>&,
+ const std::vector<hp::DoFHandler<deal_II_dimension>::cell_iterator>&);
+
+template void
 DoFRenumbering::downstream_dg<hp::DoFHandler<deal_II_dimension> >
 (hp::DoFHandler<deal_II_dimension>&,
  const Point<deal_II_dimension>&);
 
-template
-void
+template void
 DoFRenumbering::compute_downstream_dg<hp::DoFHandler<deal_II_dimension> >
 (std::vector<unsigned int>&,std::vector<unsigned int>&,
  const hp::DoFHandler<deal_II_dimension>&,
@@ -1451,6 +1683,17 @@ DoFRenumbering::compute_downstream_dg<hp::DoFHandler<deal_II_dimension> >
  const hp::DoFHandler<deal_II_dimension>&,
  const Point<deal_II_dimension>&);
 
+template void
+DoFRenumbering::downstream<hp::DoFHandler<deal_II_dimension> >
+(hp::DoFHandler<deal_II_dimension>&,
+ const Point<deal_II_dimension>&);
+
+template void
+DoFRenumbering::compute_downstream<hp::DoFHandler<deal_II_dimension> >
+(std::vector<unsigned int>&,std::vector<unsigned int>&,
+ const hp::DoFHandler<deal_II_dimension>&,
+ const Point<deal_II_dimension>&);
+
 template
 void
 DoFRenumbering::clockwise_dg<hp::DoFHandler<deal_II_dimension> >
@@ -1469,8 +1712,12 @@ DoFRenumbering::compute_clockwise_dg<hp::DoFHandler<deal_II_dimension> >
 
 template
 void DoFRenumbering::downstream_dg<deal_II_dimension>
-(MGDoFHandler<deal_II_dimension>&,
- const unsigned int,
+(MGDoFHandler<deal_II_dimension>&, const unsigned int,
+ const Point<deal_II_dimension>&);
+
+template
+void DoFRenumbering::downstream<deal_II_dimension>
+(MGDoFHandler<deal_II_dimension>&, const unsigned int,
  const Point<deal_II_dimension>&);
 
 template

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.