From: David Wells Date: Mon, 9 Sep 2024 20:08:38 +0000 (-0400) Subject: TriaAccessor::face(): consolidate implementations via if constexpr. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F17687%2Fhead;p=dealii.git TriaAccessor::face(): consolidate implementations via if constexpr. --- diff --git a/include/deal.II/grid/tria_accessor.templates.h b/include/deal.II/grid/tria_accessor.templates.h index fb0f3cadbf..2d959a6570 100644 --- a/include/deal.II/grid/tria_accessor.templates.h +++ b/include/deal.II/grid/tria_accessor.templates.h @@ -3342,49 +3342,6 @@ inline CellAccessor::CellAccessor( -namespace internal -{ - namespace CellAccessorImplementation - { - template - inline dealii::TriaIterator> - get_face(const dealii::CellAccessor<1, spacedim> &cell, - const unsigned int i) - { - dealii::TriaAccessor<0, 1, spacedim> a( - &cell.get_triangulation(), - ((i == 0) && cell.at_boundary(0) ? - dealii::TriaAccessor<0, 1, spacedim>::left_vertex : - ((i == 1) && cell.at_boundary(1) ? - dealii::TriaAccessor<0, 1, spacedim>::right_vertex : - dealii::TriaAccessor<0, 1, spacedim>::interior_vertex)), - dealii::internal::TriaAccessorImplementation::Implementation:: - vertex_index(cell, i)); - return dealii::TriaIterator>(a); - } - - - template - inline dealii::TriaIterator> - get_face(const dealii::CellAccessor<2, spacedim> &cell, - const unsigned int i) - { - return cell.line(i); - } - - - template - inline dealii::TriaIterator> - get_face(const dealii::CellAccessor<3, spacedim> &cell, - const unsigned int i) - { - return cell.quad(i); - } - } // namespace CellAccessorImplementation -} // namespace internal - - - template inline TriaIterator> CellAccessor::child(const unsigned int i) const @@ -3423,7 +3380,28 @@ inline TriaIterator> CellAccessor::face(const unsigned int i) const { AssertIndexRange(i, this->n_faces()); - return dealii::internal::CellAccessorImplementation::get_face(*this, i); + 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)); + return TriaIterator>(a); + } + else if constexpr (dim == 2) + return this->line(i); + else if constexpr (dim == 3) + return this->quad(i); + else + { + Assert(false, ExcNotImplemented()); + return {}; + } } @@ -3456,8 +3434,7 @@ CellAccessor::face_iterators() const face_iterators(this->n_faces()); for (const unsigned int i : this->face_indices()) - face_iterators[i] = - dealii::internal::CellAccessorImplementation::get_face(*this, i); + face_iterators[i] = this->face(i); return face_iterators; }