*/
template <int dim>
static void
- compute_Cuthill_McKee (std::vector<unsigned int>&,
+ compute_Cuthill_McKee (std::vector<unsigned int>& new_dof_indices,
const DoFHandler<dim>&,
const bool reversed_numbering = false,
const bool use_constraints = false,
*/
template <int dim, class ITERATOR, class ENDITERATOR>
static unsigned int
- compute_component_wise (std::vector<unsigned int>& new_index,
+ compute_component_wise (std::vector<unsigned int>& new_dof_indices,
ITERATOR& start,
const ENDITERATOR& end,
const std::vector<unsigned int> &target_component);
*/
template <int dim>
static void
- compute_cell_wise_dg (std::vector<unsigned int>&,
+ compute_cell_wise_dg (std::vector<unsigned int> &new_dof_indices,
const DoFHandler<dim> &dof_handler,
const std::vector<typename DoFHandler<dim>::cell_iterator> &cell_order);
*/
template <int dim>
static void
- compute_downstream_dg (std::vector<unsigned int>&,
+ compute_downstream_dg (std::vector<unsigned int> &new_dof_indices,
const DoFHandler<dim> &dof_handler,
const Point<dim> &direction);
*/
template <int dim>
static void
- compute_clockwise_dg (std::vector<unsigned int>&,
+ compute_clockwise_dg (std::vector<unsigned int> &new_dof_indices,
const DoFHandler<dim> &dof_handler,
const Point<dim> ¢er,
const bool counter);
*/
template <int dim>
static void
- compute_sort_selected_dofs_back (std::vector<unsigned int>&,
+ compute_sort_selected_dofs_back (std::vector<unsigned int> &new_dof_indices,
const DoFHandler<dim>&,
const std::vector<bool>&selected_dofs);
*/
template <int dim>
static void
- compute_random (std::vector<unsigned int>&,
+ compute_random (std::vector<unsigned int> &new_dof_indices,
const DoFHandler<dim> &dof_handler);
+ /**
+ * Renumber the degrees of freedom such
+ * that they are associated with the
+ * subdomain id of the cells they are
+ * living on, i.e. first all degrees of
+ * freedom that belong to cells with
+ * subdomain zero, then all with
+ * subdomain one, etc. This is useful
+ * when doing parallel computations after
+ * assigning subdomain ids using a
+ * partitioner (see the
+ * @p{GridTools::partition_triangulation}
+ * function for this).
+ *
+ * Note that degrees of freedom
+ * associated with faces, edges, and
+ * vertices may be associated with
+ * multiple subdomains if they are
+ * sitting on partition boundaries. It
+ * would therefore be undefined with
+ * which subdomain they have to be
+ * associated. For this, we use what we
+ * get from the
+ * @p{DoFTools::get_subdomain_association}
+ * function.
+ *
+ * The algorithm is stable, i.e. if
+ * two dofs i,j have @p{i<j} and belong
+ * to the same subdomain, then they
+ * will be in this order also after
+ * reordering.
+ */
+ template <int dim>
+ static void
+ subdomain_wise (DoFHandler<dim> &dof_handler);
+
+ /**
+ * Computes the renumbering vector needed
+ * by the @p{subdomain_wise}
+ * function. Does not perform the
+ * renumbering on the @p{DoFHandler} dofs
+ * but returns the renumbering vector.
+ */
+ template <int dim>
+ static void
+ compute_subdomain_wise (std::vector<unsigned int> &new_dof_indices,
+ const DoFHandler<dim> &dof_handler);
+
/**
* Exception
*/
* freedom are not mutually
* exclusive for different
* subdomain ids.
+ *
+ * If you want to get a unique
+ * association of degree of freedom with
+ * subdomains, use the
+ * @p{get_subdomain_association}
+ * function.
*/
template <int dim>
static void
extract_subdomain_dofs (const DoFHandler<dim> &dof_handler,
const unsigned int subdomain_id,
std::vector<bool> &selected_dofs);
+
+ /**
+ * For each DoF, return in the output
+ * array to which subdomain (as given by
+ * the @p{cell->subdomain_id()} function)
+ * it belongs. The output array is
+ * supposed to have the right size
+ * already when calling this function.
+ *
+ * Note that degrees of freedom
+ * associated with faces, edges, and
+ * vertices may be associated with
+ * multiple subdomains if they are
+ * sitting on partition boundaries. In
+ * these cases, we put them into one of
+ * the associated partitions in an
+ * undefined way. This may sometimes lead
+ * to different numbers of degrees of
+ * freedom in partitions, even if the
+ * number of cells is perfectly
+ * equidistributed. While this is
+ * regrettable, it is not a problem in
+ * practice since the number of degrees
+ * of freedom on partition boundaries is
+ * asymptotically vanishing as we refine
+ * the mesh as long as the number of
+ * partitions is kept constant.
+ */
+ template <int dim>
+ static void
+ get_subdomain_association (const DoFHandler<dim> &dof_handler,
+ std::vector<unsigned int> &subdomain);
/**
* Count how many degrees of
// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors
+// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
template <int dim>
void
-DoFRenumbering::component_wise (
- DoFHandler<dim> &dof_handler,
- const std::vector<unsigned int> &component_order_arg)
+DoFRenumbering::
+component_wise (DoFHandler<dim> &dof_handler,
+ const std::vector<unsigned int> &component_order_arg)
{
std::vector<unsigned int> renumbering (dof_handler.n_dofs(),
DoFHandler<dim>::invalid_dof_index);
template <int dim>
-void DoFRenumbering::component_wise (
- MGDoFHandler<dim>& dof_handler,
- unsigned int level,
- const std::vector<unsigned int> &component_order_arg)
+void DoFRenumbering::
+component_wise (MGDoFHandler<dim> &dof_handler,
+ const unsigned int level,
+ const std::vector<unsigned int> &component_order_arg)
{
std::vector<unsigned int> renumbering (dof_handler.n_dofs(level),
DoFHandler<dim>::invalid_dof_index);
template <int dim, class ITERATOR, class ENDITERATOR>
unsigned int
-DoFRenumbering::compute_component_wise (
- std::vector<unsigned int>& new_indices,
- ITERATOR& start,
- const ENDITERATOR& end,
- const std::vector<unsigned int> &component_order_arg)
+DoFRenumbering::
+compute_component_wise (std::vector<unsigned int>& new_indices,
+ ITERATOR& start,
+ const ENDITERATOR& end,
+ const std::vector<unsigned int> &component_order_arg)
{
// Modify for hp
const FiniteElement<dim>& fe = start->get_fe();
template <int dim>
void
-DoFRenumbering::compute_random (
- std::vector<unsigned int>& new_indices,
- const DoFHandler<dim> &dof_handler)
+DoFRenumbering::compute_random (std::vector<unsigned int> &new_indices,
+ const DoFHandler<dim> &dof_handler)
{
const unsigned int n_dofs = dof_handler.n_dofs();
Assert(new_indices.size() == n_dofs,
+template <int dim>
+void
+DoFRenumbering::subdomain_wise (DoFHandler<dim> &dof_handler)
+{
+ std::vector<unsigned int> renumbering(dof_handler.n_dofs(),
+ DoFHandler<dim>::invalid_dof_index);
+ compute_subdomain_wise(renumbering, dof_handler);
+
+ dof_handler.renumber_dofs(renumbering);
+}
+
+
+template <int dim>
+void
+DoFRenumbering::
+compute_subdomain_wise (std::vector<unsigned int> &new_dof_indices,
+ const DoFHandler<dim> &dof_handler)
+{
+ const unsigned int n_dofs = dof_handler.n_dofs();
+ Assert (new_dof_indices.size() == n_dofs,
+ ExcDimensionMismatch (new_dof_indices.size(), n_dofs));
+
+ // first get the association of each dof
+ // with a subdomain and determine the total
+ // number of subdomain ids used
+ std::vector<unsigned int> subdomain_association (n_dofs);
+ DoFTools::get_subdomain_association (dof_handler,
+ subdomain_association);
+ 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
+ // that the algorithm is stable, i.e. if
+ // two dofs i,j have i<j and belong to the
+ // same subdomain, then they will be in
+ // this order also after reordering
+ std::fill (new_dof_indices.begin(), new_dof_indices.end(),
+ deal_II_numbers::invalid_unsigned_int);
+ unsigned int next_free_index = 0;
+ for (unsigned int subdomain=0; subdomain<n_subdomains; ++subdomain)
+ for (unsigned int i=0; i<n_dofs; ++i)
+ if (subdomain_association[i] == subdomain)
+ {
+ Assert (new_dof_indices[i] == deal_II_numbers::invalid_unsigned_int,
+ ExcInternalError());
+ new_dof_indices[i] = next_free_index;
+ ++next_free_index;
+ }
+
+ // we should have numbered all dofs
+ Assert (next_free_index = n_dofs, ExcInternalError());
+ Assert (std::find (new_dof_indices.begin(), new_dof_indices.end(),
+ deal_II_numbers::invalid_unsigned_int)
+ == new_dof_indices.end(),
+ ExcInternalError());
+}
+
+
+
// explicit instantiations
template
void DoFRenumbering::random<deal_II_dimension>
(DoFHandler<deal_II_dimension> &);
+template
+void DoFRenumbering::subdomain_wise<deal_II_dimension>
+(DoFHandler<deal_II_dimension> &);
+
template
void
DoFRenumbering::compute_random<deal_II_dimension>
(std::vector<unsigned int>&,
const DoFHandler<deal_II_dimension> &);
+template
+void
+DoFRenumbering::compute_subdomain_wise<deal_II_dimension>
+(std::vector<unsigned int>&,
+ const DoFHandler<deal_II_dimension> &);
+
template
void DoFRenumbering::Cuthill_McKee<deal_II_dimension>
(MGDoFHandler<deal_II_dimension>&,
// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors
+// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
{
Assert(selected_dofs.size() == dof_handler.n_dofs(),
ExcDimensionMismatch(selected_dofs.size(), dof_handler.n_dofs()));
- // preset all values by false
+
+ // preset all values by false
std::fill_n (selected_dofs.begin(), dof_handler.n_dofs(), false);
const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell;
std::vector<unsigned int> local_dof_indices (dofs_per_cell);
- // this function is similar to the
- // make_sparsity_pattern function,
- // see there for more information
-
typename DoFHandler<dim>::active_cell_iterator
cell = dof_handler.begin_active(),
endc = dof_handler.end();
+template <int dim>
+void
+DoFTools::
+get_subdomain_association (const DoFHandler<dim> &dof_handler,
+ std::vector<unsigned int> &subdomain_association)
+{
+ Assert(subdomain_association.size() == dof_handler.n_dofs(),
+ ExcDimensionMismatch(subdomain_association.size(),
+ dof_handler.n_dofs()));
+
+ // preset all values by an invalid value
+ std::fill_n (subdomain_association.begin(), dof_handler.n_dofs(),
+ deal_II_numbers::invalid_unsigned_int);
+
+ const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell;
+ std::vector<unsigned int> local_dof_indices (dofs_per_cell);
+
+ // this function is similar to the
+ // make_sparsity_pattern function,
+ // see there for more information
+
+ typename DoFHandler<dim>::active_cell_iterator
+ cell = dof_handler.begin_active(),
+ endc = dof_handler.end();
+ for (; cell!=endc; ++cell)
+ {
+ const unsigned int subdomain_id = cell->subdomain_id();
+ cell->get_dof_indices (local_dof_indices);
+
+ // set subdomain ids. if dofs already
+ // have their values set then they must
+ // be on partition interfaces. don't
+ // worry about that, just overwrite it
+ for (unsigned int i=0; i<dofs_per_cell; ++i)
+ subdomain_association[local_dof_indices[i]] = subdomain_id;
+ };
+
+ Assert (std::find (subdomain_association.begin(),
+ subdomain_association.end(),
+ deal_II_numbers::invalid_unsigned_int)
+ == subdomain_association.end(),
+ ExcInternalError());
+}
+
+
+
template <int dim>
void
DoFTools::
const unsigned int subdomain_id,
std::vector<bool> &selected_dofs);
+template
+void
+DoFTools::get_subdomain_association<deal_II_dimension>
+(const DoFHandler<deal_II_dimension> &dof_handler,
+ std::vector<unsigned int> &subdomain_association);
+
template
void