]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement a second DoFTools::extract_boundary_dofs function that also works for paral...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 26 Apr 2012 08:38:34 +0000 (08:38 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 26 Apr 2012 08:38:34 +0000 (08:38 +0000)
git-svn-id: https://svn.dealii.org/trunk@25449 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/dofs/dof_tools.h
deal.II/source/dofs/dof_tools.cc

index 7f0b65932f80007cb0fc45ef53972291798230ed..fbbccfd87ebfdabc4555507376414b40fa1bb624 100644 (file)
@@ -175,6 +175,12 @@ enabled due to a missing include file in file
 <h3>Specific improvements</h3>
 
 <ol>
+<li> New: There is now a second DoFTools::extract_boundary_dofs
+function that also works for parallel::distributed::Triangulation
+triangulations.
+<br>
+(Wolfgang Bangerth, 2012/04/26)
+
 <li> New: The FullMatrix::extract_submatrix_from, FullMatrix::scatter_matrix_to,
 FullMatrix::set functions are new.
 <br>
index 633d040ec494f1285ec8d3c0e1d8ba846c16a78f..27014eeb90648f38a7ae757eb5ffcef87a8f4c32 100644 (file)
@@ -1069,68 +1069,84 @@ namespace DoFTools
                       std::vector<bool>       &selected_dofs,
                       const bool               blocks = false);
 
-                                   /**
-                                    * Extract all degrees of freedom
-                                    * which are at the boundary and
-                                    * belong to specified components
-                                    * of the solution. The function
-                                    * returns its results in the
-                                    * last non-default-valued
-                                    * parameter which contains
-                                    * @p true if a degree of
-                                    * freedom is at the boundary and
-                                    * belongs to one of the selected
-                                    * components, and @p false
-                                    * otherwise.
-                                    *
-                                    * By specifying the
-                                    * @p boundary_indicator
-                                    * variable, you can select which
-                                    * boundary indicators the faces
-                                    * have to have on which the
-                                    * degrees of freedom are located
-                                    * that shall be extracted. If it
-                                    * is an empty list, then all
-                                    * boundary indicators are
-                                    * accepted.
-                                    *
-                                    * The size of @p component_select
-                                    * shall equal the number of
-                                    * components in the finite
-                                    * element used by @p dof. The
-                                    * size of @p selected_dofs shall
-                                    * equal
-                                    * <tt>dof_handler.n_dofs()</tt>. Previous
-                                    * contents of this array or
-                                    * overwritten.
-                                    *
-                                    * Using the usual convention, if
-                                    * a shape function is non-zero
-                                    * in more than one component
-                                    * (i.e. it is non-primitive),
-                                    * then the element in the
-                                    * component mask is used that
-                                    * corresponds to the first
-                                    * non-zero components. Elements
-                                    * in the mask corresponding to
-                                    * later components are ignored.
-                                    *
-                                    * @note This function will not work
-                                    * for DoFHandler objects that are
-                                    * built on a
-                                    * parallel::distributed::Triangulation
-                                    * object. The reasons is that the
-                                    * output argument @p selected_dofs
-                                    * has to have a length equal to <i>all</i>
-                                    * global degrees of freedom.
-                                    * Consequently, this does not scale
-                                    * to very large problems. If you
-                                    * need the functionality of this
-                                    * function for parallel triangulations,
-                                    * then you need to use the other
-                                    * DoFTools::extract_boundary_dofs
-                                    * function.
-                                    */
+  /**
+   * Extract all degrees of freedom
+   * which are at the boundary and
+   * belong to specified components
+   * of the solution. The function
+   * returns its results in the
+   * last non-default-valued
+   * parameter which contains
+   * @p true if a degree of
+   * freedom is at the boundary and
+   * belongs to one of the selected
+   * components, and @p false
+   * otherwise.
+   *
+   * By specifying the
+   * @p boundary_indicator
+   * variable, you can select which
+   * boundary indicators the faces
+   * have to have on which the
+   * degrees of freedom are located
+   * that shall be extracted. If it
+   * is an empty list, then all
+   * boundary indicators are
+   * accepted.
+   *
+   * The size of @p component_select
+   * shall equal the number of
+   * components in the finite
+   * element used by @p dof. The
+   * size of @p selected_dofs shall
+   * equal
+   * <tt>dof_handler.n_dofs()</tt>. Previous
+   * contents of this array or
+   * overwritten.
+   *
+   * Using the usual convention, if
+   * a shape function is non-zero
+   * in more than one component
+   * (i.e. it is non-primitive),
+   * then the element in the
+   * component mask is used that
+   * corresponds to the first
+   * non-zero components. Elements
+   * in the mask corresponding to
+   * later components are ignored.
+   *
+   * @note This function will not work
+   * for DoFHandler objects that are
+   * built on a
+   * parallel::distributed::Triangulation
+   * object. The reasons is that the
+   * output argument @p selected_dofs
+   * has to have a length equal to <i>all</i>
+   * global degrees of freedom.
+   * Consequently, this does not scale
+   * to very large problems. If you
+   * need the functionality of this
+   * function for parallel triangulations,
+   * then you need to use the other
+   * DoFTools::extract_boundary_dofs
+   * function.
+   *
+   * @param dof_handler The object that describes which degrees of freedom
+   *          live on which cell
+   * @param component_select A mask denoting the vector components of the
+   *          finite element that should be considered.
+   * @param selected_dofs The IndexSet object that is returned and that
+   *          will contain the indices of degrees of freedom that are
+   *          located on the boundary (and correspond to the selected
+   *          vector components and boundary indicators, depending on
+   *          the values of the @p component_select and @p boundary_indicators
+   *          arguments).
+   * @param boundary_indicators If empty, this function extracts the
+   *          indices of the degrees of freedom for all parts of the boundary.
+   *          If it is a non-empty list, then the function only considers
+   *          boundary faces with the boundary indicators listed in this
+   *          argument.
+   */
   template <class DH>
   void
   extract_boundary_dofs (const DH                   &dof_handler,
@@ -1138,6 +1154,41 @@ namespace DoFTools
                          std::vector<bool>          &selected_dofs,
                          const std::set<types::boundary_id_t> &boundary_indicators = std::set<types::boundary_id_t>());
 
+  /**
+   * This function does the same as the previous one but it
+   * returns its result as an IndexSet rather than a std::vector@<bool@>.
+   * Thus, it can also be called for DoFHandler objects that are
+   * defined on parallel::distributed::Triangulation objects.
+   *
+   * @note If the DoFHandler object is indeed defined on a
+   * parallel::distributed::Triangulation, then the @p selected_dofs
+   * index set will contain only those degrees of freedom on the
+   * boundary that belong to the locally relevant set (see
+   * @ref GlossLocallyRelevantDof "locally relevant DoFs").
+   *
+   * @param dof_handler The object that describes which degrees of freedom
+   *          live on which cell
+   * @param component_select A mask denoting the vector components of the
+   *          finite element that should be considered.
+   * @param selected_dofs The IndexSet object that is returned and that
+   *          will contain the indices of degrees of freedom that are
+   *          located on the boundary (and correspond to the selected
+   *          vector components and boundary indicators, depending on
+   *          the values of the @p component_select and @p boundary_indicators
+   *          arguments).
+   * @param boundary_indicators If empty, this function extracts the
+   *          indices of the degrees of freedom for all parts of the boundary.
+   *          If it is a non-empty list, then the function only considers
+   *          boundary faces with the boundary indicators listed in this
+   *          argument.
+   */
+  template <class DH>
+  void
+  extract_boundary_dofs (const DH                   &dof_handler,
+                         const std::vector<bool>    &component_select,
+                         IndexSet                    &selected_dofs,
+                         const std::set<types::boundary_id_t> &boundary_indicators = std::set<types::boundary_id_t>());
+
                                    /**
                                     * This function is similar to
                                     * the extract_boundary_dofs()
index f865b18bf24b1f03233f6487cd620172d41b0018..b5eb53fe81ed0bc46808492bfced7cc88fa31dd4 100644 (file)
@@ -3752,9 +3752,6 @@ namespace DoFTools
                          std::vector<bool>             &selected_dofs,
                          const std::set<types::boundary_id_t> &boundary_indicators)
   {
-    AssertDimension (component_select.size(), n_components(dof_handler));
-    Assert (boundary_indicators.find (types::internal_face_boundary_id) == boundary_indicators.end(),
-            ExcInvalidBoundaryIndicator());
     Assert ((dynamic_cast<const parallel::distributed::Triangulation<DH::dimension,DH::space_dimension>*>
              (&dof_handler.get_tria())
              ==
@@ -3762,8 +3759,35 @@ namespace DoFTools
             ExcMessage ("This function can not be used with distributed triangulations."
                         "See the documentation for more information."));
 
+    IndexSet indices;
+    extract_boundary_dofs (dof_handler, component_select,
+                           indices, boundary_indicators);
+
+    // clear and reset array by default values
+    selected_dofs.clear ();
+    selected_dofs.resize (dof_handler.n_dofs(), false);
+
+    // then convert the values computed above to the binary vector
+    indices.fill_binary_vector(selected_dofs);
+  }
+
+
+  template <class DH>
+  void
+  extract_boundary_dofs (const DH                      &dof_handler,
+                         const std::vector<bool>       &component_select,
+                         IndexSet             &selected_dofs,
+                         const std::set<types::boundary_id_t> &boundary_indicators)
+  {
+    AssertDimension (component_select.size(), n_components(dof_handler));
+    Assert (boundary_indicators.find (types::internal_face_boundary_id) == boundary_indicators.end(),
+            ExcInvalidBoundaryIndicator());
     const unsigned int dim=DH::dimension;
 
+    // first reset output argument
+    selected_dofs.clear ();
+    selected_dofs.set_size(dof_handler.n_dofs());
+
                                      // let's see whether we have to
                                      // check for certain boundary
                                      // indicators or whether we can
@@ -3777,10 +3801,6 @@ namespace DoFTools
       = (component_select != std::vector<bool>(component_select.size(),
                                                true));
 
-                                     // clear and reset array by default
-                                     // values
-    selected_dofs.clear ();
-    selected_dofs.resize (dof_handler.n_dofs(), false);
     std::vector<unsigned int> face_dof_indices;
     face_dof_indices.reserve (max_dofs_per_face(dof_handler));
 
@@ -3799,79 +3819,83 @@ namespace DoFTools
                                      // sooner or later
     for (typename DH::active_cell_iterator cell=dof_handler.begin_active();
          cell!=dof_handler.end(); ++cell)
-      for (unsigned int face=0;
-           face<GeometryInfo<DH::dimension>::faces_per_cell; ++face)
-        if (cell->at_boundary(face))
-          if (! check_boundary_indicator ||
-              (boundary_indicators.find (cell->face(face)->boundary_indicator())
-               != boundary_indicators.end()))
-            {
-              const FiniteElement<DH::dimension, DH::space_dimension> &fe = cell->get_fe();
+      // only work on cells that are either locally owned or at
+      // least ghost cells
+      if (cell->is_artificial() == false)
+        for (unsigned int face=0;
+            face<GeometryInfo<DH::dimension>::faces_per_cell; ++face)
+          if (cell->at_boundary(face))
+            if (! check_boundary_indicator ||
+                (boundary_indicators.find (cell->face(face)->boundary_indicator())
+                    != boundary_indicators.end()))
+              {
+                const FiniteElement<DH::dimension, DH::space_dimension> &fe = cell->get_fe();
 
-              const unsigned int dofs_per_face = fe.dofs_per_face;
-              face_dof_indices.resize (dofs_per_face);
-              cell->face(face)->get_dof_indices (face_dof_indices,
-                                                 cell->active_fe_index());
+                const unsigned int dofs_per_face = fe.dofs_per_face;
+                face_dof_indices.resize (dofs_per_face);
+                cell->face(face)->get_dof_indices (face_dof_indices,
+                                                   cell->active_fe_index());
 
-              for (unsigned int i=0; i<fe.dofs_per_face; ++i)
-                if (!check_vector_component)
-                  selected_dofs[face_dof_indices[i]] = true;
-                else
-                                                   // check for
-                                                   // component is
-                                                   // required. somewhat
-                                                   // tricky as usual
-                                                   // for the case that
-                                                   // the shape function
-                                                   // is non-primitive,
-                                                   // but use usual
-                                                   // convention (see
-                                                   // docs)
-                  {
-                                                     // first get at the
-                                                     // cell-global
-                                                     // number of a face
-                                                     // dof, to ask the
-                                                     // fe certain
-                                                     // questions
-                    const unsigned int cell_index
+                for (unsigned int i=0; i<fe.dofs_per_face; ++i)
+                  if (!check_vector_component)
+                    selected_dofs.add_index (face_dof_indices[i]);
+                  else
+                    // check for
+                    // component is
+                    // required. somewhat
+                    // tricky as usual
+                    // for the case that
+                    // the shape function
+                    // is non-primitive,
+                    // but use usual
+                    // convention (see
+                    // docs)
+                    {
+                      // first get at the
+                      // cell-global
+                      // number of a face
+                      // dof, to ask the
+                      // fe certain
+                      // questions
+                      const unsigned int cell_index
                       = (dim == 1 ?
-                         i
-                         :
-                         (dim == 2 ?
-                          (i<2*fe.dofs_per_vertex ? i : i+2*fe.dofs_per_vertex)
+                          i
                           :
-                          (dim == 3 ?
-                           (i<4*fe.dofs_per_vertex ?
-                            i
-                            :
-                            (i<4*fe.dofs_per_vertex+4*fe.dofs_per_line ?
-                             i+4*fe.dofs_per_vertex
-                             :
-                             i+4*fe.dofs_per_vertex+8*fe.dofs_per_line))
-                           :
-                           numbers::invalid_unsigned_int)));
-                    if (fe.is_primitive (cell_index))
-                      selected_dofs[face_dof_indices[i]]
-                        = (component_select[fe.face_system_to_component_index(i).first]
-                           == true);
-                    else // not primitive
-                      {
-                        const unsigned int first_nonzero_comp
+                          (dim == 2 ?
+                              (i<2*fe.dofs_per_vertex ? i : i+2*fe.dofs_per_vertex)
+                              :
+                              (dim == 3 ?
+                                  (i<4*fe.dofs_per_vertex ?
+                                      i
+                                      :
+                                      (i<4*fe.dofs_per_vertex+4*fe.dofs_per_line ?
+                                          i+4*fe.dofs_per_vertex
+                                          :
+                                          i+4*fe.dofs_per_vertex+8*fe.dofs_per_line))
+                                          :
+                                          numbers::invalid_unsigned_int)));
+                      if (fe.is_primitive (cell_index))
+                        {
+                          if (component_select[fe.face_system_to_component_index(i).first]
+                                               == true)
+                            selected_dofs.add_index (face_dof_indices[i]);
+                        }
+                      else // not primitive
+                        {
+                          const unsigned int first_nonzero_comp
                           = (std::find (fe.get_nonzero_components(cell_index).begin(),
                                         fe.get_nonzero_components(cell_index).end(),
                                         true)
-                             -
-                             fe.get_nonzero_components(cell_index).begin());
-                        Assert (first_nonzero_comp < fe.n_components(),
-                                ExcInternalError());
-
-                        selected_dofs[face_dof_indices[i]]
-                          = (component_select[first_nonzero_comp]
-                             == true);
-                      }
-                  }
-            }
+                          -
+                          fe.get_nonzero_components(cell_index).begin());
+                          Assert (first_nonzero_comp < fe.n_components(),
+                                  ExcInternalError());
+
+                          if (component_select[first_nonzero_comp] == true)
+                            selected_dofs.add_index (face_dof_indices[i]);
+                        }
+                    }
+              }
   }
 
 

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.