From: Jean-Paul Pelteret Date: Sat, 27 Jan 2018 09:37:20 +0000 (+0100) Subject: Ensure parallel::Tria can't be used with DoFRenumbering::subdomain_wise X-Git-Tag: v9.0.0-rc1~508^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1630706bd391e04ab1e6715e73f9a06878c1bc3c;p=dealii.git Ensure parallel::Tria can't be used with DoFRenumbering::subdomain_wise The documentation is also updated to reflect this point. Fixes 5763. --- diff --git a/include/deal.II/dofs/dof_renumbering.h b/include/deal.II/dofs/dof_renumbering.h index 706ba23995..9abbd21f02 100644 --- a/include/deal.II/dofs/dof_renumbering.h +++ b/include/deal.II/dofs/dof_renumbering.h @@ -1138,9 +1138,13 @@ namespace DoFRenumbering * 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 - * GridTools::partition_triangulation function for this). + * one, etc. This is useful when doing parallel computations with a standard + * Triangulation after assigning subdomain ids using a partitioner (see the + * GridTools::partition_triangulation function for this). Calling this + * function is unnecessary when using a parallel::shared::Triangulation or + * parallel::distributed::Triangulation, as the degrees of freedom are + * already enumerated according to the MPI process id. Therefore, if the + * underlying triangulation is of this type then an error will be thrown. * * Note that degrees of freedom associated with faces, edges, and vertices * may be associated with multiple subdomains if they are sitting on diff --git a/source/dofs/dof_renumbering.cc b/source/dofs/dof_renumbering.cc index 585739eb10..68bffe826d 100644 --- a/source/dofs/dof_renumbering.cc +++ b/source/dofs/dof_renumbering.cc @@ -2031,6 +2031,9 @@ namespace DoFRenumbering void subdomain_wise (DoFHandlerType &dof_handler) { + Assert((!dynamic_cast*> (&dof_handler.get_triangulation())), + ExcMessage("Parallel triangulations are already enumerated according to their MPI process id.")); + std::vector renumbering(dof_handler.n_dofs(), numbers::invalid_dof_index); compute_subdomain_wise(renumbering, dof_handler);