]> https://gitweb.dealii.org/ - dealii.git/commitdiff
implement comments 4103/head
authorDenis Davydov <davydden@gmail.com>
Sat, 15 Apr 2017 08:35:17 +0000 (10:35 +0200)
committerDenis Davydov <davydden@gmail.com>
Sat, 15 Apr 2017 08:35:17 +0000 (10:35 +0200)
16 files changed:
doc/doxygen/images/extract_dofs_with_support_contained_within.png [moved from doc/doxygen/images/extract_dofs_with_support_on.png with 100% similarity]
doc/news/changes/minor/20170324DenisDavydov
include/deal.II/dofs/dof_tools.h
source/dofs/dof_tools.cc
source/dofs/dof_tools.inst.in
tests/dofs/dof_tools_extract_dofs_with_support_contained_within_01.cc [moved from tests/dofs/dof_tools_extract_dofs_with_support_on_01.cc with 97% similarity]
tests/dofs/dof_tools_extract_dofs_with_support_contained_within_01.output [moved from tests/dofs/dof_tools_extract_dofs_with_support_on_01.output with 100% similarity]
tests/dofs/dof_tools_extract_dofs_with_support_contained_within_02.cc [moved from tests/dofs/dof_tools_extract_dofs_with_support_on_02.cc with 94% similarity]
tests/dofs/dof_tools_extract_dofs_with_support_contained_within_02.mpirun=1.output [moved from tests/dofs/dof_tools_extract_dofs_with_support_on_02.mpirun=1.output with 100% similarity]
tests/dofs/dof_tools_extract_dofs_with_support_contained_within_02.mpirun=3.output [moved from tests/dofs/dof_tools_extract_dofs_with_support_on_02.mpirun=3.output with 100% similarity]
tests/dofs/dof_tools_extract_dofs_with_support_contained_within_02.mpirun=5.output [moved from tests/dofs/dof_tools_extract_dofs_with_support_on_02.mpirun=5.output with 100% similarity]
tests/dofs/dof_tools_extract_dofs_with_support_contained_within_03.cc [moved from tests/dofs/dof_tools_extract_dofs_with_support_on_03.cc with 93% similarity]
tests/dofs/dof_tools_extract_dofs_with_support_contained_within_03.output [moved from tests/dofs/dof_tools_extract_dofs_with_support_on_03.output with 100% similarity]
tests/dofs/dof_tools_extract_dofs_with_support_contained_within_04.cc [moved from tests/dofs/dof_tools_extract_dofs_with_support_on_04.cc with 98% similarity]
tests/dofs/dof_tools_extract_dofs_with_support_contained_within_04.mpirun=2.output [moved from tests/dofs/dof_tools_extract_dofs_with_support_on_04.mpirun=2.output with 100% similarity]
tests/dofs/dof_tools_extract_dofs_with_support_contained_within_04.output [moved from tests/dofs/dof_tools_extract_dofs_with_support_on_04.output with 100% similarity]

index 462f1a659c6071b5a19b65d29187e4117a4eabe0..212711ae1d3cfb98a947ba6c7675c9ce9f6d4dff 100644 (file)
@@ -1,4 +1,5 @@
-New: DoFTools::extract_dofs_with_support_on() returns a set of degrees of
-freedom that have support on a subset of locally owned cells.
+New: DoFTools::extract_dofs_with_support_contained_within() returns a set of 
+degrees of freedom whose support is entirely contained within the cells for 
+which the predicate returns true.
 <br>
 (Denis Davydov, 2017/03/24)
index d0eead02416be909a1feb7178e920c0bdffd4cdd..2ae871e631ef86007dcf4e1a6d45c19daa098f05 100644 (file)
@@ -1437,18 +1437,18 @@ namespace DoFTools
                                          const std::set<types::boundary_id> &boundary_ids = std::set<types::boundary_id>());
 
   /**
-   * Extract all degrees of freedom (DoFs) that have support only within cells,
-   * for which @p predicate is <code>true</code>. The result is returned as
-   * an IndexSet.
+   * Extract all indices of shape functions such that their support is entirely
+   * contained within the cells for which the @p predicate is <code>true</code>.
+   * The result is returned as an IndexSet.
    *
    * Consider the following FE space where predicate returns <code>true</code>
-   * for all cells on the left side:
+   * for all cells on the left half of the domain:
    *
-   * @image html extract_dofs_with_support_on.png
+   * @image html extract_dofs_with_support_contained_within.png
    *
-   * This functions will return the union of all DoFs on those cells minus
+   * This functions will return the union of all DoF indices on those cells minus
    * DoF 11, 13, 2 and 0; the result will be <code>[9,10], 12, [14,38]</code>.
-   * On the image above the returned DoFs are separated from the rest by the
+   * In the image above the returned DoFs are separated from the rest by the
    * red line
    *
    * Essentially, the question this functions answers is the following:
@@ -1460,14 +1460,14 @@ namespace DoFTools
    *
    * In case of parallel::distributed::Triangulation @p predicate will be called
    * only for locally owned and ghost cells. The resulting index set may contain
-   * DoFs that are associated to the locally owned or ghost cells, but are not
+   * DoFs that are associated with the locally owned or ghost cells, but are not
    * owned by the current MPI core.
    */
   template <typename DoFHandlerType>
   IndexSet
-  extract_dofs_with_support_on(const DoFHandlerType &dof_handler,
-                               const std::function< bool(const typename DoFHandlerType::active_cell_iterator &)> &predicate,
-                               const ConstraintMatrix &cm = ConstraintMatrix());
+  extract_dofs_with_support_contained_within(const DoFHandlerType &dof_handler,
+                                             const std::function< bool(const typename DoFHandlerType::active_cell_iterator &)> &predicate,
+                                             const ConstraintMatrix &cm = ConstraintMatrix());
 
   /**
    * Extract a vector that represents the constant modes of the DoFHandler for
index 796335ef7fce87e4d25e6958965521d3a73812c9..899169dc8ac6ebc1c7b277c4fac7d71e9262921e 100644 (file)
@@ -767,9 +767,9 @@ namespace DoFTools
 
   template <typename DoFHandlerType>
   IndexSet
-  extract_dofs_with_support_on(const DoFHandlerType &dof_handler,
-                               const std::function< bool(const typename DoFHandlerType::active_cell_iterator &)> &predicate,
-                               const ConstraintMatrix &cm)
+  extract_dofs_with_support_contained_within(const DoFHandlerType &dof_handler,
+                                             const std::function< bool(const typename DoFHandlerType::active_cell_iterator &)> &predicate,
+                                             const ConstraintMatrix &cm)
   {
     const std::function< bool(const typename DoFHandlerType::active_cell_iterator &)> predicate_local
     = [=] (const typename DoFHandlerType::active_cell_iterator & cell) -> bool
@@ -796,7 +796,7 @@ namespace DoFTools
         }
 
     // Get halo layer and accumulate its DoFs
-    std::set<types::global_dof_index> halo_dofs;
+    std::set<types::global_dof_index> dofs_with_support_on_halo_cells;
 
     const std::vector<typename DoFHandlerType::active_cell_iterator> halo_cells =
       GridTools::compute_active_cell_halo_layer(dof_handler, predicate_local);
@@ -816,55 +816,51 @@ namespace DoFTools
         const unsigned int dofs_per_cell = (*it)->get_fe().dofs_per_cell;
         local_dof_indices.resize (dofs_per_cell);
         (*it)->get_dof_indices (local_dof_indices);
-        halo_dofs.insert(local_dof_indices.begin(),
-                         local_dof_indices.end());
+        dofs_with_support_on_halo_cells.insert(local_dof_indices.begin(),
+                                               local_dof_indices.end());
       }
 
     // A complication coming from the fact that DoFs living on locally
-    // owned cells which satisfy predicate may participate in constrains for
+    // owned cells which satisfy predicate may participate in constraints for
     // DoFs outside of this region.
     if (cm.n_constraints() > 0)
       {
-        const std::vector<std::pair<types::global_dof_index,double> > *line_ptr;
-
-        // Remove DoFs that are part of constrains for DoFs outside
+        // Remove DoFs that are part of constraints for DoFs outside
         // of predicate. Since we will subtract halo_dofs from predicate_dofs,
         // achieve this by extending halo_dofs with DoFs to which
         // halo_dofs are constrained.
         std::set<types::global_dof_index> extra_halo;
-        for (std::set<types::global_dof_index>::iterator it = halo_dofs.begin();
-             it != halo_dofs.end(); ++it)
-          {
-            line_ptr = cm.get_constraint_entries(*it);
-            // if halo DoF is constrained, add all DoFs to which it's constrained.
-            if (line_ptr!=NULL)
-              {
-                const unsigned int line_size = line_ptr->size();
-                for (unsigned int j=0; j<line_size; ++j)
-                  extra_halo.insert((*line_ptr)[j].first);
-              }
-          }
+        for (std::set<types::global_dof_index>::iterator it = dofs_with_support_on_halo_cells.begin();
+             it != dofs_with_support_on_halo_cells.end(); ++it)
+          // if halo DoF is constrained, add all DoFs to which it's constrained
+          // because after resolving constraints, the support of the DoFs that
+          // constrain the current DoF will extend to the halo cells.
+          if (const std::vector<std::pair<types::global_dof_index,double> > *line_ptr = cm.get_constraint_entries(*it))
+            {
+              const unsigned int line_size = line_ptr->size();
+              for (unsigned int j=0; j<line_size; ++j)
+                extra_halo.insert((*line_ptr)[j].first);
+            }
 
-        halo_dofs.insert(extra_halo.begin(),
-                         extra_halo.end());
+        dofs_with_support_on_halo_cells.insert(extra_halo.begin(),
+                                               extra_halo.end());
       }
 
     // Rework our std::set's into IndexSet and subtract halo DoFs from the
     // predicate_dofs:
     IndexSet support_set(dof_handler.n_dofs());
-    support_set.add_indices(predicate_dofs.begin(), predicate_dofs.end());
+    support_set.add_indices(predicate_dofs.begin(),
+                            predicate_dofs.end());
     support_set.compress();
 
     IndexSet halo_set(dof_handler.n_dofs());
-    halo_set.add_indices(halo_dofs.begin(), halo_dofs.end());
+    halo_set.add_indices(dofs_with_support_on_halo_cells.begin(),
+                         dofs_with_support_on_halo_cells.end());
     halo_set.compress();
 
     support_set.subtract_set(halo_set);
 
     // we intentionally do not want to limit the output to locally owned DoFs.
-    // In other words, the output may contain DoFs that are on the locally
-    // owned cells, but are not owned by this MPI core.
-
     return support_set;
   }
 
index 08d536ae16c1f06b7667460824f7da67a0b3d5cc..6bf1ac92936eaaeca5a72d538b53726b84deeda6 100644 (file)
@@ -162,7 +162,7 @@ for (deal_II_dimension : DIMENSIONS)
 
     template
     IndexSet
-    DoFTools::extract_dofs_with_support_on<DoFHandler<deal_II_dimension> >
+    DoFTools::extract_dofs_with_support_contained_within<DoFHandler<deal_II_dimension> >
     (const DoFHandler<deal_II_dimension> &,
      const std::function< bool(const typename DoFHandler<deal_II_dimension>::active_cell_iterator &)> &,
      const ConstraintMatrix&);
similarity index 97%
rename from tests/dofs/dof_tools_extract_dofs_with_support_on_01.cc
rename to tests/dofs/dof_tools_extract_dofs_with_support_contained_within_01.cc
index 741020cd8c35b9aca9d6641a5af1232e5aaf3bfd..a6cbc171f6ad4c8baaefdb8bf5bb0d82612da242 100644 (file)
@@ -90,7 +90,7 @@ void test (const unsigned int flag)
   ConstraintMatrix cm;
   DoFTools::make_hanging_node_constraints(dh,cm);
 
-  IndexSet support = DoFTools::extract_dofs_with_support_on(dh, pred_d<dim>, cm);
+  IndexSet support = DoFTools::extract_dofs_with_support_contained_within(dh, pred_d<dim>, cm);
   support.print(deallog);
 
   // print grid and DoFs for visual inspection
similarity index 94%
rename from tests/dofs/dof_tools_extract_dofs_with_support_on_02.cc
rename to tests/dofs/dof_tools_extract_dofs_with_support_contained_within_02.cc
index d74ba742fb1f83994eb4b14cd5c498d3dc8c3bb4..9c2431a50ad8ed84bb673ffbf77f214fe6bbcf31 100644 (file)
@@ -14,9 +14,9 @@
 // ---------------------------------------------------------------------
 //
 
-// similar to dof_tools_23, but for the p::d::Tria. As an MPI test, make sure
-// that the total number of shape functions that are non-zero within the domain
-// is the same.
+// test DoFTools::extract_dofs_with_support_contained_within() for the p::d::Tria.
+// As an MPI test, make sure that the total number of shape functions that are
+// non-zero within the domain is the same.
 
 #include "../tests.h"
 #include <deal.II/distributed/tria.h>
@@ -102,7 +102,7 @@ void test (const unsigned int flag)
   DoFTools::make_hanging_node_constraints(dh,cm);
   cm.close ();
 
-  const IndexSet support = DoFTools::extract_dofs_with_support_on(dh, pred_d<dim>, cm);
+  const IndexSet support = DoFTools::extract_dofs_with_support_contained_within(dh, pred_d<dim>, cm);
   const IndexSet support_local = support & dh.locally_owned_dofs();
 
   deallog << support.n_elements() << std::endl;
similarity index 93%
rename from tests/dofs/dof_tools_extract_dofs_with_support_on_03.cc
rename to tests/dofs/dof_tools_extract_dofs_with_support_contained_within_03.cc
index e19e4ec03662e37e9a3120a5e52fb2adcbc26b6c..91d2131322b61ffb1f3e3a0fa4699231f6eccbff 100644 (file)
@@ -14,7 +14,7 @@
 // ---------------------------------------------------------------------
 //
 
-// same as dof_tools_32 test DoFTools::extract_dofs_with_support_on() but
+// test DoFTools::extract_dofs_with_support_contained_within() but
 // for a slightly different configuration
 
 #include "../tests.h"
@@ -96,8 +96,8 @@ void test (const bool left = true)
   DoFTools::make_hanging_node_constraints(dh,cm);
 
   IndexSet support = left ?
-                     DoFTools::extract_dofs_with_support_on(dh, pred_left<dim>, cm) :
-                     DoFTools::extract_dofs_with_support_on(dh, pred_right<dim>, cm);
+                     DoFTools::extract_dofs_with_support_contained_within(dh, pred_left<dim>, cm) :
+                     DoFTools::extract_dofs_with_support_contained_within(dh, pred_right<dim>, cm);
   support.print(deallog);
 
   // print grid and DoFs for visual inspection
similarity index 98%
rename from tests/dofs/dof_tools_extract_dofs_with_support_on_04.cc
rename to tests/dofs/dof_tools_extract_dofs_with_support_contained_within_04.cc
index 7126ab055897245c4adeef3f61ceb7f9918ce563..38bd6fa4545cf2cdafd04035b17d213980796bd8 100644 (file)
@@ -14,8 +14,8 @@
 // ---------------------------------------------------------------------
 //
 
-// test DoFTools::extract_dofs_with_support_on() for calculation of the RHS
-// function in parallel
+// test DoFTools::extract_dofs_with_support_contained_within() for calculation
+// of the RHS function in parallel
 
 #include "../tests.h"
 #include <deal.II/base/logstream.h>
@@ -111,7 +111,7 @@ void test ()
   std::vector<types::global_dof_index> local_dof_indices(fe.dofs_per_cell);
 
   // get support on the predicate
-  IndexSet support = DoFTools::extract_dofs_with_support_on(dh, pred_d<dim>, cm);
+  IndexSet support = DoFTools::extract_dofs_with_support_contained_within(dh, pred_d<dim>, cm);
   IndexSet local_support = support & locally_owned_set;
 
   // rhs vectors:

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.