]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Use std::set for boundary_indicators and rename set_zero_boundary_dofs 2847/head
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Wed, 20 Jul 2016 09:07:17 +0000 (11:07 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sat, 23 Jul 2016 19:05:09 +0000 (21:05 +0200)
doc/news/changes.h
include/deal.II/multigrid/mg_constrained_dofs.h
include/deal.II/multigrid/mg_tools.h
source/multigrid/mg_tools.cc
source/multigrid/mg_tools.inst.in
tests/multigrid/constrained_dofs_03.cc

index df133e88fc506c85812b8cd20f96c5d7d49a58ec..bc3106010ed0733abcf11b7b7eb8d05e2c5b1c41 100644 (file)
@@ -360,9 +360,10 @@ inconvenience this causes.
 <h3>Specific improvements</h3>
 
 <ol>
- <li> Improved: Replace by a std::vector for initializing
- the constrained boundary DoFs, since the function values
- where not used. Allow for non-primitive FiniteElements.   
+ <li> Improved: Allow for initializing the constrained
+ boundary DoFs in MGConstrainedDoFs using a std::set
+ instead of a FunctionMap whose function values were not used.
+ Allow for non-primitive FiniteElements.
  <br>
  (Daniel Arndt, 2016/07/20)
  </li>
index a2706716a1ad63aa306fa08b2c88be0d691538c0..68103f605a8c9650191a0473cfa7ba9303e4e354 100644 (file)
@@ -61,6 +61,8 @@ public:
    * This function internally calls the initialize() function above and the
    * constrains degrees on the external boundary of the domain by calling
    * MGTools::make_boundary_list() with the given second and third argument.
+   *
+   * @deprecated Use initialize() followed by set_zero_boundary_dofs() instead
    */
   template <int dim, int spacedim>
   void initialize(const DoFHandler<dim,spacedim> &dof,
@@ -68,16 +70,19 @@ public:
                   const ComponentMask &component_mask = ComponentMask()) DEAL_II_DEPRECATED;
 
   /**
-   * Fill the internal data structures with information about interface and boundary dofs.
+   * Fill the internal data structures with information
+   * about Dirichlet boundary dofs.
    *
-   * This function internally calls the initialize() function above and the
-   * constrains degrees on the external boundary of the domain by calling
-   * MGTools::make_boundary_list() with the given second and third argument.
+   * The initialize() function has to be called before
+   * to set hanging node constraints.
+   *
+   * This function can be called multiple times to allow considering
+   * different sets of boundary_ids for different components.
    */
   template <int dim, int spacedim>
-  void set_zero_boundary_dofs(const DoFHandler<dim,spacedim> &dof,
-                              const std::vector<types::boundary_id> &boundary_indicators,
-                              const ComponentMask &component_mask = ComponentMask());
+  void make_zero_boundary_constraints(const DoFHandler<dim,spacedim> &dof,
+                                      const std::set<types::boundary_id> &boundary_ids,
+                                      const ComponentMask &component_mask = ComponentMask());
 
   /**
    * Reset the data structures.
@@ -140,9 +145,11 @@ void
 MGConstrainedDoFs::initialize(const DoFHandler<dim,spacedim> &dof)
 {
   boundary_indices.clear();
+  refinement_edge_indices.clear();
 
   const unsigned int nlevels = dof.get_triangulation().n_global_levels();
 
+  //At this point refinement_edge_indices is empty.
   refinement_edge_indices.resize(nlevels);
   for (unsigned int l=0; l<nlevels; ++l)
     refinement_edge_indices[l] = IndexSet(dof.n_dofs(l));
@@ -163,6 +170,7 @@ MGConstrainedDoFs::initialize(const DoFHandler<dim,spacedim> &dof,
   // allocate an IndexSet for each global level. Contents will be
   // overwritten inside make_boundary_list.
   const unsigned int n_levels = dof.get_triangulation().n_global_levels();
+  //At this point boundary_indices is empty.
   boundary_indices.resize(n_levels);
 
   MGTools::make_boundary_list (dof,
@@ -175,17 +183,19 @@ MGConstrainedDoFs::initialize(const DoFHandler<dim,spacedim> &dof,
 template <int dim, int spacedim>
 inline
 void
-MGConstrainedDoFs::set_zero_boundary_dofs(const DoFHandler<dim,spacedim> &dof,
-                                          const std::vector<types::boundary_id> &boundary_indicators,
-                                          const ComponentMask &component_mask)
+MGConstrainedDoFs::make_zero_boundary_constraints(const DoFHandler<dim,spacedim> &dof,
+                                                  const std::set<types::boundary_id> &boundary_ids,
+                                                  const ComponentMask &component_mask)
 {
   // allocate an IndexSet for each global level. Contents will be
   // overwritten inside make_boundary_list.
   const unsigned int n_levels = dof.get_triangulation().n_global_levels();
+  Assert(boundary_indices.size() == 0 || boundary_indices.size() == n_levels,
+         ExcInternalError());
   boundary_indices.resize(n_levels);
 
   MGTools::make_boundary_list (dof,
-                               boundary_indicators,
+                               boundary_ids,
                                boundary_indices,
                                component_mask);
 }
index 816ef7500cdd56ff53ad99db7e875bce56e567f3..dfc7b8b67d8b7a4e1f934daf21564adbb04867b1 100644 (file)
@@ -178,6 +178,9 @@ namespace MGTools
    * indices of degrees of freedom for each level that are at the part of the
    * boundary identified by the function_map argument. Its length has to match
    * the number of levels in the dof handler object.
+   *
+   * Previous content in @p boundary_indices is not overwritten,
+   * but added to.
    */
   template <int dim, int spacedim>
   void
@@ -189,6 +192,9 @@ namespace MGTools
   /**
    * The same function as above, but return an IndexSet rather than a
    * std::set<unsigned int> on each level.
+   *
+   * Previous content in @p boundary_indices is not overwritten,
+   * but added to.
    */
   template <int dim, int spacedim>
   void
@@ -199,12 +205,16 @@ namespace MGTools
 
   /**
   * The same function as above, but return an IndexSet rather than a
-  * std::set<unsigned int> on each level.
+  * std::set<unsigned int> on each level and use a std::set of
+  * boundary_ids as input.
+  *
+  * Previous content in @p boundary_indices is not overwritten,
+  * but added to.
   */
   template <int dim, int spacedim>
   void
   make_boundary_list (const DoFHandler<dim,spacedim>      &mg_dof,
-                      const std::vector<types::boundary_id> &boundary_indicators,
+                      const std::set<types::boundary_id> &boundary_ids,
                       std::vector<IndexSet>                 &boundary_indices,
                       const ComponentMask               &component_mask = ComponentMask());
 
index a65f2cb6e3f0a2aa4d9cd48c23d94bf40b9d628d..387d5903d1b271e7d771763af2641d403800ef1b 100644 (file)
@@ -1179,20 +1179,6 @@ namespace MGTools
   }
 
 
-
-  template <>
-  void
-  make_boundary_list(
-    const DoFHandler<1,2> &,
-    const FunctionMap<1>::type &,
-    std::vector<std::set<types::global_dof_index> > &,
-    const ComponentMask &)
-  {
-    Assert(false, ExcNotImplemented());
-  }
-
-
-
   template <int dim, int spacedim>
   void
   make_boundary_list(
@@ -1205,20 +1191,17 @@ namespace MGTools
             ExcDimensionMismatch (boundary_indices.size(),
                                   dof.get_triangulation().n_global_levels()));
 
-    std::vector<types::boundary_id> boundary_indicators;
+    std::set<types::boundary_id> boundary_ids;
     typename std::map<types::boundary_id, const Function<dim>* >::const_iterator it
       = function_map.begin();
     for (; it!=function_map.end(); ++it)
-      boundary_indicators.push_back(it->first);
+      boundary_ids.insert(it->first);
 
     std::vector<IndexSet>  boundary_indexset;
-    make_boundary_list(dof, boundary_indicators, boundary_indexset, component_mask);
+    make_boundary_list(dof, boundary_ids, boundary_indexset, component_mask);
     for (unsigned int i=0; i<dof.get_triangulation().n_global_levels(); ++i)
-      {
-        IndexSet::ElementIterator it_element = boundary_indexset[i].begin();
-        for (; it_element != boundary_indexset[i].end(); ++it_element)
-          boundary_indices[i].insert(*it_element);
-      }
+      boundary_indices[i].insert(boundary_indexset[i].begin(),
+                                 boundary_indexset[i].end());
   }
 
 
@@ -1233,12 +1216,12 @@ namespace MGTools
             ExcDimensionMismatch (boundary_indices.size(),
                                   dof.get_triangulation().n_global_levels()));
 
-    std::vector<types::boundary_id> boundary_indicators;
+    std::set<types::boundary_id> boundary_ids;
     typename std::map<types::boundary_id, const Function<dim>* >::const_iterator it = function_map.begin();
     for (; it!=function_map.end(); ++it)
-      boundary_indicators.push_back(it->first);
+      boundary_ids.insert(it->first);
 
-    make_boundary_list (dof, boundary_indicators, boundary_indices, component_mask);
+    make_boundary_list (dof, boundary_ids, boundary_indices, component_mask);
   }
 
 
@@ -1246,20 +1229,7 @@ namespace MGTools
   void
   make_boundary_list(
     const DoFHandler<1,1> &,
-    const std::vector<types::boundary_id> &,
-    std::vector<IndexSet> &,
-    const ComponentMask &)
-  {
-    Assert(false, ExcNotImplemented());
-  }
-
-
-
-  template <>
-  void
-  make_boundary_list(
-    const DoFHandler<1,2> &,
-    const std::vector<types::boundary_id> &,
+    const std::set<types::boundary_id> &,
     std::vector<IndexSet> &,
     const ComponentMask &)
   {
@@ -1270,14 +1240,14 @@ namespace MGTools
   template <int dim, int spacedim>
   void
   make_boundary_list(const DoFHandler<dim,spacedim> &dof,
-                     const std::vector<types::boundary_id> &boundary_indicators,
+                     const std::set<types::boundary_id> &boundary_ids,
                      std::vector<IndexSet> &boundary_indices,
                      const ComponentMask &component_mask)
   {
     boundary_indices.resize(dof.get_triangulation().n_global_levels());
 
-    // if for whatever reason we were passed an empty vector, return immediately
-    if (boundary_indicators.size() == 0)
+    // if for whatever reason we were passed an empty set, return immediately
+    if (boundary_ids.size() == 0)
       return;
 
     for (unsigned int i=0; i<dof.get_triangulation().n_global_levels(); ++i)
@@ -1293,7 +1263,7 @@ namespace MGTools
                local_dofs.end (),
                DoFHandler<dim,spacedim>::invalid_dof_index);
 
-    // First, deal with the simpler case when we have to identif all boundary dofs
+    // First, deal with the simpler case when we have to identify all boundary dofs
     if (component_mask.n_selected_components(n_components) == n_components)
       {
         typename DoFHandler<dim,spacedim>::cell_iterator
@@ -1316,7 +1286,7 @@ namespace MGTools
                   face = cell->face(face_no);
                   const types::boundary_id bi = face->boundary_id();
                   // Face is listed in boundary map
-                  if (std::find(boundary_indicators.begin(), boundary_indicators.end(), bi) != boundary_indicators.end())
+                  if (boundary_ids.find(bi) != boundary_ids.end())
                     {
                       face->get_mg_dof_indices(level, local_dofs);
                       boundary_indices[level].add_indices(local_dofs.begin(), local_dofs.end());
@@ -1346,8 +1316,7 @@ namespace MGTools
 
                 typename DoFHandler<dim,spacedim>::face_iterator face = cell->face(face_no);
                 const types::boundary_id boundary_component = face->boundary_id();
-                if (std::find(boundary_indicators.begin(), boundary_indicators.end(), boundary_component)
-                    != boundary_indicators.end())
+                if (boundary_ids.find(boundary_component) != boundary_ids.end())
                   // we want to constrain this boundary
                   {
 
@@ -1378,9 +1347,6 @@ namespace MGTools
                     face->get_mg_dof_indices (level, local_dofs);
                     if (fe_is_system)
                       {
-                        // enter those dofs into the list that match the component
-                        // signature. avoid the usual complication that we can't
-                        // just use *_system_to_component_index for non-primitive FEs
                         for (unsigned int i=0; i<local_dofs.size(); ++i)
                           {
                             unsigned int component = numbers::invalid_unsigned_int;
@@ -1406,8 +1372,7 @@ namespace MGTools
                           }
                       }
                     else
-                      for (unsigned int i=0; i<local_dofs.size(); ++i)
-                        boundary_indices[level].add_indices(local_dofs.begin(), local_dofs.end());
+                      boundary_indices[level].add_indices(local_dofs.begin(), local_dofs.end());
                   }
               }
       }
index 26a1b645458c6e51a90b545f3fab74cad80733f6..2876c42cda4c69d738c973e7f83585237926fde4 100644 (file)
@@ -103,10 +103,10 @@ for (deal_II_dimension : DIMENSIONS)
        const ComponentMask &);
 
     template void  make_boundary_list(
-       const DoFHandler<deal_II_dimension> &dof,
-       const std::vector<types::boundary_id> &boundary_indicators,
-        std::vector<IndexSet> &boundary_indices,
-        const ComponentMask &component_mask);
+        const DoFHandler<deal_II_dimension> &,
+        const std::set<types::boundary_id> &,
+        std::vector<IndexSet> &,
+        const ComponentMask &);
 
 
       template
index 6f55029b9688f11d866cc40f3e5345f8dd649750..4ce74154d83fad1ac10670b6745b28750b51cc25 100644 (file)
@@ -55,10 +55,12 @@ void check_fe(FiniteElement<dim> &fe, ComponentMask &component_mask)
   dofh.distribute_dofs(fe);
   dofh.distribute_mg_dofs(fe);
 
-  MGConstrainedDoFs                    mg_constrained_dofs;
-  std::vector<types::boundary_id>      boundary_indicators(1,0);
+  MGConstrainedDoFs                 mg_constrained_dofs;
+  std::set<types::boundary_id>      boundary_indicators;
+
+  boundary_indicators.insert(0);
   mg_constrained_dofs.initialize(dofh);
-  mg_constrained_dofs.set_zero_boundary_dofs(dofh, boundary_indicators, component_mask);
+  mg_constrained_dofs.make_zero_boundary_constraints(dofh, boundary_indicators, component_mask);
 
   const unsigned int n_levels = tr.n_global_levels();
   for (unsigned int level = 0; level < n_levels; ++level)

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.