From 7158719c5e7cfaf4fe55403b13191e92b93017d3 Mon Sep 17 00:00:00 2001 From: Bruno Blais Date: Mon, 2 Dec 2019 18:45:48 -0500 Subject: [PATCH] Refactored map_dofs_to_support_points to allow for the use of a mask and removed the additional function that had been introduced --- include/deal.II/dofs/dof_tools.h | 69 +++++--------------------------- source/dofs/dof_tools.cc | 56 ++++++++++++++++++++------ source/dofs/dof_tools.inst.in | 26 +++++++----- 3 files changed, 69 insertions(+), 82 deletions(-) diff --git a/include/deal.II/dofs/dof_tools.h b/include/deal.II/dofs/dof_tools.h index 25a07560c2..fcfe83f181 100644 --- a/include/deal.II/dofs/dof_tools.h +++ b/include/deal.II/dofs/dof_tools.h @@ -2310,7 +2310,8 @@ namespace DoFTools void map_dofs_to_support_points(const Mapping & mapping, const DoFHandler &dof_handler, - std::vector> & support_points); + std::vector> & support_points, + const ComponentMask &mask = ComponentMask()); /** * Same as the previous function but for the hp case. @@ -2321,7 +2322,8 @@ namespace DoFTools map_dofs_to_support_points( const dealii::hp::MappingCollection &mapping, const hp::DoFHandler & dof_handler, - std::vector> & support_points); + std::vector> & support_points, + const ComponentMask & mask = ComponentMask()); /** * This function is a version of the above map_dofs_to_support_points @@ -2349,70 +2351,16 @@ namespace DoFTools * @param support_points A map that for every locally relevant DoF index * contains the corresponding location in real space coordinates. Previous * content of this object is deleted in this function. + * @param component_mask An optional component mask that restricts the + * components from which the support points are extracted */ - template - void - map_dofs_to_support_points( - const Mapping & mapping, - const DoFHandler & dof_handler, - std::map> &support_points); - - // Same as above but with support for a component mask - template void map_dofs_to_support_points( const Mapping & mapping, const DoFHandler & dof_handler, std::map> &support_points, - const ComponentMask & mask) - { - const auto &tria = dof_handler.get_triangulation(); - const auto &fe = dof_handler.get_fe(); - - // Check if the FE has support points - Assert(fe.has_support_points(), - typename FiniteElement::ExcFEHasNoSupportPoints()); - - // Create quadrature at support point - Quadrature quad(fe.get_unit_support_points()); - - // Now loop over all cells and enquire the support points on each - // of these. we use dummy quadrature formulas where the quadrature - // points are located at the unit support points to enquire the - // location of the support points in real space. - // - // The weights of the quadrature rule have been set to invalid - // values by the used constructor. - // - // The mask is used to assess if a component should be kept - // or not - FEValues fe_values(mapping, - fe, - quad, - update_quadrature_points); - - std::vector local_dof_indices; - for (const auto &cell : dof_handler.active_cell_iterators()) - // only work on locally relevant cells - if (cell->is_artificial() == false) - { - fe_values.reinit(cell); - - local_dof_indices.resize(cell->get_fe().dofs_per_cell); - cell->get_dof_indices(local_dof_indices); - - const std::vector> &points = - fe_values.get_quadrature_points(); - for (unsigned int i = 0; i < cell->get_fe().dofs_per_cell; ++i) - { - unsigned int dof_comp = fe.system_to_component_index(i).first; - // insert the values into the map - if (mask[dof_comp]) - support_points[local_dof_indices[i]] = points[i]; - } - } - } + const ComponentMask & mask = ComponentMask()); /** * Same as the previous function but for the hp case. @@ -2422,7 +2370,8 @@ namespace DoFTools map_dofs_to_support_points( const dealii::hp::MappingCollection &mapping, const hp::DoFHandler & dof_handler, - std::map> &support_points); + std::map> &support_points, + const ComponentMask & mask = ComponentMask()); /** diff --git a/source/dofs/dof_tools.cc b/source/dofs/dof_tools.cc index be74d9c8f2..1ffb78db1e 100644 --- a/source/dofs/dof_tools.cc +++ b/source/dofs/dof_tools.cc @@ -2077,7 +2077,8 @@ namespace DoFTools DoFHandlerType::space_dimension> &mapping, const DoFHandlerType & dof_handler, std::map> &support_points) + Point> &support_points, + const ComponentMask & in_mask) { const unsigned int dim = DoFHandlerType::dimension; const unsigned int spacedim = DoFHandlerType::space_dimension; @@ -2096,6 +2097,12 @@ namespace DoFTools fe_collection[fe_index].get_unit_support_points())); } + // Take care of components + const ComponentMask mask = + (in_mask.size() == 0 ? + ComponentMask(fe_collection.n_components(), true) : + in_mask); + // Now loop over all cells and enquire the support points on each // of these. we use dummy quadrature formulas where the quadrature // points are located at the unit support points to enquire the @@ -2126,8 +2133,14 @@ namespace DoFTools const std::vector> &points = fe_values.get_quadrature_points(); for (unsigned int i = 0; i < cell->get_fe().dofs_per_cell; ++i) - // insert the values into the map - support_points[local_dof_indices[i]] = points[i]; + { + unsigned int dof_comp = + cell->get_fe().system_to_component_index(i).first; + + // insert the values into the map if it is a valid component + if (mask[dof_comp]) + support_points[local_dof_indices[i]] = points[i]; + } } } @@ -2138,13 +2151,17 @@ namespace DoFTools const hp::MappingCollection &mapping, const DoFHandlerType & dof_handler, - std::vector> &support_points) + std::vector> &support_points, + const ComponentMask & mask) { // get the data in the form of the map as above std::map> x_support_points; - map_dofs_to_support_points(mapping, dof_handler, x_support_points); + map_dofs_to_support_points(mapping, + dof_handler, + x_support_points, + mask); // now convert from the map to the linear vector. make sure every // entry really appeared in the map @@ -2152,6 +2169,7 @@ namespace DoFTools { Assert(x_support_points.find(i) != x_support_points.end(), ExcInternalError()); + support_points[i] = x_support_points[i]; } } @@ -2162,7 +2180,8 @@ namespace DoFTools void map_dofs_to_support_points(const Mapping & mapping, const DoFHandler &dof_handler, - std::vector> & support_points) + std::vector> & support_points, + const ComponentMask & mask) { AssertDimension(support_points.size(), dof_handler.n_dofs()); Assert((dynamic_cast< @@ -2178,7 +2197,8 @@ namespace DoFTools internal::map_dofs_to_support_points(mapping_collection, dof_handler, - support_points); + support_points, + mask); } @@ -2187,7 +2207,8 @@ namespace DoFTools map_dofs_to_support_points( const hp::MappingCollection &mapping, const hp::DoFHandler & dof_handler, - std::vector> & support_points) + std::vector> & support_points, + const ComponentMask & mask) { AssertDimension(support_points.size(), dof_handler.n_dofs()); Assert((dynamic_cast< @@ -2199,7 +2220,10 @@ namespace DoFTools // Let the internal function do all the work, just make sure that it // gets a MappingCollection - internal::map_dofs_to_support_points(mapping, dof_handler, support_points); + internal::map_dofs_to_support_points(mapping, + dof_handler, + support_points, + mask); } @@ -2208,7 +2232,8 @@ namespace DoFTools map_dofs_to_support_points( const Mapping & mapping, const DoFHandler & dof_handler, - std::map> &support_points) + std::map> &support_points, + const ComponentMask & mask) { support_points.clear(); @@ -2218,7 +2243,8 @@ namespace DoFTools internal::map_dofs_to_support_points(mapping_collection, dof_handler, - support_points); + support_points, + mask); } @@ -2227,13 +2253,17 @@ namespace DoFTools map_dofs_to_support_points( const hp::MappingCollection & mapping, const hp::DoFHandler & dof_handler, - std::map> &support_points) + std::map> &support_points, + const ComponentMask & mask) { support_points.clear(); // Let the internal function do all the work, just make sure that it // gets a MappingCollection - internal::map_dofs_to_support_points(mapping, dof_handler, support_points); + internal::map_dofs_to_support_points(mapping, + dof_handler, + support_points, + mask); } template diff --git a/source/dofs/dof_tools.inst.in b/source/dofs/dof_tools.inst.in index bcde5579b4..f87fe38a54 100644 --- a/source/dofs/dof_tools.inst.in +++ b/source/dofs/dof_tools.inst.in @@ -489,25 +489,29 @@ for (deal_II_dimension : DIMENSIONS) template void DoFTools::map_dofs_to_support_points( const Mapping &, const DoFHandler &, - std::vector> &); + std::vector> &, + const ComponentMask &); template void DoFTools::map_dofs_to_support_points( const hp::MappingCollection &, const hp::DoFHandler &, - std::vector> &); + std::vector> &, + const ComponentMask &); template void DoFTools::map_dofs_to_support_points( const Mapping &, const DoFHandler &, - std::map> &); + std::map> &, + const ComponentMask &); template void DoFTools::map_dofs_to_support_points( const hp::MappingCollection &, const hp::DoFHandler &, - std::map> &); + std::map> &, + const ComponentMask &); #if deal_II_dimension < 3 @@ -515,25 +519,29 @@ for (deal_II_dimension : DIMENSIONS) deal_II_dimension + 1>( const Mapping &, const DoFHandler &, - std::vector> &); + std::vector> &, + const ComponentMask &); template void DoFTools::map_dofs_to_support_points( const hp::MappingCollection &, const hp::DoFHandler &, - std::vector> &); + std::vector> &, + const ComponentMask &); template void DoFTools::map_dofs_to_support_points( const Mapping &, const DoFHandler &, - std::map> &); + std::map> &, + const ComponentMask &); template void DoFTools::map_dofs_to_support_points( const hp::MappingCollection &, const hp::DoFHandler &, - std::map> &); + std::map> &, + const ComponentMask &); template void DoFTools::count_dofs_per_block< @@ -554,7 +562,7 @@ for (deal_II_dimension : DIMENSIONS) #if deal_II_dimension == 3 template void DoFTools::map_dofs_to_support_points<1, 3>( - const Mapping<1, 3> &, const DoFHandler<1, 3> &, std::vector> &); + const Mapping<1, 3> &, const DoFHandler<1, 3> &, std::vector> &,const ComponentMask &); template void DoFTools::count_dofs_per_block>( const DoFHandler<1, 3> &, -- 2.39.5