From: Timo Heister 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-Tag: v8.0.0~3590 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=da548ba2ac66833e807b4907fb1429cd1e42f074;p=dealii.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.

Specific improvements

    +
  1. 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. +
    +(Timo Heister, 2011/08/24) +
  2. Changed: Triangulation::get_boundary_indicators() returned wrong data for dim=1.
    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 &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 + void + hierarchical (DoFHandler &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 + unsigned int + compute_hierarchical_recursive ( + unsigned int next_free, + std::vector& new_indices, + const iterator & cell, + const IndexSet & locally_owned) + { + if (cell->has_children()) + { + //recursion + for (unsigned int c = 0;c < GeometryInfo::max_children_per_cell; ++c) + next_free = compute_hierarchical_recursive ( + 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 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::invalid_dof_index) + { + new_indices[idx] = locally_owned.nth_index_in_set (next_free); + next_free++; + } + } + } + } + } + return next_free; + } + } + + + + template + void + hierarchical (DoFHandler &dof_handler) + { + std::vector renumbering (dof_handler.n_locally_owned_dofs(), + DoFHandler::invalid_dof_index); + + typename DoFHandler::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 (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 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&, const std::vector&); + template + void hierarchical + (DoFHandler&); + + // DG renumbering for DoFHandler template