]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Deprecate some versions of DoFTools::extract_dofs().
authorWolfgang Bangerth <bangerth@colostate.edu>
Fri, 20 Dec 2019 00:23:20 +0000 (17:23 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Sat, 21 Dec 2019 05:19:49 +0000 (22:19 -0700)
Add other versions that *return* their results, rather than doing so via
output arguments.

include/deal.II/dofs/dof_tools.h
source/dofs/dof_tools.cc
source/dofs/dof_tools.inst.in

index 2fe5012db38da3a12fa9066c78b6b674ac8d15c3..2ac374a080b50e12adf4aed0afebd851d9a95114 100644 (file)
@@ -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 <int dim, int spacedim>
-  void
+  DEAL_II_DEPRECATED void
   extract_hanging_node_dofs(const DoFHandler<dim, spacedim> &dof_handler,
                             std::vector<bool> &              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 <typename DoFHandlerType>
-  void
+  DEAL_II_DEPRECATED void
   extract_dofs(const DoFHandlerType &dof_handler,
                const ComponentMask & component_mask,
                std::vector<bool> &   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 <em>one</em> shape vector component of this element is
+   * flagged in @p component_mask (see
+   * @ref GlossComponentMask),
+   * then this is equivalent to selecting <em>all</em> 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 <typename DoFHandlerType>
+  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 <typename DoFHandlerType>
-  void
+  DEAL_II_DEPRECATED void
   extract_dofs(const DoFHandlerType &dof_handler,
                const BlockMask &     block_mask,
                std::vector<bool> &   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 <typename DoFHandlerType>
+  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.
index d785bba776e5e42f5d6cd9325ffedb382e7a7183..da3cb6402e02b58e677ce8cf4105826d0fd4c33d 100644 (file)
@@ -227,9 +227,7 @@ namespace DoFTools
 
       // then loop over all cells and do the work
       std::vector<types::global_dof_index> 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 <typename DoFHandlerType>
+  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<std::vector<bool>> 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<unsigned char> local_component_association =
+          internal::get_local_component_association(fe_collection[f],
+                                                    component_mask);
+
+        // Check which dofs were selected
+        std::vector<bool> 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<types::global_dof_index> 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 <typename DoFHandlerType>
   void
   extract_dofs(const DoFHandlerType &dof,
@@ -456,6 +518,17 @@ namespace DoFTools
 
 
 
+  template <typename DoFHandlerType>
+  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<DoFHandlerType>(
+      dof, dof.get_fe_collection().component_mask(block_mask));
+  }
+
+
+
   template <typename DoFHandlerType>
   std::vector<IndexSet>
   locally_owned_dofs_per_component(const DoFHandlerType &dof,
index c3d0d14f60c6531c47c19948bf85421d18ad7343..569ade4fe0eab25b0580a37498d37c7b1c8198fa 100644 (file)
@@ -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<DoFHandler<deal_II_dimension, deal_II_space_dimension>>(
         const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
@@ -728,6 +730,30 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : DIMENSIONS)
         const BlockMask &,
         std::vector<bool> &);
 
+
+      // New versions that return information as an IndexSet
+      template IndexSet
+      extract_dofs<DoFHandler<deal_II_dimension, deal_II_space_dimension>>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const ComponentMask &);
+
+      template IndexSet
+      extract_dofs<DoFHandler<deal_II_dimension, deal_II_space_dimension>>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const BlockMask &);
+
+      template IndexSet
+      extract_dofs<hp::DoFHandler<deal_II_dimension, deal_II_space_dimension>>(
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const ComponentMask &);
+
+      template IndexSet
+      extract_dofs<hp::DoFHandler<deal_II_dimension, deal_II_space_dimension>>(
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const BlockMask &);
+
+
+
       template void
       make_cell_patches<deal_II_dimension, deal_II_space_dimension>(
         SparsityPattern &,

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.