From 714a73730896608362bb411237a7743dad40312f Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Tue, 10 Mar 2020 13:53:52 +0100 Subject: [PATCH] Simplify ParallelDistributed::distribute_dofs() by using get/set_dof_indices --- source/dofs/dof_handler_policy.cc | 384 +----------------------------- 1 file changed, 13 insertions(+), 371 deletions(-) diff --git a/source/dofs/dof_handler_policy.cc b/source/dofs/dof_handler_policy.cc index 39b4c5ce34..2494008672 100644 --- a/source/dofs/dof_handler_policy.cc +++ b/source/dofs/dof_handler_policy.cc @@ -219,373 +219,6 @@ namespace internal { /* -------------- distribute_dofs functionality ------------- */ - /** - * Distribute dofs on the given cell, with new dofs starting with index - * @p next_free_dof. Return the next unused index number. - * - * This function is refactored from the main @p distribute_dofs function since - * it can not be implemented dimension independent. - */ - template - static types::global_dof_index - distribute_dofs_on_cell( - const DoFHandler<1, spacedim> &dof_handler, - const typename DoFHandler<1, spacedim>::active_cell_iterator &cell, - types::global_dof_index next_free_dof) - { - // distribute dofs of vertices - if (dof_handler.get_fe().dofs_per_vertex > 0) - for (const unsigned int v : GeometryInfo<1>::vertex_indices()) - { - if (cell->vertex_dof_index(v, 0) == numbers::invalid_dof_index) - for (unsigned int d = 0; - d < dof_handler.get_fe().dofs_per_vertex; - ++d) - { - Assert((cell->vertex_dof_index(v, d) == - numbers::invalid_dof_index), - ExcInternalError()); - cell->set_vertex_dof_index(v, d, next_free_dof++); - } - else - for (unsigned int d = 0; - d < dof_handler.get_fe().dofs_per_vertex; - ++d) - Assert((cell->vertex_dof_index(v, d) != - numbers::invalid_dof_index), - ExcInternalError()); - } - - // dofs of line - for (unsigned int d = 0; d < dof_handler.get_fe().dofs_per_line; ++d) - cell->set_dof_index(d, next_free_dof++); - - return next_free_dof; - } - - - - template - static types::global_dof_index - distribute_dofs_on_cell( - const DoFHandler<2, spacedim> &dof_handler, - const typename DoFHandler<2, spacedim>::active_cell_iterator &cell, - types::global_dof_index next_free_dof) - { - if (dof_handler.get_fe().dofs_per_vertex > 0) - // number dofs on vertices - for (const unsigned int vertex : GeometryInfo<2>::vertex_indices()) - // check whether dofs for this vertex have been distributed - // (checking the first dof should be good enough) - if (cell->vertex_dof_index(vertex, 0) == - numbers::invalid_dof_index) - for (unsigned int d = 0; - d < dof_handler.get_fe().dofs_per_vertex; - ++d) - cell->set_vertex_dof_index(vertex, d, next_free_dof++); - - // for the four sides - if (dof_handler.get_fe().dofs_per_line > 0) - for (const unsigned int side : GeometryInfo<2>::face_indices()) - { - const typename DoFHandler<2, spacedim>::line_iterator line = - cell->line(side); - - // distribute dofs if necessary: check whether line dof is - // already numbered (checking the first dof should be good - // enough) - if (line->dof_index(0) == numbers::invalid_dof_index) - // if not: distribute dofs - for (unsigned int d = 0; - d < dof_handler.get_fe().dofs_per_line; - ++d) - line->set_dof_index(d, next_free_dof++); - } - - - // dofs of quad - if (dof_handler.get_fe().dofs_per_quad > 0) - for (unsigned int d = 0; d < dof_handler.get_fe().dofs_per_quad; - ++d) - cell->set_dof_index(d, next_free_dof++); - - return next_free_dof; - } - - - - template - static types::global_dof_index - distribute_dofs_on_cell( - const DoFHandler<3, spacedim> &dof_handler, - const typename DoFHandler<3, spacedim>::active_cell_iterator &cell, - types::global_dof_index next_free_dof) - { - if (dof_handler.get_fe().dofs_per_vertex > 0) - // number dofs on vertices - for (const unsigned int vertex : GeometryInfo<3>::vertex_indices()) - // check whether dofs for this vertex have been distributed - // (checking the first dof should be good enough) - if (cell->vertex_dof_index(vertex, 0) == - numbers::invalid_dof_index) - for (unsigned int d = 0; - d < dof_handler.get_fe().dofs_per_vertex; - ++d) - cell->set_vertex_dof_index(vertex, d, next_free_dof++); - - // for the lines - if (dof_handler.get_fe().dofs_per_line > 0) - for (unsigned int l = 0; l < GeometryInfo<3>::lines_per_cell; ++l) - { - const typename DoFHandler<3, spacedim>::line_iterator line = - cell->line(l); - - // distribute dofs if necessary: check whether line dof is - // already numbered (checking the first dof should be good - // enough) - if (line->dof_index(0) == numbers::invalid_dof_index) - // if not: distribute dofs - for (unsigned int d = 0; - d < dof_handler.get_fe().dofs_per_line; - ++d) - line->set_dof_index(d, next_free_dof++); - } - - // for the quads - if (dof_handler.get_fe().dofs_per_quad > 0) - for (unsigned int q = 0; q < GeometryInfo<3>::quads_per_cell; ++q) - { - const typename DoFHandler<3, spacedim>::quad_iterator quad = - cell->quad(q); - - // distribute dofs if necessary: check whether line dof is - // already numbered (checking the first dof should be good - // enough) - if (quad->dof_index(0) == numbers::invalid_dof_index) - // if not: distribute dofs - for (unsigned int d = 0; - d < dof_handler.get_fe().dofs_per_quad; - ++d) - quad->set_dof_index( - dof_handler.get_fe() - .adjust_quad_dof_index_for_face_orientation( - d, - cell->face_orientation(q), - cell->face_flip(q), - cell->face_rotation(q)), - next_free_dof++); - } - - - // dofs of hex - if (dof_handler.get_fe().dofs_per_hex > 0) - for (unsigned int d = 0; d < dof_handler.get_fe().dofs_per_hex; ++d) - cell->set_dof_index(d, next_free_dof++); - - return next_free_dof; - } - - - - // same for the hp::DoFHandler - template - static types::global_dof_index - distribute_dofs_on_cell( - const hp::DoFHandler<1, spacedim> &, - const typename hp::DoFHandler<1, spacedim>::active_cell_iterator - & cell, - types::global_dof_index next_free_dof) - { - const unsigned int dim = 1; - - const FiniteElement &fe = cell->get_fe(); - const unsigned int fe_index = cell->active_fe_index(); - - // number dofs on vertices. to do so, check whether dofs for - // this vertex have been distributed and for the present fe - // (only check the first dof), and if this isn't the case - // distribute new ones there - if (fe.dofs_per_vertex > 0) - for (const unsigned int vertex : GeometryInfo<1>::vertex_indices()) - if (cell->vertex_dof_index(vertex, 0, fe_index) == - numbers::invalid_dof_index) - for (unsigned int d = 0; d < fe.dofs_per_vertex; - ++d, ++next_free_dof) - cell->set_vertex_dof_index(vertex, - d, - next_free_dof, - fe_index); - - // finally for the line. this one shouldn't be numbered yet - if (fe.dofs_per_line > 0) - { - Assert((cell->dof_index(0, fe_index) == - numbers::invalid_dof_index), - ExcInternalError()); - - for (unsigned int d = 0; d < fe.dofs_per_line; - ++d, ++next_free_dof) - cell->set_dof_index(d, next_free_dof, fe_index); - } - - // note that this cell has been processed - cell->set_user_flag(); - - return next_free_dof; - } - - - - template - static types::global_dof_index - distribute_dofs_on_cell( - const hp::DoFHandler<2, spacedim> &, - const typename hp::DoFHandler<2, spacedim>::active_cell_iterator - & cell, - types::global_dof_index next_free_dof) - { - const unsigned int dim = 2; - - const FiniteElement &fe = cell->get_fe(); - const unsigned int fe_index = cell->active_fe_index(); - - // number dofs on vertices. to do so, check whether dofs for - // this vertex have been distributed and for the present fe - // (only check the first dof), and if this isn't the case - // distribute new ones there - if (fe.dofs_per_vertex > 0) - for (const unsigned int vertex : GeometryInfo<2>::vertex_indices()) - if (cell->vertex_dof_index(vertex, 0, fe_index) == - numbers::invalid_dof_index) - for (unsigned int d = 0; d < fe.dofs_per_vertex; - ++d, ++next_free_dof) - cell->set_vertex_dof_index(vertex, - d, - next_free_dof, - fe_index); - - // next the sides. do the same as above: check whether the - // line is already numbered for the present fe_index, and if - // not do it - if (fe.dofs_per_line > 0) - for (unsigned int l = 0; l < GeometryInfo<2>::lines_per_cell; ++l) - { - const typename hp::DoFHandler::line_iterator - line = cell->line(l); - - if (line->dof_index(0, fe_index) == numbers::invalid_dof_index) - for (unsigned int d = 0; d < fe.dofs_per_line; - ++d, ++next_free_dof) - line->set_dof_index(d, next_free_dof, fe_index); - } - - - // finally for the quad. this one shouldn't be numbered yet - if (fe.dofs_per_quad > 0) - { - Assert((cell->dof_index(0, fe_index) == - numbers::invalid_dof_index), - ExcInternalError()); - - for (unsigned int d = 0; d < fe.dofs_per_quad; - ++d, ++next_free_dof) - cell->set_dof_index(d, next_free_dof, fe_index); - } - - // note that this cell has been processed - cell->set_user_flag(); - - return next_free_dof; - } - - - - template - static types::global_dof_index - distribute_dofs_on_cell( - const hp::DoFHandler<3, spacedim> &, - const typename hp::DoFHandler<3, spacedim>::active_cell_iterator - & cell, - types::global_dof_index next_free_dof) - { - const unsigned int dim = 3; - - const FiniteElement &fe = cell->get_fe(); - const unsigned int fe_index = cell->active_fe_index(); - - // number dofs on vertices. to do so, check whether dofs for - // this vertex have been distributed and for the present fe - // (only check the first dof), and if this isn't the case - // distribute new ones there - if (fe.dofs_per_vertex > 0) - for (const unsigned int vertex : GeometryInfo<3>::vertex_indices()) - if (cell->vertex_dof_index(vertex, 0, fe_index) == - numbers::invalid_dof_index) - for (unsigned int d = 0; d < fe.dofs_per_vertex; - ++d, ++next_free_dof) - cell->set_vertex_dof_index(vertex, - d, - next_free_dof, - fe_index); - - // next the eight lines. do the same as above: check whether - // the line is already numbered for the present fe_index, - // and if not do it - if (fe.dofs_per_line > 0) - for (unsigned int l = 0; l < GeometryInfo<3>::lines_per_cell; ++l) - { - const typename hp::DoFHandler::line_iterator - line = cell->line(l); - - if (line->dof_index(0, fe_index) == numbers::invalid_dof_index) - for (unsigned int d = 0; d < fe.dofs_per_line; - ++d, ++next_free_dof) - line->set_dof_index(d, next_free_dof, fe_index); - } - - // same for quads - if (fe.dofs_per_quad > 0) - for (unsigned int q = 0; q < GeometryInfo<3>::quads_per_cell; ++q) - { - const typename hp::DoFHandler::quad_iterator - quad = cell->quad(q); - - if (quad->dof_index(0, fe_index) == numbers::invalid_dof_index) - for (unsigned int d = 0; d < fe.dofs_per_quad; - ++d, ++next_free_dof) - quad->set_dof_index( - fe.adjust_quad_dof_index_for_face_orientation( - d, - cell->face_orientation(q), - cell->face_flip(q), - cell->face_rotation(q)), - next_free_dof, - fe_index); - } - - - // finally for the hex. this one shouldn't be numbered yet - // because there is no other cell from which we could already - // have gotten to numbering anything on this hex (=cell) - if (fe.dofs_per_hex > 0) - { - Assert((cell->dof_index(0, fe_index) == - numbers::invalid_dof_index), - ExcInternalError()); - - for (unsigned int d = 0; d < fe.dofs_per_hex; - ++d, ++next_free_dof) - cell->set_dof_index(d, next_free_dof, fe_index); - } - - // note that this cell has been processed - cell->set_user_flag(); - - return next_free_dof; - } - - - /** * Compute identities between DoFs located on vertices. Called from * distribute_dofs(). @@ -2034,10 +1667,19 @@ namespace internal if (!cell->is_artificial()) if ((subdomain_id == numbers::invalid_subdomain_id) || (cell->subdomain_id() == subdomain_id)) - next_free_dof = - Implementation::distribute_dofs_on_cell(dof_handler, - cell, - next_free_dof); + { + std::vector dof_indices( + cell->get_fe().dofs_per_cell); + + internal::DoFAccessorImplementation::get_dof_indices( + *cell, dof_indices, cell->active_fe_index()); + + for (auto &dof_index : dof_indices) + if (dof_index == numbers::invalid_dof_index) + dof_index = next_free_dof++; + + cell->set_dof_indices(dof_indices); + } update_all_active_cell_dof_indices_caches(dof_handler); -- 2.39.5