From 085080d58ac536ae0e80f854dbf8ff92c9b9ad7c Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Thu, 19 Dec 2019 17:23:20 -0700 Subject: [PATCH] Deprecate some versions of DoFTools::extract_dofs(). Add other versions that *return* their results, rather than doing so via output arguments. --- include/deal.II/dofs/dof_tools.h | 85 ++++++++++++++++++++++++++++++-- source/dofs/dof_tools.cc | 79 +++++++++++++++++++++++++++-- source/dofs/dof_tools.inst.in | 26 ++++++++++ 3 files changed, 184 insertions(+), 6 deletions(-) diff --git a/include/deal.II/dofs/dof_tools.h b/include/deal.II/dofs/dof_tools.h index 2fe5012db3..2ac374a080 100644 --- a/include/deal.II/dofs/dof_tools.h +++ b/include/deal.II/dofs/dof_tools.h @@ -1323,9 +1323,12 @@ namespace DoFTools * for large parallel computations -- in fact, it may be too large to store on * each processor. In that case, you may want to choose the variant of this * function that returns an IndexSet object. + * + * @deprecated For the reason stated above, this function is deprecated in + * favor of the following function. */ template - void + DEAL_II_DEPRECATED void extract_hanging_node_dofs(const DoFHandler &dof_handler, std::vector & selected_dofs); @@ -1372,13 +1375,56 @@ namespace DoFTools * vector just holds the locally owned extracted degrees of freedom, which * first have to be mapped to the global degrees of freedom, to correspond * with them. + * + * @deprecated This function is difficult to use in parallel contexts because + * it returns a vector that needs to be indexed based on the position + * of a degree of freedom within the set of locally owned degrees of + * freedom. It is consequently deprecated in favor of the following + * function. */ template - void + DEAL_II_DEPRECATED void extract_dofs(const DoFHandlerType &dof_handler, const ComponentMask & component_mask, std::vector & selected_dofs); + /** + * Extract the (locally owned) indices of the degrees of freedom belonging to + * certain vector components of a vector-valued finite element. The + * @p component_mask defines which components or blocks of an FESystem or + * vector-valued element are to be extracted + * from the DoFHandler @p dof. The entries in the output object then + * correspond to degrees of freedom belonging to these + * components. + * + * If the finite element under consideration is not primitive, i.e., some or + * all of its shape functions are non-zero in more than one vector component + * (which holds, for example, for FE_Nedelec or FE_RaviartThomas elements), + * then shape functions cannot be associated with a single vector component. + * In this case, if one shape vector component of this element is + * flagged in @p component_mask (see + * @ref GlossComponentMask), + * then this is equivalent to selecting all vector components + * corresponding to this non-primitive base element. + * + * @param[in] dof_handler The DoFHandler whose enumerated degrees of freedom + * are to be filtered by this function. + * @param[in] component_mask A mask that states which components you want + * to select. The size of this mask must be compatible with the number of + * components in the FiniteElement used by the @p dof_handler. See + * @ref GlossComponentMask "the glossary entry on component masks" + * for more information. + * @return An IndexSet object that will contain exactly those entries that + * (i) correspond to degrees of freedom selected by the mask above, and + * (ii) are locally owned. The size of the index set is equal to the global + * number of degrees of freedom. Note that the resulting object is always + * a subset of what DoFHandler::locally_owned_dofs() returns. + */ + template + IndexSet + extract_dofs(const DoFHandlerType &dof_handler, + const ComponentMask & component_mask); + /** * This function is the equivalent to the DoFTools::extract_dofs() functions * above except that the selection of which degrees of freedom to extract is @@ -1402,13 +1448,46 @@ namespace DoFTools * of this array must equal DoFHandler::n_locally_owned_dofs(), which for * sequential computations of course equals DoFHandler::n_dofs(). The * previous contents of this array are overwritten. + * + * @deprecated This function is difficult to use in parallel contexts because + * it returns a vector that needs to be indexed based on the position + * of a degree of freedom within the set of locally owned degrees of + * freedom. It is consequently deprecated in favor of the following + * function. */ template - void + DEAL_II_DEPRECATED void extract_dofs(const DoFHandlerType &dof_handler, const BlockMask & block_mask, std::vector & selected_dofs); + /** + * This function is the equivalent to the DoFTools::extract_dofs() functions + * above except that the selection of which degrees of freedom to extract is + * not done based on components (see + * @ref GlossComponent) + * but instead based on whether they are part of a particular block (see + * @ref GlossBlock). + * Consequently, the second argument is not a ComponentMask but a BlockMask + * object. + * + * @param[in] dof_handler The DoFHandler whose enumerated degrees of freedom + * are to be filtered by this function. + * @param[in] block_mask A mask that states which blocks you want + * to select. The size of this mask must be compatible with the number of + * blocks in the FiniteElement used by the @p dof_handler. See + * @ref GlossBlockMask "the glossary entry on block masks" + * for more information. + * @return An IndexSet object that will contain exactly those entries that + * (i) correspond to degrees of freedom selected by the mask above, and + * (ii) are locally owned. The size of the index set is equal to the global + * number of degrees of freedom. Note that the resulting object is always + * a subset of what DoFHandler::locally_owned_dofs() returns. + */ + template + IndexSet + extract_dofs(const DoFHandlerType &dof_handler, const BlockMask &block_mask); + /** * Do the same thing as the corresponding extract_dofs() function for one * level of a multi-grid DoF numbering. diff --git a/source/dofs/dof_tools.cc b/source/dofs/dof_tools.cc index d785bba776..da3cb6402e 100644 --- a/source/dofs/dof_tools.cc +++ b/source/dofs/dof_tools.cc @@ -227,9 +227,7 @@ namespace DoFTools // then loop over all cells and do the work std::vector indices; - for (typename DoFHandlerType::active_cell_iterator c = dof.begin_active(); - c != dof.end(); - ++c) + for (const auto &c : dof.active_cell_iterators()) if (c->is_locally_owned()) { const unsigned int fe_index = c->active_fe_index(); @@ -443,6 +441,70 @@ namespace DoFTools + template + IndexSet + extract_dofs(const DoFHandlerType &dof, const ComponentMask &component_mask) + { + Assert(component_mask.represents_n_components( + dof.get_fe_collection().n_components()), + ExcMessage( + "The given component mask is not sized correctly to represent the " + "components of the given finite element.")); + + // Two special cases: no component is selected, and all components are + // selected; both rather stupid, but easy to catch + if (component_mask.n_selected_components( + dof.get_fe_collection().n_components()) == 0) + return IndexSet(dof.n_dofs()); + else if (component_mask.n_selected_components( + dof.get_fe_collection().n_components()) == + dof.get_fe_collection().n_components()) + return dof.locally_owned_dofs(); + + // If we have to deal with a partial component mask, then + // set up a table for the degrees of freedom on the reference cell (for + // each element in the collection) to see whether it belongs to one of the + // desired components. + const auto & fe_collection = dof.get_fe_collection(); + std::vector> local_selected_dofs; + Assert(fe_collection.n_components() < 256, ExcNotImplemented()); + for (unsigned int f = 0; f < fe_collection.size(); ++f) + { + // First get the component for each shape function + const std::vector local_component_association = + internal::get_local_component_association(fe_collection[f], + component_mask); + + // Check which dofs were selected + std::vector this_selected_dofs(fe_collection[f].dofs_per_cell); + for (unsigned int i = 0; i < fe_collection[f].dofs_per_cell; ++i) + this_selected_dofs[i] = + component_mask[local_component_association[i]]; + + local_selected_dofs.emplace_back(std::move(this_selected_dofs)); + } + + // Then loop over all cells and do the work + IndexSet selected_dofs(dof.n_dofs()); + std::vector local_dof_indices; + for (const auto &c : dof.active_cell_iterators()) + if (c->is_locally_owned()) + { + const unsigned int fe_index = c->active_fe_index(); + const unsigned int dofs_per_cell = c->get_fe().dofs_per_cell; + local_dof_indices.resize(dofs_per_cell); + c->get_dof_indices(local_dof_indices); + for (unsigned int i = 0; i < dofs_per_cell; ++i) + if (dof.locally_owned_dofs().is_element(local_dof_indices[i]) && + local_selected_dofs[fe_index][i]) + selected_dofs.add_index(local_dof_indices[i]); + } + + return selected_dofs; + } + + + template void extract_dofs(const DoFHandlerType &dof, @@ -456,6 +518,17 @@ namespace DoFTools + template + IndexSet + extract_dofs(const DoFHandlerType &dof, const BlockMask &block_mask) + { + // simply forward to the function that works based on a component mask + return extract_dofs( + dof, dof.get_fe_collection().component_mask(block_mask)); + } + + + template std::vector locally_owned_dofs_per_component(const DoFHandlerType &dof, diff --git a/source/dofs/dof_tools.inst.in b/source/dofs/dof_tools.inst.in index c3d0d14f60..569ade4fe0 100644 --- a/source/dofs/dof_tools.inst.in +++ b/source/dofs/dof_tools.inst.in @@ -704,6 +704,8 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : DIMENSIONS) deal_II_space_dimension>::active_cell_iterator> &patch); + // Deprecated versions that return information through the last + // argument template void extract_dofs>( const DoFHandler &, @@ -728,6 +730,30 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : DIMENSIONS) const BlockMask &, std::vector &); + + // New versions that return information as an IndexSet + template IndexSet + extract_dofs>( + const DoFHandler &, + const ComponentMask &); + + template IndexSet + extract_dofs>( + const DoFHandler &, + const BlockMask &); + + template IndexSet + extract_dofs>( + const hp::DoFHandler &, + const ComponentMask &); + + template IndexSet + extract_dofs>( + const hp::DoFHandler &, + const BlockMask &); + + + template void make_cell_patches( SparsityPattern &, -- 2.39.5