From: heister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Date: Wed, 24 Aug 2011 16:36:32 +0000 (+0000)
Subject: Added GridRefinement:hierarchical() to reorder the degrees of freedom by going throug... 
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=f9deadee13e67d26a4e549c69f0b65546ecd7127;p=dealii-svn.git

Added GridRefinement:hierarchical() to reorder the degrees of freedom by going through the cells in hierarchical order. This ensures consistent DoF numbering in parallel computations.

git-svn-id: https://svn.dealii.org/trunk@24178 0785d39b-7218-0410-832d-ea1e28bc413d
---

diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h
index c167c3b976..084f3a35a8 100644
--- a/deal.II/doc/news/changes.h
+++ b/deal.II/doc/news/changes.h
@@ -281,6 +281,12 @@ and DoF handlers embedded in higher dimensional space have been added.
 <h3>Specific improvements</h3>
 
 <ol>
+<li> New: Added GridRefinement:hierarchical() to reorder the degrees of freedom
+by going through the cells in hierarchical order. This ensures consistent
+DoF numbering in parallel computations.
+<br>
+(Timo Heister, 2011/08/24)
+
 <li> Changed: Triangulation<dim>::get_boundary_indicators() returned
 wrong data for dim=1.
 <br>
diff --git a/deal.II/include/deal.II/dofs/dof_renumbering.h b/deal.II/include/deal.II/dofs/dof_renumbering.h
index 925a57ca2a..7b8500f2ae 100644
--- a/deal.II/include/deal.II/dofs/dof_renumbering.h
+++ b/deal.II/include/deal.II/dofs/dof_renumbering.h
@@ -790,6 +790,17 @@ namespace DoFRenumbering
 			  const ENDITERATOR& end,
 			  const std::vector<unsigned int> &target_component);
 
+				   /**
+					* Renumber the degrees cell by cell in hierarchical order
+					* (also known as z-order). The main usage is that this
+					* guarantees the same ordering independent of the
+					* number of processors involved in a parallel
+					* distributed computation.
+					*/
+  template <int dim>
+  void
+  hierarchical (DoFHandler<dim> &dof_handler);
+  
 				   /**
 				    * Cell-wise renumbering for DG
 				    * elements.  This function takes
diff --git a/deal.II/source/dofs/dof_renumbering.cc b/deal.II/source/dofs/dof_renumbering.cc
index 6affc7c7ac..7af3b620e2 100644
--- a/deal.II/source/dofs/dof_renumbering.cc
+++ b/deal.II/source/dofs/dof_renumbering.cc
@@ -899,6 +899,96 @@ namespace DoFRenumbering
     return next_free_index;
   }
 
+  namespace
+  {
+	// helper function for hierarchical()
+    template <int dim, class iterator>
+    unsigned int
+    compute_hierarchical_recursive (
+      unsigned int next_free,
+      std::vector<unsigned int>& new_indices,
+      const iterator & cell,
+      const IndexSet & locally_owned)
+    {
+      if (cell->has_children())
+        {
+		  //recursion
+          for (unsigned int c = 0;c < GeometryInfo<dim>::max_children_per_cell; ++c)
+            next_free = compute_hierarchical_recursive<dim> (
+                          next_free,
+                          new_indices,
+                          cell->child (c),
+                          locally_owned);
+        }
+      else
+        {
+          if (!cell->is_artificial() && !cell->is_ghost())
+            {
+              const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell;
+              std::vector<unsigned int> local_dof_indices (dofs_per_cell);
+              cell->get_dof_indices (local_dof_indices);
+
+              for (unsigned int i = 0;i < dofs_per_cell;++i)
+                {
+                  if (locally_owned.is_element (local_dof_indices[i]))
+                    {
+                      // this is a locally owned DoF, assign new number if not assigned a number yet
+                      unsigned int idx = locally_owned.index_within_set (local_dof_indices[i]);
+                      if (new_indices[idx] == DoFHandler<dim>::invalid_dof_index)
+                        {
+                          new_indices[idx] = locally_owned.nth_index_in_set (next_free);
+                          next_free++;
+                        }
+                    }
+                }
+            }
+        }
+      return next_free;
+    }
+  }
+
+
+
+  template <int dim>
+  void
+  hierarchical (DoFHandler<dim> &dof_handler)
+  {
+    std::vector<unsigned int> renumbering (dof_handler.n_locally_owned_dofs(),
+					   DoFHandler<dim>::invalid_dof_index);
+
+	typename DoFHandler<dim>::cell_iterator cell;
+
+	unsigned int next_free = 0;
+	const IndexSet locally_owned = dof_handler.locally_owned_dofs();
+
+    for(cell = dof_handler.begin(0); cell != dof_handler.end(0); ++cell)
+	  next_free = compute_hierarchical_recursive<dim> (next_free,
+													   renumbering,
+													   cell,
+													   locally_owned);
+
+				     // verify that the last numbered
+				     // degree of freedom is either
+				     // equal to the number of degrees
+				     // of freedom in total (the
+				     // sequential case) or in the
+				     // distributed case at least
+				     // makes sense
+    Assert ((next_free == dof_handler.n_locally_owned_dofs())
+	    ||
+	    ((dof_handler.n_locally_owned_dofs() < dof_handler.n_dofs())
+	     &&
+	     (next_free <= dof_handler.n_dofs())),
+	    ExcRenumberingIncomplete());
+
+	// make sure that all local DoFs got new numbers assigned
+	Assert (std::find (renumbering.begin(), renumbering.end(),
+		       numbers::invalid_unsigned_int)
+	    == renumbering.end(),
+	    ExcInternalError());
+	dof_handler.renumber_dofs(renumbering);
+  }
+
 
 
   template <class DH>
diff --git a/deal.II/source/dofs/dof_renumbering.inst.in b/deal.II/source/dofs/dof_renumbering.inst.in
index c633fa06e2..1db0cd32fe 100644
--- a/deal.II/source/dofs/dof_renumbering.inst.in
+++ b/deal.II/source/dofs/dof_renumbering.inst.in
@@ -171,6 +171,11 @@ for (deal_II_dimension : DIMENSIONS)
   (MGDoFHandler<deal_II_dimension>&,
    const std::vector<unsigned int>&);
 
+  template
+  void hierarchical<deal_II_dimension>
+  (DoFHandler<deal_II_dimension>&);
+
+
 // DG renumbering for DoFHandler
 
   template