]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Simplify MGConstrainedDoFs
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Mon, 18 Jul 2016 15:37:12 +0000 (17:37 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Wed, 20 Jul 2016 00:22:49 +0000 (02:22 +0200)
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

index 6a839a7378f9235dfe7c3272904eff13225774a3..a2706716a1ad63aa306fa08b2c88be0d691538c0 100644 (file)
@@ -65,7 +65,19 @@ public:
   template <int dim, int spacedim>
   void initialize(const DoFHandler<dim,spacedim> &dof,
                   const typename FunctionMap<dim>::type &function_map,
-                  const ComponentMask &component_mask = ComponentMask());
+                  const ComponentMask &component_mask = ComponentMask()) DEAL_II_DEPRECATED;
+
+  /**
+   * Fill the internal data structures with information about interface and 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.
+   */
+  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());
 
   /**
    * Reset the data structures.
@@ -160,6 +172,25 @@ 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)
+{
+  // 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();
+  boundary_indices.resize(n_levels);
+
+  MGTools::make_boundary_list (dof,
+                               boundary_indicators,
+                               boundary_indices,
+                               component_mask);
+}
+
+
 inline
 void
 MGConstrainedDoFs::clear()
index b7c148f74c8b5d6c72509c68b6fd52fc3c98a830..dae7f34d926739e9e2100b7cb8cccb2952b45e2b 100644 (file)
@@ -197,6 +197,17 @@ namespace MGTools
                       std::vector<IndexSet>                 &boundary_indices,
                       const ComponentMask               &component_mask = ComponentMask());
 
+   /**
+   * The same function as above, but return an IndexSet rather than a
+   * std::set<unsigned int> on each level.
+   */
+  template <int dim, int spacedim>
+  void
+  make_boundary_list (const DoFHandler<dim,spacedim>      &mg_dof,
+                      const std::vector<types::boundary_id> &boundary_indicators,
+                      std::vector<IndexSet>                 &boundary_indices,
+                      const ComponentMask               &component_mask = ComponentMask());
+
   /**
    * For each level in a multigrid hierarchy, produce an IndexSet that
    * indicates which of the degrees of freedom are along interfaces of this
index bea3e5779bd558b343b58c6c3c1fdf12e3e2a739..183114aa0fc099028175c5c6601b0653a828334b 100644 (file)
@@ -1432,6 +1432,37 @@ namespace MGTools
   }
 
 
+  template <int dim, int spacedim>
+  void
+  make_boundary_list(const DoFHandler<dim,spacedim> &dof,
+                     const std::vector<types::boundary_id> &boundary_indicators,
+                     std::vector<IndexSet> &boundary_indices,
+                     const ComponentMask &component_mask)
+  {
+    Assert (boundary_indices.size() == dof.get_triangulation().n_global_levels(),
+            ExcDimensionMismatch (boundary_indices.size(),
+                                  dof.get_triangulation().n_global_levels()));
+
+    std::vector<std::set<types::global_dof_index> >
+    my_boundary_indices (dof.get_triangulation().n_global_levels());
+
+    //TODO get rid off this
+    typename FunctionMap<dim>::type homogen_bc;
+    ZeroFunction<dim> zero_function (dim);
+    for (unsigned int i=0;i<boundary_indicators.size();++i)
+       homogen_bc[i] = &zero_function;
+
+    make_boundary_list (dof, homogen_bc, my_boundary_indices, component_mask);
+    for (unsigned int i=0; i<dof.get_triangulation().n_global_levels(); ++i)
+      {
+        if (boundary_indices[i].size()!=dof.n_dofs(i))
+          boundary_indices[i] = IndexSet (dof.n_dofs(i));
+        boundary_indices[i].add_indices (my_boundary_indices[i].begin(),
+                                         my_boundary_indices[i].end());
+      }
+  }
+
+
   template <int dim, int spacedim>
   void
   extract_non_interface_dofs (const DoFHandler<dim,spacedim> &mg_dof_handler,
index 329b7ce229ac35ff6bda396fbd787de7cb09b709..26a1b645458c6e51a90b545f3fab74cad80733f6 100644 (file)
@@ -102,6 +102,13 @@ for (deal_II_dimension : DIMENSIONS)
        std::vector<IndexSet>&,
        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);
+
+
       template
        void
        extract_inner_interface_dofs (const DoFHandler<deal_II_dimension> &mg_dof_handler,

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.