From 6446e49ec22ba05c8eaedc523e71f3f067d8b22b Mon Sep 17 00:00:00 2001 From: Spencer Patty Date: Tue, 9 Feb 2016 11:59:43 -0600 Subject: [PATCH] change template parameter to DoFHandlerType and make various changes as recommended in pull request --- include/deal.II/grid/grid_tools.h | 51 +++++++++++++++---------------- source/grid/grid_tools.cc | 44 +++++++++++++------------- 2 files changed, 47 insertions(+), 48 deletions(-) diff --git a/include/deal.II/grid/grid_tools.h b/include/deal.II/grid/grid_tools.h index 1dc398aff2..7a9a985c7e 100644 --- a/include/deal.II/grid/grid_tools.h +++ b/include/deal.II/grid/grid_tools.h @@ -1081,44 +1081,43 @@ namespace GridTools /** * This function runs through the degrees of freedom defined by the - * Container and for each dof constructs a vector of active_cell_iterators - * representing the cells of support of the nodal basis element - * at that degree of freedom. This function was designed for the + * DoFHandlerType and for each dof constructs a vector of active_cell_iterators + * representing the cells of support of the associated basis element + * at that degree of freedom. This function was originally designed for the * implementation of local projections, for instance the Clement interpolant, * in conjunction with other local patch functions like * GridTools::build_triangulation_from_patch. * - * Containers built on top of Triangulation or - * parallel:distributed::Triangulation are supported and handled + * DoFHandlerType's built on top of Triangulation or + * parallel:distributed::Triangulation are supported and handled * appropriately. * - * It is necessary that the finite element underlying the MeshType used has - * degrees of freedom that are logically associated to a vertex, line, quad, - * or hex. This includes both nodal degrees of freedom and other modal types - * of dofs that are associated to an edge, etc. The result is the patch of - * cells representing the support of the basis element associated to the - * degree of freedom. For instance using an FE_Q finite element, we obtain - * the standard patch of cells touching the degree of freedom and then others - * that take care of possible hanging node constraints. Using a FE_DGQ finite - * element, the degrees of freedom are logically considered to be "interior" to - * the cells so the patch would be the cell on which the degree of freedom is + * The result is the patch of cells representing the support of the basis + * element associated to the degree of freedom. For instance using an FE_Q + * finite element, we obtain the standard patch of cells touching the degree + * of freedom and then add other cells that take care of possible hanging node + * constraints. Using a FE_DGQ finite element, the degrees of freedom are + * logically considered to be "interior" to the cells so the patch would + * consist exclusively of the single cell on which the degree of freedom is * located. * - * @tparam MeshType The MeshType should be a DoFHandler or hp::DoFHandler. - * @param[in] dof_handler The MeshType which could be built on a Triangulation - * or a parallel::distributed::Triangulation and should be using a finite - * element that has degrees of freedom that are logically associated to a vertex, - * line, quad, or hex. - * @param[out] dof_to_support_patch_map A map from the global_dof_index of DoFs - * on locally relevant cells to vectors containing MeshType::active_cell_iterators of - * cells in support of basis function at that degree of freedom. + * @tparam DoFHandlerType The DoFHandlerType should be a DoFHandler or + * hp::DoFHandler. + * @param[in] dof_handler The DoFHandlerType which could be built on a + * Triangulation or a parallel::distributed::Triangulation with a finite + * element that has degrees of freedom that are logically associated to a + * vertex, line, quad, or hex. + * @param[out] dof_to_support_patch_map A map from the global_dof_index of + * degrees of freedom on locally relevant cells to vectors containing + * DoFHandlerType::active_cell_iterators of cells in the support of the basis + * function at that degree of freedom. * * @author Spencer Patty, 2016 * */ - template - std::map< types::global_dof_index,std::vector > - get_dof_to_support_patch_map(MeshType &dof_handler); + template + std::map< types::global_dof_index,std::vector > + get_dof_to_support_patch_map(DoFHandlerType &dof_handler); /*@}*/ diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index 336c66159d..eb61d2164c 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -3049,9 +3049,9 @@ next_cell: - template - std::map< types::global_dof_index,std::vector > - get_dof_to_support_patch_map(MeshType &dof_handler) + template + std::map< types::global_dof_index,std::vector > + get_dof_to_support_patch_map(DoFHandlerType &dof_handler) { // This is the map from global_dof_index to @@ -3064,7 +3064,7 @@ next_cell: // constructed, we will copy to a map of vectors // since that is the prefered output for other // functions. - std::map< types::global_dof_index,std::set > dof_to_set_of_cells_map; + std::map< types::global_dof_index,std::set > dof_to_set_of_cells_map; std::vector local_dof_indices; std::vector local_face_dof_indices; @@ -3077,17 +3077,17 @@ next_cell: // in 3d, we need pointers from active lines to the // active parent lines, so we construct it as needed. - std::map lines_to_parent_lines_map; - if (MeshType::dimension == 3) + std::map lines_to_parent_lines_map; + if (DoFHandlerType::dimension == 3) { // save user flags as they will be modified and then later restored dof_handler.get_triangulation().save_user_flags(user_flags); - const_cast &>(dof_handler.get_triangulation()).clear_user_flags (); + const_cast &>(dof_handler.get_triangulation()).clear_user_flags (); - typename MeshType::active_cell_iterator cell = dof_handler.begin_active(), - endc = dof_handler.end(); + typename DoFHandlerType::active_cell_iterator cell = dof_handler.begin_active(), + endc = dof_handler.end(); for (; cell!=endc; ++cell) { // We only want lines that are locally_relevant @@ -3096,7 +3096,7 @@ next_cell: // few and we don't have to use them. if (cell->is_artificial() == false) { - for (unsigned int l=0; l::lines_per_cell; ++l) + for (unsigned int l=0; l::lines_per_cell; ++l) if (cell->line(l)->has_children()) for (unsigned int c=0; cline(l)->n_children(); ++c) { @@ -3116,8 +3116,8 @@ next_cell: // which it is a part, mainly the ones that must // be added on account of adaptivity hanging node // constraints. - typename MeshType::active_cell_iterator cell = dof_handler.begin_active(), - endc = dof_handler.end(); + typename DoFHandlerType::active_cell_iterator cell = dof_handler.begin_active(), + endc = dof_handler.end(); for (; cell!=endc; ++cell) { // Need to loop through all cells that could @@ -3142,7 +3142,7 @@ next_cell: // face (or line). // Take care of dofs on neighbor faces - for (unsigned int f=0; f::faces_per_cell; ++f) + for (unsigned int f=0; f::faces_per_cell; ++f) { if (cell->face(f)->has_children()) { @@ -3225,9 +3225,9 @@ next_cell: // if cell's line has an active parent, then // distribute cell to dofs on parent line // and dofs on all children of parent line. - if (MeshType::dimension == 3) + if (DoFHandlerType::dimension == 3) { - for (unsigned int l=0; l::lines_per_cell; ++l) + for (unsigned int l=0; l::lines_per_cell; ++l) { if (cell->line(l)->has_children()) { @@ -3252,7 +3252,7 @@ next_cell: // children else if (cell->line(l)->user_flag_set() == true) { - typename MeshType::line_iterator parent_line = lines_to_parent_lines_map[cell->line(l)]; + typename DoFHandlerType::line_iterator parent_line = lines_to_parent_lines_map[cell->line(l)]; Assert (parent_line->has_children(), ExcInternalError() ); // dofs_per_line returns number of dofs @@ -3281,24 +3281,24 @@ next_cell: } } // for lines l - }// if MeshType::dimension == 3 + }// if DoFHandlerType::dimension == 3 }// if cell->is_locally_owned() }// for cells - if (MeshType::dimension == 3) + if (DoFHandlerType::dimension == 3) { // finally, restore user flags that were changed above // to when we constructed the pointers to parent of lines // Since dof_handler is const, we must leave it unchanged. - const_cast &>(dof_handler.get_triangulation()).load_user_flags (user_flags); + const_cast &>(dof_handler.get_triangulation()).load_user_flags (user_flags); } // Finally, we copy map of sets to - // map of vectors using assign() - std::map< types::global_dof_index, std::vector > dof_to_cell_patches; + // map of vectors using the std::vector::assign() function + std::map< types::global_dof_index, std::vector > dof_to_cell_patches; - typename std::map >::iterator + typename std::map >::iterator it = dof_to_set_of_cells_map.begin(), it_end = dof_to_set_of_cells_map.end(); for ( ; it!=it_end; ++it) -- 2.39.5