From f3cf8e75e918d89aea1d237f48fcfdbc4dddfa01 Mon Sep 17 00:00:00 2001 From: =?utf8?q?Micha=C5=82=20Wichrowski?= Date: Fri, 4 Jul 2025 15:35:25 +0200 Subject: [PATCH] Renumber lexycographic + face patch dof numbering --- doc/news/changes/minor/20250502Wichrowski | 4 + doc/news/changes/minor/20250701Wichrowski | 2 + include/deal.II/dofs/dof_renumbering.h | 32 ++++++++ include/deal.II/fe/fe_tools.h | 28 +++++++ include/deal.II/fe/fe_tools.templates.h | 98 +++++++++++++++++++++++ source/dofs/dof_renumbering.cc | 66 +++++++++++++++ source/dofs/dof_renumbering.inst.in | 3 + source/fe/fe_tools.inst.in | 6 ++ tests/dofs/dof_renumbering_11.cc | 93 +++++++++++++++++++++ tests/dofs/dof_renumbering_11.output | 65 +++++++++++++++ tests/fe/cell_face_lex.cc | 91 +++++++++++++++++++++ tests/fe/cell_face_lex.output | 40 +++++++++ 12 files changed, 528 insertions(+) create mode 100644 doc/news/changes/minor/20250502Wichrowski create mode 100644 doc/news/changes/minor/20250701Wichrowski create mode 100644 tests/dofs/dof_renumbering_11.cc create mode 100644 tests/dofs/dof_renumbering_11.output create mode 100644 tests/fe/cell_face_lex.cc create mode 100644 tests/fe/cell_face_lex.output diff --git a/doc/news/changes/minor/20250502Wichrowski b/doc/news/changes/minor/20250502Wichrowski new file mode 100644 index 0000000000..796bc13e55 --- /dev/null +++ b/doc/news/changes/minor/20250502Wichrowski @@ -0,0 +1,4 @@ +New: A new function FETools::cell_to_face_lexicographic() to generate a +mapping from cell-local DoFs to lexicographic ordering of DoFs on two adjacent cells +has been added. +(Michał Wichrowski, 2025/05/03) diff --git a/doc/news/changes/minor/20250701Wichrowski b/doc/news/changes/minor/20250701Wichrowski new file mode 100644 index 0000000000..0f11fec48a --- /dev/null +++ b/doc/news/changes/minor/20250701Wichrowski @@ -0,0 +1,2 @@ +New: The ordering strategy DoFRenumbering::lexicographic() has been added. +(Michał Wichrowski, 2025/07/01) diff --git a/include/deal.II/dofs/dof_renumbering.h b/include/deal.II/dofs/dof_renumbering.h index a2c5a9ca0b..39a9a36dfa 100644 --- a/include/deal.II/dofs/dof_renumbering.h +++ b/include/deal.II/dofs/dof_renumbering.h @@ -1273,6 +1273,38 @@ namespace DoFRenumbering std::vector &new_dof_indices, const DoFHandler &dof_handler); + /** + * Reorder the degrees of freedom in lexicographic order of support points. + * The tolerance is used to determine whether two support points are equal. + * + * This is particularly useful for finite element spaces that are defined on + * tensor product grids, such as the ones generated by + * GridGenerator::hyper_rectangle() or GridGenerator::hyper_cube() or + * their subdivided variants. + * With this ordering certain higher-dimension finite element matrices can be + * expressed as Kronecker products of lower-dimensional matrices. + * + * @warning parallel::distributed::Triangulation objects are not supported. + * The finite element must have support points. + */ + template + void + lexicographic(DoFHandler &dof_handler, const double tolerance = 1e-12); + + /** + * Compute the renumbering vector needed by the lexicographic() function. + * Does not perform the renumbering on the @p DoFHandler dofs but returns the + * renumbering vector. + * The tolerance is treated as an absolute tolerance for comparing support + * point coordinates. To use a relative tolerance, multiply a relative value + * by dof_handler.begin(0)->diameter() or a similar characteristic length. + */ + template + void + compute_lexicographic(std::vector &new_dof_indices, + const DoFHandler &handler, + const double tolerance = 1e-12); + /** * @} */ diff --git a/include/deal.II/fe/fe_tools.h b/include/deal.II/fe/fe_tools.h index 69ed37fe1c..c34cec1ab6 100644 --- a/include/deal.II/fe/fe_tools.h +++ b/include/deal.II/fe/fe_tools.h @@ -931,6 +931,34 @@ namespace FETools std::vector lexicographic_to_hierarchic_numbering(unsigned int degree); + /** + * Given the polynomial degree, direction, and numbering options, this + * function returns a pair of vectors mapping cell DoFs to face patch DoFs + * in lexicographic order. + * + * This function works for scalar finite elements, specifically those derived + * from FE_Q and FE_DGQ. It computes the mapping from the DoFs of + * two adjacent cells sharing a face perpendicular to @p direction in + * reference space, to the DoFs on the face patch. The result is a pair of + * vectors: the first for the DoFs on the first cell, the second for the + * DoFs on the second cell. The @p cell_hierarchical_numbering flag + * determines whether hierarchical or lexicographic numbering is assumed for + * the cell DoFs. The @p is_continuous flag indicates whether the finite + * element space is continuous (i.e., FE_Q). + * + * @param degree Polynomial degree of the finite element. + * @param direction Axis perpendicular to the considered face. + * @param cell_hierarchical_numbering If true, use hierarchical numbering for cell DoFs. + * @param is_continuous If true, assumes the finite element space is continuous (FE_Q). + * @return A pair of vectors mapping cell DoFs to face patch DoFs in lexicographic order. + */ + template + std::pair, std::vector> + cell_to_face_patch(const unsigned int °ree, + const unsigned int &direction, + const bool &cell_hierarchical_numbering, + const bool &is_continuous); + /** * A namespace that contains functions that help setting up internal * data structures when implementing FiniteElement which are build diff --git a/include/deal.II/fe/fe_tools.templates.h b/include/deal.II/fe/fe_tools.templates.h index 420be5f268..65167f6c53 100644 --- a/include/deal.II/fe/fe_tools.templates.h +++ b/include/deal.II/fe/fe_tools.templates.h @@ -3243,6 +3243,104 @@ namespace FETools return Utilities::invert_permutation( hierarchic_to_lexicographic_numbering(degree)); } + + + + template + std::pair, std::vector> + cell_to_face_patch(const unsigned int °ree, + const unsigned int &direction, + const bool &cell_hierarchical_numbering, + const bool &is_continuous) + { + AssertThrow(dim > 1, ExcInvalidFEDimension(1, dim)); + AssertThrow(direction < dim, ExcInvalidFEDimension(direction, dim)); + + // n_1d: DoFs per 1D edge (parallel to face) + const unsigned int n_1d = degree + 1; + // n_ort: DoFs in the direction orthogonal to the face in the combined face + // space For discontinuous: 2*n (n for cell0 face, n for cell1 face) For + // continuous: 2*n - 1 (shared DoFs) + const unsigned int n_ort_1d = 2 * n_1d - (is_continuous ? 1 : 0); + + const unsigned int dofs_per_cell = Utilities::fixed_power(n_1d); + const unsigned int dofs_per_face_patch = + Utilities::fixed_power(n_1d) * n_ort_1d; + + std::array face_patch_dofs_in_direction; + for (unsigned int d = 0; d < dim; ++d) + face_patch_dofs_in_direction[d] = n_1d; + face_patch_dofs_in_direction[direction] = n_ort_1d; + + // Initialize maps for cell0 and cell1 DoFs to face DoFs + std::array, 2> results; + results[0].resize(dofs_per_cell, numbers::invalid_unsigned_int); + results[1].resize(dofs_per_cell, numbers::invalid_unsigned_int); + + // Compute with lexicographic first + for (unsigned int i = 0; i < dofs_per_face_patch; ++i) + { + // Decompose i into multiindex + std::array face_multiindex; + unsigned int temp = i; + for (unsigned int d = 0; d < dim; ++d) + { + face_multiindex[d] = temp % face_patch_dofs_in_direction[d]; + temp /= face_patch_dofs_in_direction[d]; + } + + std::array cell_multiindex = face_multiindex; + + // compute corresponding cell + unsigned int cell = face_multiindex[direction] / n_1d; + cell_multiindex[direction] = face_multiindex[direction] % n_1d; + if (is_continuous) + cell_multiindex[direction] += cell; + + { + unsigned int cell_index = 0; + unsigned int stride = 1; + for (unsigned int d = 0; d < dim; ++d) + { + cell_index += cell_multiindex[d] * stride; + stride *= n_1d; + } + results[cell][cell_index] = i; + } + + if (face_multiindex[direction] == n_1d - 1 && is_continuous) + { + // fill the other cell + unsigned other_cell = 1; + cell_multiindex[direction] = 0; + unsigned int cell_index = 0; + unsigned int stride = 1; + for (unsigned int d = 0; d < dim; ++d) + { + cell_index += cell_multiindex[d] * stride; + stride *= n_1d; + } + results[other_cell][cell_index] = i; + } + } + // Now we have the DoFs in lexicographic order + if (!cell_hierarchical_numbering) + return std::make_pair(results[0], results[1]); + + // Renumbering required for hierarchical numbering + std::vector lex_to_hie = + lexicographic_to_hierarchic_numbering(degree); + + for (auto &result : results) + { + std::vector renumbered_result = result; + for (unsigned int i = 0; i < dofs_per_cell; ++i) + renumbered_result[lex_to_hie[i]] = result[i]; + result = renumbered_result; + } + return std::make_pair(results[0], results[1]); + } + } // namespace FETools diff --git a/source/dofs/dof_renumbering.cc b/source/dofs/dof_renumbering.cc index a0b0f3240e..8600dbf6d5 100644 --- a/source/dofs/dof_renumbering.cc +++ b/source/dofs/dof_renumbering.cc @@ -2442,6 +2442,72 @@ namespace DoFRenumbering + template + void + lexicographic(DoFHandler &dof_handler, const double tolerance) + { + std::vector renumbering; + compute_lexicographic(renumbering, dof_handler, tolerance); + dof_handler.renumber_dofs(renumbering); + } + + + template + void + compute_lexicographic(std::vector &new_dof_indices, + const DoFHandler &dof_handler, + const double tolerance) + { + Assert(dynamic_cast *>( + &dof_handler.get_triangulation()) == nullptr, + ExcMessage( + "Lexicographic renumbering is not implemented for distributed " + "triangulations.")); + + std::map> dof_location_map; + DoFTools::map_dofs_to_support_points(MappingQ1(), + dof_handler, + dof_location_map); + std::vector>> + dof_location_vector; + + dof_location_vector.reserve(dof_location_map.size()); + for (const auto &s : dof_location_map) + dof_location_vector.push_back(s); + + std::sort( + dof_location_vector.begin(), + dof_location_vector.end(), + [&](const std::pair> &p1, + const std::pair> &p2) -> bool { + for (int i = dim - 1; i >= 0; --i) + { + const double diff = p1.second(i) - p2.second(i); + // Check if p1 is significantly smaller than p2 in this dimension + if (diff < -tolerance) + return true; + // Check if p1 is significantly larger than p2 in this dimension + if (diff > tolerance) + return false; + // Otherwise, coordinates are considered equal within tolerance for + // this dimension, continue to the next dimension. + } + // All coordinates are within tolerance, points are considered + // equivalent. Use the original DoF index as a tie-breaker to preserve + // relative order. + return p1.first < p2.first; + }); + + new_dof_indices.resize(dof_location_vector.size(), + numbers::invalid_dof_index); + + for (types::global_dof_index dof = 0; dof < dof_location_vector.size(); + ++dof) + new_dof_indices[dof_location_vector[dof].first] = dof; + } + + + template &); + + template void + lexicographic(DoFHandler &, const double); \} #endif } diff --git a/source/fe/fe_tools.inst.in b/source/fe/fe_tools.inst.in index 7c9e101c04..7dfec11c67 100644 --- a/source/fe/fe_tools.inst.in +++ b/source/fe/fe_tools.inst.in @@ -262,6 +262,12 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) template std::vector lexicographic_to_hierarchic_numbering(unsigned int); + template std::pair, std::vector> + cell_to_face_patch(const unsigned int &, + const unsigned int &, + const bool &, + const bool &); + #endif \} } diff --git a/tests/dofs/dof_renumbering_11.cc b/tests/dofs/dof_renumbering_11.cc new file mode 100644 index 0000000000..bbf3c9e23d --- /dev/null +++ b/tests/dofs/dof_renumbering_11.cc @@ -0,0 +1,93 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2019 - 2023 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 lexicographic renumbering on a patch + +#include +#include +#include + +#include +#include +#include +#include + +#include +#include +#include + +#include +#include +#include +#include + +#include + +#include "../tests.h" + +template +void +test(const unsigned int fe_degree, unsigned int n_components) +{ + Triangulation tria; + GridGenerator::hyper_cube(tria, 0, 1, true); + tria.refine_global(1); + + deallog << "Testing dim=" << dim << ", degree=" << fe_degree + << ", components=" << n_components << std::endl; + + FESystem fe(FE_Q(fe_degree), n_components); + DoFHandler dof_handler(tria); + + const unsigned int n_dof_1d = (2 * fe_degree + 1); + + dof_handler.distribute_dofs(fe); + DoFRenumbering::lexicographic(dof_handler); + DoFRenumbering::component_wise(dof_handler); + + std::vector local_dof_indices(fe.dofs_per_cell); + + std::map> support_points; + MappingQ1 mapping; + DoFTools::map_dofs_to_support_points(mapping, dof_handler, support_points); + deallog << "Support points:" << std::endl; + for (unsigned int i = 0; i < dof_handler.n_dofs(); ++i) + { + deallog << " " << i << " (" << support_points[i] << ")"; + if ((i + 1) % n_dof_1d == 0) + deallog << std::endl; + } + + deallog << std::endl; +} + + +int +main(int argc, char **argv) +{ + initlog(); + + deallog.precision(2); + + test<2>(1, 1); + test<2>(2, 1); + + test<2>(1, 2); + test<2>(2, 2); + + test<3>(2, 1); + + return 0; +} diff --git a/tests/dofs/dof_renumbering_11.output b/tests/dofs/dof_renumbering_11.output new file mode 100644 index 0000000000..c3a2bc6498 --- /dev/null +++ b/tests/dofs/dof_renumbering_11.output @@ -0,0 +1,65 @@ + +DEAL::Testing dim=2, degree=1, components=1 +DEAL::Support points: +DEAL:: 0 (0.0 0.0) 1 (0.50 0.0) 2 (1.0 0.0) +DEAL:: 3 (0.0 0.50) 4 (0.50 0.50) 5 (1.0 0.50) +DEAL:: 6 (0.0 1.0) 7 (0.50 1.0) 8 (1.0 1.0) +DEAL:: +DEAL::Testing dim=2, degree=2, components=1 +DEAL::Support points: +DEAL:: 0 (0.0 0.0) 1 (0.25 0.0) 2 (0.50 0.0) 3 (0.75 0.0) 4 (1.0 0.0) +DEAL:: 5 (0.0 0.25) 6 (0.25 0.25) 7 (0.50 0.25) 8 (0.75 0.25) 9 (1.0 0.25) +DEAL:: 10 (0.0 0.50) 11 (0.25 0.50) 12 (0.50 0.50) 13 (0.75 0.50) 14 (1.0 0.50) +DEAL:: 15 (0.0 0.75) 16 (0.25 0.75) 17 (0.50 0.75) 18 (0.75 0.75) 19 (1.0 0.75) +DEAL:: 20 (0.0 1.0) 21 (0.25 1.0) 22 (0.50 1.0) 23 (0.75 1.0) 24 (1.0 1.0) +DEAL:: +DEAL::Testing dim=2, degree=1, components=2 +DEAL::Support points: +DEAL:: 0 (0.0 0.0) 1 (0.50 0.0) 2 (1.0 0.0) +DEAL:: 3 (0.0 0.50) 4 (0.50 0.50) 5 (1.0 0.50) +DEAL:: 6 (0.0 1.0) 7 (0.50 1.0) 8 (1.0 1.0) +DEAL:: 9 (0.0 0.0) 10 (0.50 0.0) 11 (1.0 0.0) +DEAL:: 12 (0.0 0.50) 13 (0.50 0.50) 14 (1.0 0.50) +DEAL:: 15 (0.0 1.0) 16 (0.50 1.0) 17 (1.0 1.0) +DEAL:: +DEAL::Testing dim=2, degree=2, components=2 +DEAL::Support points: +DEAL:: 0 (0.0 0.0) 1 (0.25 0.0) 2 (0.50 0.0) 3 (0.75 0.0) 4 (1.0 0.0) +DEAL:: 5 (0.0 0.25) 6 (0.25 0.25) 7 (0.50 0.25) 8 (0.75 0.25) 9 (1.0 0.25) +DEAL:: 10 (0.0 0.50) 11 (0.25 0.50) 12 (0.50 0.50) 13 (0.75 0.50) 14 (1.0 0.50) +DEAL:: 15 (0.0 0.75) 16 (0.25 0.75) 17 (0.50 0.75) 18 (0.75 0.75) 19 (1.0 0.75) +DEAL:: 20 (0.0 1.0) 21 (0.25 1.0) 22 (0.50 1.0) 23 (0.75 1.0) 24 (1.0 1.0) +DEAL:: 25 (0.0 0.0) 26 (0.25 0.0) 27 (0.50 0.0) 28 (0.75 0.0) 29 (1.0 0.0) +DEAL:: 30 (0.0 0.25) 31 (0.25 0.25) 32 (0.50 0.25) 33 (0.75 0.25) 34 (1.0 0.25) +DEAL:: 35 (0.0 0.50) 36 (0.25 0.50) 37 (0.50 0.50) 38 (0.75 0.50) 39 (1.0 0.50) +DEAL:: 40 (0.0 0.75) 41 (0.25 0.75) 42 (0.50 0.75) 43 (0.75 0.75) 44 (1.0 0.75) +DEAL:: 45 (0.0 1.0) 46 (0.25 1.0) 47 (0.50 1.0) 48 (0.75 1.0) 49 (1.0 1.0) +DEAL:: +DEAL::Testing dim=3, degree=2, components=1 +DEAL::Support points: +DEAL:: 0 (0.0 0.0 0.0) 1 (0.25 0.0 0.0) 2 (0.50 0.0 0.0) 3 (0.75 0.0 0.0) 4 (1.0 0.0 0.0) +DEAL:: 5 (0.0 0.25 0.0) 6 (0.25 0.25 0.0) 7 (0.50 0.25 0.0) 8 (0.75 0.25 0.0) 9 (1.0 0.25 0.0) +DEAL:: 10 (0.0 0.50 0.0) 11 (0.25 0.50 0.0) 12 (0.50 0.50 0.0) 13 (0.75 0.50 0.0) 14 (1.0 0.50 0.0) +DEAL:: 15 (0.0 0.75 0.0) 16 (0.25 0.75 0.0) 17 (0.50 0.75 0.0) 18 (0.75 0.75 0.0) 19 (1.0 0.75 0.0) +DEAL:: 20 (0.0 1.0 0.0) 21 (0.25 1.0 0.0) 22 (0.50 1.0 0.0) 23 (0.75 1.0 0.0) 24 (1.0 1.0 0.0) +DEAL:: 25 (0.0 0.0 0.25) 26 (0.25 0.0 0.25) 27 (0.50 0.0 0.25) 28 (0.75 0.0 0.25) 29 (1.0 0.0 0.25) +DEAL:: 30 (0.0 0.25 0.25) 31 (0.25 0.25 0.25) 32 (0.50 0.25 0.25) 33 (0.75 0.25 0.25) 34 (1.0 0.25 0.25) +DEAL:: 35 (0.0 0.50 0.25) 36 (0.25 0.50 0.25) 37 (0.50 0.50 0.25) 38 (0.75 0.50 0.25) 39 (1.0 0.50 0.25) +DEAL:: 40 (0.0 0.75 0.25) 41 (0.25 0.75 0.25) 42 (0.50 0.75 0.25) 43 (0.75 0.75 0.25) 44 (1.0 0.75 0.25) +DEAL:: 45 (0.0 1.0 0.25) 46 (0.25 1.0 0.25) 47 (0.50 1.0 0.25) 48 (0.75 1.0 0.25) 49 (1.0 1.0 0.25) +DEAL:: 50 (0.0 0.0 0.50) 51 (0.25 0.0 0.50) 52 (0.50 0.0 0.50) 53 (0.75 0.0 0.50) 54 (1.0 0.0 0.50) +DEAL:: 55 (0.0 0.25 0.50) 56 (0.25 0.25 0.50) 57 (0.50 0.25 0.50) 58 (0.75 0.25 0.50) 59 (1.0 0.25 0.50) +DEAL:: 60 (0.0 0.50 0.50) 61 (0.25 0.50 0.50) 62 (0.50 0.50 0.50) 63 (0.75 0.50 0.50) 64 (1.0 0.50 0.50) +DEAL:: 65 (0.0 0.75 0.50) 66 (0.25 0.75 0.50) 67 (0.50 0.75 0.50) 68 (0.75 0.75 0.50) 69 (1.0 0.75 0.50) +DEAL:: 70 (0.0 1.0 0.50) 71 (0.25 1.0 0.50) 72 (0.50 1.0 0.50) 73 (0.75 1.0 0.50) 74 (1.0 1.0 0.50) +DEAL:: 75 (0.0 0.0 0.75) 76 (0.25 0.0 0.75) 77 (0.50 0.0 0.75) 78 (0.75 0.0 0.75) 79 (1.0 0.0 0.75) +DEAL:: 80 (0.0 0.25 0.75) 81 (0.25 0.25 0.75) 82 (0.50 0.25 0.75) 83 (0.75 0.25 0.75) 84 (1.0 0.25 0.75) +DEAL:: 85 (0.0 0.50 0.75) 86 (0.25 0.50 0.75) 87 (0.50 0.50 0.75) 88 (0.75 0.50 0.75) 89 (1.0 0.50 0.75) +DEAL:: 90 (0.0 0.75 0.75) 91 (0.25 0.75 0.75) 92 (0.50 0.75 0.75) 93 (0.75 0.75 0.75) 94 (1.0 0.75 0.75) +DEAL:: 95 (0.0 1.0 0.75) 96 (0.25 1.0 0.75) 97 (0.50 1.0 0.75) 98 (0.75 1.0 0.75) 99 (1.0 1.0 0.75) +DEAL:: 100 (0.0 0.0 1.0) 101 (0.25 0.0 1.0) 102 (0.50 0.0 1.0) 103 (0.75 0.0 1.0) 104 (1.0 0.0 1.0) +DEAL:: 105 (0.0 0.25 1.0) 106 (0.25 0.25 1.0) 107 (0.50 0.25 1.0) 108 (0.75 0.25 1.0) 109 (1.0 0.25 1.0) +DEAL:: 110 (0.0 0.50 1.0) 111 (0.25 0.50 1.0) 112 (0.50 0.50 1.0) 113 (0.75 0.50 1.0) 114 (1.0 0.50 1.0) +DEAL:: 115 (0.0 0.75 1.0) 116 (0.25 0.75 1.0) 117 (0.50 0.75 1.0) 118 (0.75 0.75 1.0) 119 (1.0 0.75 1.0) +DEAL:: 120 (0.0 1.0 1.0) 121 (0.25 1.0 1.0) 122 (0.50 1.0 1.0) 123 (0.75 1.0 1.0) 124 (1.0 1.0 1.0) +DEAL:: diff --git a/tests/fe/cell_face_lex.cc b/tests/fe/cell_face_lex.cc new file mode 100644 index 0000000000..aa394561a7 --- /dev/null +++ b/tests/fe/cell_face_lex.cc @@ -0,0 +1,91 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2009 - 2021 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. +// +// ------------------------------------------------------------------------ + + + +// since early 2009, the FEValues objects try to be more efficient by only +// recomputing things like gradients of shape functions if the cell on which +// we are is not a translation of the previous one. in this series of tests we +// make sure that this actually works the way it's supposed to be; in +// particular, if we create a mesh of two identical cells but one has a curved +// boundary, then they are the same if we use a Q1 mapping, but not a Q2 +// mapping. so we test that the mass matrix we get from each of these cells is +// actually different in the latter case, but the same in the former +// +// Check cell to face lexicographic ordering + + +#include + +#include "../tests.h" + + + +template +void +test(unsigned int degree, + unsigned int direction, + bool cell_hierarchical_numbering, + bool is_continuous) +{ + std::pair, std::vector> result = + FETools::cell_to_face_patch(degree, + direction, + cell_hierarchical_numbering, + is_continuous); + + deallog << "Results for degree = " << degree << ", dim = " << dim + << ", direction = " << direction + << ", is_continuous = " << (is_continuous ? "true" : "false") + << ", cell_hierarchical_numbering = " + << (cell_hierarchical_numbering ? "true" : "false") << std::endl; + + deallog << "Cell 0: "; + for (unsigned int val : result.first) + deallog << val << " "; + deallog << std::endl; + + deallog << "Cell 1: "; + for (unsigned int val : result.second) + deallog << val << " "; + deallog << std::endl; +} + + +int +main() +{ + initlog(); + deallog << std::setprecision(3); + + test<2>(2, 0, false, false); + test<2>(2, 0, false, true); + test<2>(2, 0, true, true); + + test<2>(4, 0, false, true); + + test<2>(2, 1, false, false); + test<2>(2, 1, false, true); + + + test<3>(2, 0, false, false); + test<3>(2, 0, false, true); + test<3>(2, 0, true, true); + + test<3>(2, 1, false, false); + test<3>(2, 1, false, true); + + test<3>(2, 2, false, false); + test<3>(2, 2, false, true); +} diff --git a/tests/fe/cell_face_lex.output b/tests/fe/cell_face_lex.output new file mode 100644 index 0000000000..d8264bf2f1 --- /dev/null +++ b/tests/fe/cell_face_lex.output @@ -0,0 +1,40 @@ + +DEAL::Results for degree = 2, dim = 2, direction = 0, is_continuous = false, cell_hierarchical_numbering = false +DEAL::Cell 0: 0 1 2 6 7 8 12 13 14 +DEAL::Cell 1: 3 4 5 9 10 11 15 16 17 +DEAL::Results for degree = 2, dim = 2, direction = 0, is_continuous = true, cell_hierarchical_numbering = false +DEAL::Cell 0: 0 1 2 5 6 7 10 11 12 +DEAL::Cell 1: 2 3 4 7 8 9 12 13 14 +DEAL::Results for degree = 2, dim = 2, direction = 0, is_continuous = true, cell_hierarchical_numbering = true +DEAL::Cell 0: 0 2 10 12 5 7 1 11 6 +DEAL::Cell 1: 2 4 12 14 7 9 3 13 8 +DEAL::Results for degree = 4, dim = 2, direction = 0, is_continuous = true, cell_hierarchical_numbering = false +DEAL::Cell 0: 0 1 2 3 4 9 10 11 12 13 18 19 20 21 22 27 28 29 30 31 36 37 38 39 40 +DEAL::Cell 1: 4 5 6 7 8 13 14 15 16 17 22 23 24 25 26 31 32 33 34 35 40 41 42 43 44 +DEAL::Results for degree = 2, dim = 2, direction = 1, is_continuous = false, cell_hierarchical_numbering = false +DEAL::Cell 0: 0 1 2 3 4 5 6 7 8 +DEAL::Cell 1: 9 10 11 12 13 14 15 16 17 +DEAL::Results for degree = 2, dim = 2, direction = 1, is_continuous = true, cell_hierarchical_numbering = false +DEAL::Cell 0: 0 1 2 3 4 5 6 7 8 +DEAL::Cell 1: 6 7 8 9 10 11 12 13 14 +DEAL::Results for degree = 2, dim = 3, direction = 0, is_continuous = false, cell_hierarchical_numbering = false +DEAL::Cell 0: 0 1 2 6 7 8 12 13 14 18 19 20 24 25 26 30 31 32 36 37 38 42 43 44 48 49 50 +DEAL::Cell 1: 3 4 5 9 10 11 15 16 17 21 22 23 27 28 29 33 34 35 39 40 41 45 46 47 51 52 53 +DEAL::Results for degree = 2, dim = 3, direction = 0, is_continuous = true, cell_hierarchical_numbering = false +DEAL::Cell 0: 0 1 2 5 6 7 10 11 12 15 16 17 20 21 22 25 26 27 30 31 32 35 36 37 40 41 42 +DEAL::Cell 1: 2 3 4 7 8 9 12 13 14 17 18 19 22 23 24 27 28 29 32 33 34 37 38 39 42 43 44 +DEAL::Results for degree = 2, dim = 3, direction = 0, is_continuous = true, cell_hierarchical_numbering = true +DEAL::Cell 0: 0 2 10 12 30 32 40 42 5 7 1 11 35 37 31 41 15 17 25 27 20 22 16 26 6 36 21 +DEAL::Cell 1: 2 4 12 14 32 34 42 44 7 9 3 13 37 39 33 43 17 19 27 29 22 24 18 28 8 38 23 +DEAL::Results for degree = 2, dim = 3, direction = 1, is_continuous = false, cell_hierarchical_numbering = false +DEAL::Cell 0: 0 1 2 3 4 5 6 7 8 18 19 20 21 22 23 24 25 26 36 37 38 39 40 41 42 43 44 +DEAL::Cell 1: 9 10 11 12 13 14 15 16 17 27 28 29 30 31 32 33 34 35 45 46 47 48 49 50 51 52 53 +DEAL::Results for degree = 2, dim = 3, direction = 1, is_continuous = true, cell_hierarchical_numbering = false +DEAL::Cell 0: 0 1 2 3 4 5 6 7 8 15 16 17 18 19 20 21 22 23 30 31 32 33 34 35 36 37 38 +DEAL::Cell 1: 6 7 8 9 10 11 12 13 14 21 22 23 24 25 26 27 28 29 36 37 38 39 40 41 42 43 44 +DEAL::Results for degree = 2, dim = 3, direction = 2, is_continuous = false, cell_hierarchical_numbering = false +DEAL::Cell 0: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 +DEAL::Cell 1: 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 +DEAL::Results for degree = 2, dim = 3, direction = 2, is_continuous = true, cell_hierarchical_numbering = false +DEAL::Cell 0: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 +DEAL::Cell 1: 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 -- 2.39.5