From b643fc258428d53a1ebe466892c806f1de4769d4 Mon Sep 17 00:00:00 2001 From: David Wells Date: Sun, 24 Mar 2024 22:25:42 -0400 Subject: [PATCH] Add ReferenceCell::get_inverse_combined_orientation(). --- include/deal.II/grid/reference_cell.h | 45 +++++++ source/dofs/dof_tools_constraints.cc | 5 +- source/grid/tria.cc | 23 ++-- .../reference_cell_reverse_orientation.cc | 126 ++++++++++++++++++ .../reference_cell_reverse_orientation.output | 115 ++++++++++++++++ 5 files changed, 298 insertions(+), 16 deletions(-) create mode 100644 tests/grid/reference_cell_reverse_orientation.cc create mode 100644 tests/grid/reference_cell_reverse_orientation.output diff --git a/include/deal.II/grid/reference_cell.h b/include/deal.II/grid/reference_cell.h index f62d52037c..19b760fa1b 100644 --- a/include/deal.II/grid/reference_cell.h +++ b/include/deal.II/grid/reference_cell.h @@ -711,6 +711,14 @@ public: permute_by_combined_orientation(const ArrayView &vertices, const unsigned char orientation) const; + /** + * Return the inverse orientation. This is the value such that calling + * permute_by_combined_orientation() with o and then calling it again + * with get_inverse_combined_orientation(o) is the identity operation. + */ + unsigned char + get_inverse_combined_orientation(const unsigned char orientation) const; + /** * Return a vector of faces a given @p vertex_index belongs to. */ @@ -2975,6 +2983,43 @@ ReferenceCell::permute_by_combined_orientation( +inline unsigned char +ReferenceCell::get_inverse_combined_orientation( + const unsigned char orientation) const +{ + switch (this->kind) + { + case ReferenceCells::Vertex: + // Things are always default-oriented in 1D + return orientation; + + case ReferenceCells::Line: + // the 1d orientations are the identity and a flip: i.e., the identity + // and an involutory mapping + return orientation; + + case ReferenceCells::Triangle: + { + AssertIndexRange(orientation, 6); + constexpr std::array inverses{{0, 1, 2, 5, 4, 3}}; + return inverses[orientation]; + } + case ReferenceCells::Quadrilateral: + { + AssertIndexRange(orientation, 8); + constexpr std::array inverses{ + {0, 1, 2, 7, 4, 5, 6, 3}}; + return inverses[orientation]; + } + default: + DEAL_II_NOT_IMPLEMENTED(); + } + + return std::numeric_limits::max(); +} + + + template <> unsigned int ReferenceCell::vtk_lexicographic_to_node_index<0>( diff --git a/source/dofs/dof_tools_constraints.cc b/source/dofs/dof_tools_constraints.cc index 56f32a40bc..0272036e72 100644 --- a/source/dofs/dof_tools_constraints.cc +++ b/source/dofs/dof_tools_constraints.cc @@ -3780,14 +3780,15 @@ namespace DoFTools // face_flip has to be toggled if face_rotation is true: In case of // inverted orientation, nothing has to be done. + const auto face_reference_cell = face_1->reference_cell(); internal::set_periodicity_constraints( face_1, face_2, transformation, affine_constraints, component_mask, - ::dealii::internal::combined_face_orientation( - orientation, rotation, orientation ? rotation ^ flip : flip), + face_reference_cell.get_inverse_combined_orientation( + combined_orientation), periodicity_factor); } } diff --git a/source/grid/tria.cc b/source/grid/tria.cc index bfd5252b33..4b9c2d6057 100644 --- a/source/grid/tria.cc +++ b/source/grid/tria.cc @@ -16025,21 +16025,16 @@ void Triangulation::update_periodic_face_map() it->orientation, periodic_face_map); + const auto face_reference_cell = + it->cell[0]->reference_cell().face_reference_cell(it->face_idx[0]); // for the other way, we need to invert the orientation - unsigned char inverted_orientation; - { - auto [orientation, rotation, flip] = - internal::split_face_orientation(it->orientation); - flip = orientation ? rotation ^ flip : flip; - inverted_orientation = - internal::combined_face_orientation(orientation, rotation, flip); - } - update_periodic_face_map_recursively(it->cell[1], - it->cell[0], - it->face_idx[1], - it->face_idx[0], - inverted_orientation, - periodic_face_map); + update_periodic_face_map_recursively( + it->cell[1], + it->cell[0], + it->face_idx[1], + it->face_idx[0], + face_reference_cell.get_inverse_combined_orientation(it->orientation), + periodic_face_map); } // check consistency diff --git a/tests/grid/reference_cell_reverse_orientation.cc b/tests/grid/reference_cell_reverse_orientation.cc new file mode 100644 index 0000000000..df6c2b0b38 --- /dev/null +++ b/tests/grid/reference_cell_reverse_orientation.cc @@ -0,0 +1,126 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2024 - 2024 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + + +// Test that we get the correct reversed orientations (i.e., the composition of +// permutations is the identity operation) + +#include +#include + +#include "../tests.h" + +int +main() +{ + initlog(); + + deallog << "lines" << std::endl; + { + const std::array inverse_permutations{{0u, 1u}}; + + for (unsigned char o = 0; o < 2; ++o) + { + deallog << "o = " << int(o) << std::endl; + std::array vs{{0u, 1u}}; + + const auto result1 = + ReferenceCells::Line.permute_by_combined_orientation( + make_const_array_view(vs), o); + + ArrayView view1(result1.data(), result1.size()); + const auto result2 = + ReferenceCells::Line.permute_by_combined_orientation( + view1, inverse_permutations[o]); + + for (const auto &v : result2) + deallog << " " << v << std::endl; + } + } + + deallog << "triangles" << std::endl; + { + const std::array inverse_permutations{ + {0u, 1u, 2u, 5u, 4u, 3u}}; + + for (unsigned char o = 0; o < 6; ++o) + { + deallog << "o = " << int(o) << std::endl; + std::array vs{{0u, 1u, 2u}}; + + const auto result1 = + ReferenceCells::Triangle.permute_by_combined_orientation( + make_const_array_view(vs), o); + + ArrayView view1(result1.data(), result1.size()); + const auto result2 = + ReferenceCells::Triangle.permute_by_combined_orientation( + view1, inverse_permutations[o]); + + for (const auto &v : result2) + deallog << " " << v << std::endl; + } + } + + deallog << "quadrilaterals" << std::endl; + { + const std::array inverse_permutations{ + {0u, 1u, 2u, 7u, 4u, 5u, 6u, 3u}}; + + for (unsigned char o = 0; o < 8; ++o) + { + deallog << "o = " << int(o) << std::endl; + std::array vs{{0u, 1u, 2u, 3u}}; + + const auto result1 = + ReferenceCells::Quadrilateral.permute_by_combined_orientation( + make_const_array_view(vs), o); + + ArrayView view1(result1.data(), result1.size()); + const auto result2 = + ReferenceCells::Quadrilateral.permute_by_combined_orientation( + view1, inverse_permutations[o]); + + for (const auto &v : result2) + deallog << " " << v << std::endl; + } + } + + // Verify that the manual version created the same results. + deallog << "quadrilaterals (manual)" << std::endl; + { + for (unsigned char o = 0; o < 8; ++o) + { + deallog << "o = " << int(o) << std::endl; + std::array vs{{0u, 1u, 2u, 3u}}; + + const auto result1 = + ReferenceCells::Quadrilateral.permute_by_combined_orientation( + make_const_array_view(vs), o); + + ArrayView view1(result1.data(), result1.size()); + + auto [orientation, rotation, flip] = + internal::split_face_orientation(o); + flip = orientation ? rotation ^ flip : flip; + const auto result2 = + ReferenceCells::Quadrilateral.permute_by_combined_orientation( + view1, + internal::combined_face_orientation(orientation, rotation, flip)); + + for (const auto &v : result2) + deallog << " " << v << std::endl; + } + } +} diff --git a/tests/grid/reference_cell_reverse_orientation.output b/tests/grid/reference_cell_reverse_orientation.output new file mode 100644 index 0000000000..11e5d0cb17 --- /dev/null +++ b/tests/grid/reference_cell_reverse_orientation.output @@ -0,0 +1,115 @@ + +DEAL::lines +DEAL::o = 0 +DEAL:: 0 +DEAL:: 1 +DEAL::o = 1 +DEAL:: 0 +DEAL:: 1 +DEAL::triangles +DEAL::o = 0 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL::o = 1 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL::o = 2 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL::o = 3 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL::o = 4 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL::o = 5 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL::quadrilaterals +DEAL::o = 0 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 1 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 2 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 3 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 4 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 5 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 6 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 7 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::quadrilaterals (manual) +DEAL::o = 0 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 1 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 2 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 3 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 4 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 5 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 6 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 +DEAL::o = 7 +DEAL:: 0 +DEAL:: 1 +DEAL:: 2 +DEAL:: 3 -- 2.39.5