From 366468e305d3ad8a6b6ad7968ca71096ab6b5717 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Sun, 2 Jul 2023 23:32:58 +0200 Subject: [PATCH] Introduce tria_orientation.h --- include/deal.II/grid/reference_cell.h | 137 ++++++++-------- .../deal.II/grid/tria_accessor.templates.h | 31 +++- include/deal.II/grid/tria_orientation.h | 60 +++++++ source/base/qprojector.cc | 16 +- source/fe/fe.cc | 5 +- source/fe/fe_poly_tensor.cc | 6 +- source/fe/fe_q_base.cc | 32 ++-- source/fe/fe_raviart_thomas.cc | 146 ++++++++++-------- source/fe/fe_raviart_thomas_nodal.cc | 56 +++---- 9 files changed, 293 insertions(+), 196 deletions(-) create mode 100644 include/deal.II/grid/tria_orientation.h diff --git a/include/deal.II/grid/reference_cell.h b/include/deal.II/grid/reference_cell.h index 65c962fa11..a1700420a5 100644 --- a/include/deal.II/grid/reference_cell.h +++ b/include/deal.II/grid/reference_cell.h @@ -25,6 +25,8 @@ #include #include +#include + #include #include @@ -1521,12 +1523,8 @@ ReferenceCell::child_cell_on_face( } case ReferenceCells::Quadrilateral: { - const bool face_orientation = - Utilities::get_bit(combined_face_orientation, 0); - const bool face_flip = - Utilities::get_bit(combined_face_orientation, 2); - const bool face_rotation = - Utilities::get_bit(combined_face_orientation, 1); + const auto [face_orientation, face_rotation, face_flip] = + internal::split_face_orientation(combined_face_orientation); return GeometryInfo<2>::child_cell_on_face( RefinementCase<2>(RefinementPossibilities<2>::isotropic_refinement), @@ -1545,12 +1543,8 @@ ReferenceCell::child_cell_on_face( } case ReferenceCells::Hexahedron: { - const bool face_orientation = - Utilities::get_bit(combined_face_orientation, 0); - const bool face_flip = - Utilities::get_bit(combined_face_orientation, 2); - const bool face_rotation = - Utilities::get_bit(combined_face_orientation, 1); + const auto [face_orientation, face_rotation, face_flip] = + internal::split_face_orientation(combined_face_orientation); return GeometryInfo<3>::child_cell_on_face( RefinementCase<3>(RefinementPossibilities<3>::isotropic_refinement), @@ -1776,9 +1770,10 @@ ReferenceCell::line_to_cell_vertices(const unsigned int line, inline unsigned int -ReferenceCell::face_to_cell_lines(const unsigned int face, - const unsigned int line, - const unsigned char face_orientation) const +ReferenceCell::face_to_cell_lines( + const unsigned int face, + const unsigned int line, + const unsigned char combined_face_orientation) const { AssertIndexRange(face, n_faces()); AssertIndexRange(line, face_reference_cell(face).n_lines()); @@ -1794,12 +1789,11 @@ ReferenceCell::face_to_cell_lines(const unsigned int face, } case ReferenceCells::Line: { + const auto [face_orientation, face_rotation, face_flip] = + internal::split_face_orientation(combined_face_orientation); + return GeometryInfo<1>::face_to_cell_lines( - face, - line, - Utilities::get_bit(face_orientation, 0), - Utilities::get_bit(face_orientation, 2), - Utilities::get_bit(face_orientation, 1)); + face, line, face_orientation, face_flip, face_rotation); } case ReferenceCells::Triangle: { @@ -1807,12 +1801,11 @@ ReferenceCell::face_to_cell_lines(const unsigned int face, } case ReferenceCells::Quadrilateral: { + const auto [face_orientation, face_rotation, face_flip] = + internal::split_face_orientation(combined_face_orientation); + return GeometryInfo<2>::face_to_cell_lines( - face, - line, - Utilities::get_bit(face_orientation, 0), - Utilities::get_bit(face_orientation, 2), - Utilities::get_bit(face_orientation, 1)); + face, line, face_orientation, face_flip, face_rotation); } case ReferenceCells::Tetrahedron: { @@ -1820,7 +1813,7 @@ ReferenceCell::face_to_cell_lines(const unsigned int face, {{{0, 1, 2}}, {{0, 3, 4}}, {{2, 5, 3}}, {{1, 4, 5}}}}; return table[face][standard_to_real_face_line( - line, face, face_orientation)]; + line, face, combined_face_orientation)]; } case ReferenceCells::Pyramid: { @@ -1832,7 +1825,7 @@ ReferenceCell::face_to_cell_lines(const unsigned int face, {{3, 7, 6, X}}}}; return table[face][standard_to_real_face_line( - line, face, face_orientation)]; + line, face, combined_face_orientation)]; } case ReferenceCells::Wedge: { @@ -1844,16 +1837,15 @@ ReferenceCell::face_to_cell_lines(const unsigned int face, {{8, 6, 5, 2}}}}; return table[face][standard_to_real_face_line( - line, face, face_orientation)]; + line, face, combined_face_orientation)]; } case ReferenceCells::Hexahedron: { + const auto [face_orientation, face_rotation, face_flip] = + internal::split_face_orientation(combined_face_orientation); + return GeometryInfo<3>::face_to_cell_lines( - face, - line, - Utilities::get_bit(face_orientation, 0), - Utilities::get_bit(face_orientation, 2), - Utilities::get_bit(face_orientation, 1)); + face, line, face_orientation, face_flip, face_rotation); } default: Assert(false, ExcNotImplemented()); @@ -1865,9 +1857,10 @@ ReferenceCell::face_to_cell_lines(const unsigned int face, inline unsigned int -ReferenceCell::face_to_cell_vertices(const unsigned int face, - const unsigned int vertex, - const unsigned char face_orientation) const +ReferenceCell::face_to_cell_vertices( + const unsigned int face, + const unsigned int vertex, + const unsigned char combined_face_orientation) const { AssertIndexRange(face, n_faces()); AssertIndexRange(vertex, face_reference_cell(face).n_vertices()); @@ -1881,28 +1874,29 @@ ReferenceCell::face_to_cell_vertices(const unsigned int face, } case ReferenceCells::Line: { + const auto [face_orientation, face_rotation, face_flip] = + internal::split_face_orientation(combined_face_orientation); + return GeometryInfo<1>::face_to_cell_vertices( - face, - vertex, - Utilities::get_bit(face_orientation, 0), - Utilities::get_bit(face_orientation, 2), - Utilities::get_bit(face_orientation, 1)); + face, vertex, face_orientation, face_flip, face_rotation); } case ReferenceCells::Triangle: { static constexpr ndarray table = { {{{0, 1}}, {{1, 2}}, {{2, 0}}}}; - return table[face][face_orientation != 0u ? vertex : (1 - vertex)]; + return table[face][combined_face_orientation != + reversed_combined_line_orientation() ? + vertex : + (1 - vertex)]; } case ReferenceCells::Quadrilateral: { + const auto [face_orientation, face_rotation, face_flip] = + internal::split_face_orientation(combined_face_orientation); + return GeometryInfo<2>::face_to_cell_vertices( - face, - vertex, - Utilities::get_bit(face_orientation, 0), - Utilities::get_bit(face_orientation, 2), - Utilities::get_bit(face_orientation, 1)); + face, vertex, face_orientation, face_flip, face_rotation); } case ReferenceCells::Tetrahedron: { @@ -1910,7 +1904,7 @@ ReferenceCell::face_to_cell_vertices(const unsigned int face, {{{0, 1, 2}}, {{1, 0, 3}}, {{0, 2, 3}}, {{2, 1, 3}}}}; return table[face][standard_to_real_face_vertex( - vertex, face, face_orientation)]; + vertex, face, combined_face_orientation)]; } case ReferenceCells::Pyramid: { @@ -1923,7 +1917,7 @@ ReferenceCell::face_to_cell_vertices(const unsigned int face, {{2, 3, 4, X}}}}; return table[face][standard_to_real_face_vertex( - vertex, face, face_orientation)]; + vertex, face, combined_face_orientation)]; } case ReferenceCells::Wedge: { @@ -1936,16 +1930,15 @@ ReferenceCell::face_to_cell_vertices(const unsigned int face, {{2, 0, 5, 3}}}}; return table[face][standard_to_real_face_vertex( - vertex, face, face_orientation)]; + vertex, face, combined_face_orientation)]; } case ReferenceCells::Hexahedron: { + const auto [face_orientation, face_rotation, face_flip] = + internal::split_face_orientation(combined_face_orientation); + return GeometryInfo<3>::face_to_cell_vertices( - face, - vertex, - Utilities::get_bit(face_orientation, 0), - Utilities::get_bit(face_orientation, 2), - Utilities::get_bit(face_orientation, 1)); + face, vertex, face_orientation, face_flip, face_rotation); } default: Assert(false, ExcNotImplemented()); @@ -2015,7 +2008,7 @@ inline unsigned int ReferenceCell::standard_to_real_face_line( const unsigned int line, const unsigned int face, - const unsigned char face_orientation) const + const unsigned char combined_face_orientation) const { AssertIndexRange(face, n_faces()); AssertIndexRange(line, face_reference_cell(face).n_lines()); @@ -2037,31 +2030,35 @@ ReferenceCell::standard_to_real_face_line( Assert(false, ExcNotImplemented()); break; case ReferenceCells::Tetrahedron: - return triangle_table[face_orientation][line]; + return triangle_table[combined_face_orientation][line]; case ReferenceCells::Pyramid: if (face == 0) // The quadrilateral face { - return GeometryInfo<3>::standard_to_real_face_line( - line, - Utilities::get_bit(face_orientation, 0), - Utilities::get_bit(face_orientation, 2), - Utilities::get_bit(face_orientation, 1)); + const auto [face_orientation, face_rotation, face_flip] = + internal::split_face_orientation(combined_face_orientation); + + return GeometryInfo<3>::standard_to_real_face_line(line, + face_orientation, + face_flip, + face_rotation); } else // One of the triangular faces { - return triangle_table[face_orientation][line]; + return triangle_table[combined_face_orientation][line]; } case ReferenceCells::Wedge: if (face > 1) // One of the quadrilateral faces { - return GeometryInfo<3>::standard_to_real_face_line( - line, - Utilities::get_bit(face_orientation, 0), - Utilities::get_bit(face_orientation, 2), - Utilities::get_bit(face_orientation, 1)); + const auto [face_orientation, face_rotation, face_flip] = + internal::split_face_orientation(combined_face_orientation); + + return GeometryInfo<3>::standard_to_real_face_line(line, + face_orientation, + face_flip, + face_rotation); } else // One of the triangular faces - return triangle_table[face_orientation][line]; + return triangle_table[combined_face_orientation][line]; case ReferenceCells::Hexahedron: { static constexpr ndarray table = { @@ -2073,7 +2070,7 @@ ReferenceCell::standard_to_real_face_line( {{1, 0, 3, 2}}, {{1, 0, 2, 3}}, {{2, 3, 1, 0}}}}; - return table[face_orientation][line]; + return table[combined_face_orientation][line]; } default: Assert(false, ExcNotImplemented()); diff --git a/include/deal.II/grid/tria_accessor.templates.h b/include/deal.II/grid/tria_accessor.templates.h index d063945984..9b079143c2 100644 --- a/include/deal.II/grid/tria_accessor.templates.h +++ b/include/deal.II/grid/tria_accessor.templates.h @@ -688,7 +688,32 @@ namespace internal } - inline static unsigned int + + template + inline static unsigned char + combined_face_orientation( + const TriaAccessor<1, dim, spacedim> & /*accessor*/, + const unsigned int /*face*/) + { + // There is only one way to orient a vertex + return ReferenceCell::default_combined_face_orientation(); + } + + + + template + inline static unsigned char + combined_face_orientation(const TriaAccessor<2, dim, spacedim> &accessor, + const unsigned int face) + { + return line_orientation(accessor, face) == true ? + ReferenceCell::default_combined_face_orientation() : + ReferenceCell::reversed_combined_line_orientation(); + } + + + + inline static unsigned char combined_face_orientation(const TriaAccessor<3, 3, 3> &accessor, const unsigned int face) { @@ -1363,8 +1388,8 @@ inline unsigned char TriaAccessor::combined_face_orientation( const unsigned int face) const { - return this->face_orientation(face) + 4 * this->face_flip(face) + - 2 * this->face_rotation(face); + return dealii::internal::TriaAccessorImplementation::Implementation:: + combined_face_orientation(*this, face); } diff --git a/include/deal.II/grid/tria_orientation.h b/include/deal.II/grid/tria_orientation.h new file mode 100644 index 0000000000..3047abb93d --- /dev/null +++ b/include/deal.II/grid/tria_orientation.h @@ -0,0 +1,60 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +#ifndef dealii_tria_orientation_h +#define dealii_tria_orientation_h + +#include + +#include + +#include + +DEAL_II_NAMESPACE_OPEN + +namespace internal +{ + /** + * Combine orientation flags. + */ + inline unsigned char + combined_face_orientation(const bool face_orientation, + const bool face_rotation, + const bool face_flip) + { + return (face_orientation ? 1 : 0) + (face_rotation ? 2 : 0) + + (face_flip ? 4 : 0); + } + + /** + * Split up a combined orientation flag: orientation flag, + * rotation flag, flip flag. + */ + inline std::tuple + split_face_orientation(const unsigned char combined_face_orientation) + { + return {Utilities::get_bit(combined_face_orientation, 0), + Utilities::get_bit(combined_face_orientation, 1), + Utilities::get_bit(combined_face_orientation, 2)}; + } + + +} // namespace internal + + + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/source/base/qprojector.cc b/source/base/qprojector.cc index 8d3a1261bd..8c397d507c 100644 --- a/source/base/qprojector.cc +++ b/source/base/qprojector.cc @@ -20,6 +20,7 @@ #include #include +#include DEAL_II_NAMESPACE_OPEN @@ -1330,9 +1331,10 @@ QProjector::DataSetDescriptor::face(const ReferenceCell &reference_cell, n_quadrature_points}; else if (dim == 3) { - const unsigned char orientation = (face_flip ? 4 : 0) + - (face_rotation ? 2 : 0) + - (face_orientation ? 1 : 0); + const unsigned char orientation = + internal::combined_face_orientation(face_orientation, + face_rotation, + face_flip); Assert(orientation < 6, ExcInternalError()); return {(reference_cell.n_face_orientations(face_no) * face_no + orientation) * @@ -1449,12 +1451,10 @@ QProjector::DataSetDescriptor::face( quadrature[quadrature.size() == 1 ? 0 : face_no].size()}; else if (dim == 3) { - const unsigned int orientation = (face_flip ? 4 : 0) + - (face_rotation ? 2 : 0) + - (face_orientation ? 1 : 0); - return {offset + - orientation * + internal::combined_face_orientation(face_orientation, + face_rotation, + face_flip) * quadrature[quadrature.size() == 1 ? 0 : face_no].size()}; } } diff --git a/source/fe/fe.cc b/source/fe/fe.cc index 06c10664de..d956a51fc5 100644 --- a/source/fe/fe.cc +++ b/source/fe/fe.cc @@ -684,8 +684,9 @@ FiniteElement::adjust_quad_dof_index_for_face_orientation( ExcInternalError()); return index + adjust_quad_dof_index_for_face_orientation_table[table_n]( index, - (face_orientation ? 4 : 0) + (face_flip ? 2 : 0) + - (face_rotation ? 1 : 0)); + internal::combined_face_orientation(face_orientation, + face_rotation, + face_flip)); } diff --git a/source/fe/fe_poly_tensor.cc b/source/fe/fe_poly_tensor.cc index 714642ac8a..a03ad8f818 100644 --- a/source/fe/fe_poly_tensor.cc +++ b/source/fe/fe_poly_tensor.cc @@ -30,6 +30,7 @@ #include #include +#include DEAL_II_NAMESPACE_OPEN @@ -264,8 +265,9 @@ FE_PolyTensor::adjust_quad_dof_sign_for_face_orientation( return adjust_quad_dof_sign_for_face_orientation_table [this->n_unique_2d_subobjects() == 1 ? 0 : face]( index, - 4 * static_cast(face_orientation) + 2 * static_cast(face_flip) + - static_cast(face_rotation)); + internal::combined_face_orientation(face_orientation, + face_rotation, + face_flip)); } diff --git a/source/fe/fe_q_base.cc b/source/fe/fe_q_base.cc index 68bf4dda83..e8800073cc 100644 --- a/source/fe/fe_q_base.cc +++ b/source/fe/fe_q_base.cc @@ -1061,35 +1061,35 @@ FE_Q_Base::initialize_quad_dof_index_permutation() unsigned int i = local % n, j = local / n; // face_orientation=false, face_flip=false, face_rotation=false - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 0) = + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(false, false, false)) = j + i * n - local; // face_orientation=false, face_flip=false, face_rotation=true - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 1) = + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(false, true, false)) = i + (n - 1 - j) * n - local; // face_orientation=false, face_flip=true, face_rotation=false - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 2) = + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(false, false, true)) = (n - 1 - j) + (n - 1 - i) * n - local; // face_orientation=false, face_flip=true, face_rotation=true - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 3) = + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(false, true, true)) = (n - 1 - i) + j * n - local; // face_orientation=true, face_flip=false, face_rotation=false - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 4) = 0; + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(true, false, false)) = 0; // face_orientation=true, face_flip=false, face_rotation=true - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 5) = + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(true, true, false)) = j + (n - 1 - i) * n - local; // face_orientation=true, face_flip=true, face_rotation=false - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 6) = + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(true, false, true)) = (n - 1 - i) + (n - 1 - j) * n - local; // face_orientation=true, face_flip=true, face_rotation=true - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 7) = + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(true, true, true)) = (n - 1 - j) + i * n - local; } diff --git a/source/fe/fe_raviart_thomas.cc b/source/fe/fe_raviart_thomas.cc index 442438fb59..b7c00b5bb9 100644 --- a/source/fe/fe_raviart_thomas.cc +++ b/source/fe/fe_raviart_thomas.cc @@ -341,75 +341,85 @@ FE_RaviartThomas::initialize_quad_dof_index_permutation_and_sign_change() // We have 8 cases that are all treated the same way. Note that the // corresponding case to case_no is just its binary representation. // The above example of (false | true | true) would be case_no=3 - for (unsigned int case_no = 0; case_no < 8; ++case_no) - { - // Get the binary representation of the case - const bool face_orientation = Utilities::get_bit(case_no, 2); - const bool face_flip = Utilities::get_bit(case_no, 1); - const bool face_rotation = Utilities::get_bit(case_no, 0); - - if (((!face_orientation) && (!face_rotation)) || - ((face_orientation) && (face_rotation))) - { - // We flip across the diagonal - // This is the local face dof offset - this->adjust_quad_dof_index_for_face_orientation_table[face_no]( - local_face_dof, case_no) = j + i * n - local_face_dof; - } - else - { - // Offset is zero - this->adjust_quad_dof_index_for_face_orientation_table[face_no]( - local_face_dof, case_no) = 0; - } // if face needs dof permutation - - // Get new local face_dof by adding offset - const unsigned int new_local_face_dof = - local_face_dof + - this->adjust_quad_dof_index_for_face_orientation_table[face_no]( - local_face_dof, case_no); - // compute new row and column index - i = new_local_face_dof % n; - j = new_local_face_dof / n; - - /* - * Now compute if a sign change is necessary. This is done for the - * case of face_orientation==true - */ - // flip = false, rotation=true - if (!face_flip && face_rotation) + for (const bool face_orientation : {false, true}) + for (const bool face_flip : {false, true}) + for (const bool face_rotation : {false, true}) { - this->adjust_quad_dof_sign_for_face_orientation_table[face_no]( - local_face_dof, case_no) = ((j % 2) == 1); - } - // flip = true, rotation=false - else if (face_flip && !face_rotation) - { - // This case is symmetric (although row and column may be - // switched) - this->adjust_quad_dof_sign_for_face_orientation_table[face_no]( - local_face_dof, case_no) = ((j % 2) == 1) != ((i % 2) == 1); - } - // flip = true, rotation=true - else if (face_flip && face_rotation) - { - this->adjust_quad_dof_sign_for_face_orientation_table[face_no]( - local_face_dof, case_no) = ((i % 2) == 1); - } - /* - * flip = false, rotation=false => nothing to do - */ - - /* - * If face_orientation==false the sign flip is exactly the opposite. - */ - if (!face_orientation) - this->adjust_quad_dof_sign_for_face_orientation_table[face_no]( - local_face_dof, case_no) = - !this->adjust_quad_dof_sign_for_face_orientation_table[face_no]( - local_face_dof, case_no); - } // case_no - } // local_face_dof + const auto case_no = + internal::combined_face_orientation(face_orientation, + face_rotation, + face_flip); + + if (((!face_orientation) && (!face_rotation)) || + ((face_orientation) && (face_rotation))) + { + // We flip across the diagonal + // This is the local face dof offset + this + ->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local_face_dof, case_no) = j + i * n - local_face_dof; + } + else + { + // Offset is zero + this + ->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local_face_dof, case_no) = 0; + } // if face needs dof permutation + + // Get new local face_dof by adding offset + const unsigned int new_local_face_dof = + local_face_dof + + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local_face_dof, case_no); + // compute new row and column index + i = new_local_face_dof % n; + j = new_local_face_dof / n; + + /* + * Now compute if a sign change is necessary. This is done for the + * case of face_orientation==true + */ + // flip = false, rotation=true + if (!face_flip && face_rotation) + { + this + ->adjust_quad_dof_sign_for_face_orientation_table[face_no]( + local_face_dof, case_no) = ((j % 2) == 1); + } + // flip = true, rotation=false + else if (face_flip && !face_rotation) + { + // This case is symmetric (although row and column may be + // switched) + this + ->adjust_quad_dof_sign_for_face_orientation_table[face_no]( + local_face_dof, case_no) = + ((j % 2) == 1) != ((i % 2) == 1); + } + // flip = true, rotation=true + else if (face_flip && face_rotation) + { + this + ->adjust_quad_dof_sign_for_face_orientation_table[face_no]( + local_face_dof, case_no) = ((i % 2) == 1); + } + /* + * flip = false, rotation=false => nothing to do + */ + + /* + * If face_orientation==false the sign flip is exactly the + * opposite. + */ + if (!face_orientation) + this->adjust_quad_dof_sign_for_face_orientation_table[face_no]( + local_face_dof, case_no) = + !this + ->adjust_quad_dof_sign_for_face_orientation_table[face_no]( + local_face_dof, case_no); + } // case_no + } // local_face_dof } diff --git a/source/fe/fe_raviart_thomas_nodal.cc b/source/fe/fe_raviart_thomas_nodal.cc index e59a7f6f38..4e2eccf864 100644 --- a/source/fe/fe_raviart_thomas_nodal.cc +++ b/source/fe/fe_raviart_thomas_nodal.cc @@ -170,42 +170,44 @@ FE_RaviartThomasNodal< { unsigned int i = local % n, j = local / n; - // face_orientation=false, face_flip=false, face_rotation=false - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 0) = + // face_orientation=false, face_rotation=false, face_flip=false + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(false, false, false)) = j + i * n - local; - // face_orientation=false, face_flip=false, face_rotation=true - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 1) = + // face_orientation=false, face_rotation=true, face_flip=false + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(false, true, false)) = i + (n - 1 - j) * n - local; - // face_orientation=false, face_flip=true, face_rotation=false - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 2) = + // face_orientation=false, face_rotation=false, face_flip=true + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(false, false, true)) = (n - 1 - j) + (n - 1 - i) * n - local; - // face_orientation=false, face_flip=true, face_rotation=true - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 3) = + // face_orientation=false, face_rotation=true, face_flip=true + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(false, true, true)) = (n - 1 - i) + j * n - local; - // face_orientation=true, face_flip=false, face_rotation=false - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 4) = 0; - // face_orientation=true, face_flip=false, face_rotation=true - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 5) = + // face_orientation=true, face_rotation=false, face_flip=false + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(true, false, false)) = 0; + // face_orientation=true, face_rotation=true, face_flip=false + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(true, true, false)) = j + (n - 1 - i) * n - local; - // face_orientation=true, face_flip=true, face_rotation=false - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 6) = + // face_orientation=true, face_rotation=false, face_flip=true + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(true, false, true)) = (n - 1 - i) + (n - 1 - j) * n - local; - // face_orientation=true, face_flip=true, face_rotation=true - this->adjust_quad_dof_index_for_face_orientation_table[face_no](local, - 7) = + // face_orientation=true, face_rotation=true, face_flip=true + this->adjust_quad_dof_index_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(true, true, true)) = (n - 1 - j) + i * n - local; // for face_orientation == false, we need to switch the sign - for (unsigned int i = 0; i < 4; ++i) - this->adjust_quad_dof_sign_for_face_orientation_table[face_no](local, - i) = 1; + for (const bool rotation : {false, true}) + for (const bool flip : {false, true}) + this->adjust_quad_dof_sign_for_face_orientation_table[face_no]( + local, internal::combined_face_orientation(false, rotation, flip)) = + 1; } } -- 2.39.5