From 8728789ac866dc5646b5e229f426cbfdb43828fc Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Mon, 18 Jul 2016 17:37:12 +0200 Subject: [PATCH] Simplify MGConstrainedDoFs --- .../deal.II/multigrid/mg_constrained_dofs.h | 33 ++++++++++++++++++- include/deal.II/multigrid/mg_tools.h | 11 +++++++ source/multigrid/mg_tools.cc | 31 +++++++++++++++++ source/multigrid/mg_tools.inst.in | 7 ++++ 4 files changed, 81 insertions(+), 1 deletion(-) diff --git a/include/deal.II/multigrid/mg_constrained_dofs.h b/include/deal.II/multigrid/mg_constrained_dofs.h index 6a839a7378..a2706716a1 100644 --- a/include/deal.II/multigrid/mg_constrained_dofs.h +++ b/include/deal.II/multigrid/mg_constrained_dofs.h @@ -65,7 +65,19 @@ public: template void initialize(const DoFHandler &dof, const typename FunctionMap::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 + void set_zero_boundary_dofs(const DoFHandler &dof, + const std::vector &boundary_indicators, + const ComponentMask &component_mask = ComponentMask()); /** * Reset the data structures. @@ -160,6 +172,25 @@ MGConstrainedDoFs::initialize(const DoFHandler &dof, } +template +inline +void +MGConstrainedDoFs::set_zero_boundary_dofs(const DoFHandler &dof, + const std::vector &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() diff --git a/include/deal.II/multigrid/mg_tools.h b/include/deal.II/multigrid/mg_tools.h index b7c148f74c..dae7f34d92 100644 --- a/include/deal.II/multigrid/mg_tools.h +++ b/include/deal.II/multigrid/mg_tools.h @@ -197,6 +197,17 @@ namespace MGTools std::vector &boundary_indices, const ComponentMask &component_mask = ComponentMask()); + /** + * The same function as above, but return an IndexSet rather than a + * std::set on each level. + */ + template + void + make_boundary_list (const DoFHandler &mg_dof, + const std::vector &boundary_indicators, + std::vector &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 diff --git a/source/multigrid/mg_tools.cc b/source/multigrid/mg_tools.cc index bea3e5779b..183114aa0f 100644 --- a/source/multigrid/mg_tools.cc +++ b/source/multigrid/mg_tools.cc @@ -1432,6 +1432,37 @@ namespace MGTools } + template + void + make_boundary_list(const DoFHandler &dof, + const std::vector &boundary_indicators, + std::vector &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 > + my_boundary_indices (dof.get_triangulation().n_global_levels()); + + //TODO get rid off this + typename FunctionMap::type homogen_bc; + ZeroFunction zero_function (dim); + for (unsigned int i=0;i void extract_non_interface_dofs (const DoFHandler &mg_dof_handler, diff --git a/source/multigrid/mg_tools.inst.in b/source/multigrid/mg_tools.inst.in index 329b7ce229..26a1b64545 100644 --- a/source/multigrid/mg_tools.inst.in +++ b/source/multigrid/mg_tools.inst.in @@ -102,6 +102,13 @@ for (deal_II_dimension : DIMENSIONS) std::vector&, const ComponentMask &); + template void make_boundary_list( + const DoFHandler &dof, + const std::vector &boundary_indicators, + std::vector &boundary_indices, + const ComponentMask &component_mask); + + template void extract_inner_interface_dofs (const DoFHandler &mg_dof_handler, -- 2.39.5