]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Combine the vertex_index() functions.
authorDavid Wells <drwells@email.unc.edu>
Sat, 4 Jan 2025 18:39:29 +0000 (13:39 -0500)
committerDavid Wells <drwells@email.unc.edu>
Sat, 4 Jan 2025 18:58:30 +0000 (13:58 -0500)
include/deal.II/grid/tria_accessor.h
source/grid/tria.cc

index 51c37a64593ec0021a021caa6f355409381b69d2..0906841f71dd28d4cf939dea065f4a8ad332e5da 100644 (file)
@@ -4964,61 +4964,6 @@ namespace internal
 
 
 
-      /**
-       * Implementation of the function of same name in the enclosing class.
-       */
-      template <int dim, int spacedim>
-      inline static unsigned int
-      vertex_index(const TriaAccessor<1, dim, spacedim> &accessor,
-                   const unsigned int                    corner)
-      {
-        return accessor.objects()
-          .cells[accessor.present_index * ReferenceCell::max_n_faces<1>() +
-                 corner];
-      }
-
-
-      template <int dim, int spacedim>
-      inline static unsigned int
-      vertex_index(const TriaAccessor<2, dim, spacedim> &accessor,
-                   const unsigned int                    corner)
-      {
-        const auto [line_index, vertex_index] =
-          accessor.reference_cell().standard_vertex_to_face_and_vertex_index(
-            corner);
-        const auto vertex_within_line_index =
-          accessor.reference_cell().standard_to_real_face_vertex(
-            vertex_index,
-            line_index,
-            accessor.line_orientation(line_index) == true ?
-              ReferenceCell::default_combined_face_orientation() :
-              ReferenceCell::reversed_combined_line_orientation());
-
-        return accessor.line(line_index)
-          ->vertex_index(vertex_within_line_index);
-      }
-
-
-
-      inline static unsigned int
-      vertex_index(const TriaAccessor<3, 3, 3> &accessor,
-                   const unsigned int           corner)
-      {
-        const auto [face_index, vertex_index] =
-          accessor.reference_cell().standard_vertex_to_face_and_vertex_index(
-            corner);
-        const auto vertex_within_face_index =
-          accessor.reference_cell().standard_to_real_face_vertex(
-            vertex_index,
-            face_index,
-            accessor.combined_face_orientation(face_index));
-
-        return accessor.quad(face_index)
-          ->vertex_index(vertex_within_face_index);
-      }
-
-
-
       template <int dim, int spacedim>
       static std::array<unsigned int, 1>
       get_line_indices_of_cell(const TriaAccessor<1, dim, spacedim> &)
@@ -5354,8 +5299,18 @@ TriaAccessor<structdim, dim, spacedim>::vertex_index(
 {
   AssertIndexRange(corner, this->n_vertices());
 
-  if (structdim == dim)
+  if constexpr (structdim == 1)
     {
+      // This branch needs to be first (and not combined with the structdim ==
+      // dim branch) so that we can get line vertex indices when setting up the
+      // cell vertex index cache
+      return this->objects()
+        .cells[this->present_index * ReferenceCell::max_n_faces<1>() + corner];
+    }
+  else if constexpr (structdim == dim)
+    {
+      // This branch should only be used after the cell vertex index cache is
+      // set up
       constexpr unsigned int max_vertices_per_cell = 1 << dim;
       const std::size_t      my_index =
         static_cast<std::size_t>(this->present_index) * max_vertices_per_cell;
@@ -5368,9 +5323,25 @@ TriaAccessor<structdim, dim, spacedim>::vertex_index(
       Assert(vertex_index != numbers::invalid_unsigned_int, ExcInternalError());
       return vertex_index;
     }
+  else if constexpr (structdim == 2)
+    {
+      const auto [line_index, vertex_index] =
+        this->reference_cell().standard_vertex_to_face_and_vertex_index(corner);
+      const auto vertex_within_line_index =
+        this->reference_cell().standard_to_real_face_vertex(
+          vertex_index,
+          line_index,
+          this->line_orientation(line_index) == true ?
+            ReferenceCell::default_combined_face_orientation() :
+            ReferenceCell::reversed_combined_line_orientation());
+
+      return this->line(line_index)->vertex_index(vertex_within_line_index);
+    }
   else
-    return dealii::internal::TriaAccessorImplementation::Implementation::
-      vertex_index(*this, corner);
+    {
+      DEAL_II_ASSERT_UNREACHABLE();
+      return numbers::invalid_unsigned_int;
+    }
 }
 
 
@@ -7647,15 +7618,15 @@ CellAccessor<dim, spacedim>::face(const unsigned int i) const
   AssertIndexRange(i, this->n_faces());
   if constexpr (dim == 1)
     {
-      TriaAccessor<0, 1, spacedim> a(
-        &this->get_triangulation(),
-        ((i == 0) && at_boundary(0) ?
-           TriaAccessor<0, 1, spacedim>::left_vertex :
-           ((i == 1) && at_boundary(1) ?
-              TriaAccessor<0, 1, spacedim>::right_vertex :
-              TriaAccessor<0, 1, spacedim>::interior_vertex)),
-        internal::TriaAccessorImplementation::Implementation::vertex_index(
-          *this, i));
+      using VertexKind = typename TriaAccessor<0, 1, spacedim>::VertexKind;
+      VertexKind vertex_kind = VertexKind::interior_vertex;
+      if (i == 0 && at_boundary(0))
+        vertex_kind = VertexKind::left_vertex;
+      if (i == 1 && at_boundary(1))
+        vertex_kind = VertexKind::right_vertex;
+      TriaAccessor<0, 1, spacedim> a(&this->get_triangulation(),
+                                     vertex_kind,
+                                     this->vertex_index(i));
       return TriaIterator<TriaAccessor<0, 1, spacedim>>(a);
     }
   else if constexpr (dim == 2)
index 2427d4017c2bb266b597eba56575e420c3a4c6f3..90f9b684e3ac0975ed1271941392b740019355b7 100644 (file)
@@ -16011,10 +16011,28 @@ void Triangulation<dim, spacedim>::reset_cell_vertex_indices_cache()
               for (unsigned int i = 0; i < 4; ++i)
                 cache[my_index + i] = raw_vertex_indices[i];
             }
+          else if (ref_cell == ReferenceCells::Line)
+            {
+              cache[my_index + 0] = cell->vertex_index(0);
+              cache[my_index + 1] = cell->vertex_index(1);
+            }
           else
-            for (const unsigned int i : cell->vertex_indices())
-              cache[my_index + i] = internal::TriaAccessorImplementation::
-                Implementation::vertex_index(*cell, i);
+            {
+              Assert(dim == 2 || dim == 3, ExcInternalError());
+              for (const unsigned int i : cell->vertex_indices())
+                {
+                  const auto [face_index, vertex_index] =
+                    ref_cell.standard_vertex_to_face_and_vertex_index(i);
+                  const auto vertex_within_face_index =
+                    ref_cell.standard_to_real_face_vertex(
+                      vertex_index,
+                      face_index,
+                      cell->combined_face_orientation(face_index));
+                  cache[my_index + i] =
+                    cell->face(face_index)
+                      ->vertex_index(vertex_within_face_index);
+                }
+            }
         }
     }
 }

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.