Simplify ParallelDistributed::distribute_dofs() by using get/set_dof_indices

diff --git a/source/dofs/dof_handler_policy.cc b/source/dofs/dof_handler_policy.cc
index 39b4c5ce3416da1a075962c1b7ad2a2a1c8281e7..24940086725beecd8befc2589e75d697c28da6f7 100644 (file)
--- 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 <int spacedim>
-        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 <int spacedim>
-        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 <int spacedim>
-        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 <int spacedim>
-        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<dim, spacedim> &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 <int spacedim>
-        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<dim, spacedim> &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<dim, spacedim>::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 <int spacedim>
-        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<dim, spacedim> &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<dim, spacedim>::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<dim, spacedim>::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<types::global_dof_index> 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);