--- /dev/null
+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)
--- /dev/null
+New: The ordering strategy DoFRenumbering::lexicographic() has been added.
+(Michał Wichrowski, 2025/07/01)
std::vector<types::global_dof_index> &new_dof_indices,
const DoFHandler<dim, spacedim> &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 <int dim>
+ void
+ lexicographic(DoFHandler<dim> &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 <int dim>
+ void
+ compute_lexicographic(std::vector<types::global_dof_index> &new_dof_indices,
+ const DoFHandler<dim> &handler,
+ const double tolerance = 1e-12);
+
/**
* @}
*/
std::vector<unsigned int>
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 <int dim>
+ std::pair<std::vector<unsigned int>, std::vector<unsigned int>>
+ 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
return Utilities::invert_permutation(
hierarchic_to_lexicographic_numbering<dim>(degree));
}
+
+
+
+ template <int dim>
+ std::pair<std::vector<unsigned int>, std::vector<unsigned int>>
+ 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<dim>(n_1d);
+ const unsigned int dofs_per_face_patch =
+ Utilities::fixed_power<dim - 1>(n_1d) * n_ort_1d;
+
+ std::array<unsigned int, dim> 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<std::vector<unsigned int>, 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<unsigned int, dim> 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<unsigned int, dim> 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<unsigned int> lex_to_hie =
+ lexicographic_to_hierarchic_numbering<dim>(degree);
+
+ for (auto &result : results)
+ {
+ std::vector<unsigned int> 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
+ template <int dim>
+ void
+ lexicographic(DoFHandler<dim> &dof_handler, const double tolerance)
+ {
+ std::vector<types::global_dof_index> renumbering;
+ compute_lexicographic(renumbering, dof_handler, tolerance);
+ dof_handler.renumber_dofs(renumbering);
+ }
+
+
+ template <int dim>
+ void
+ compute_lexicographic(std::vector<types::global_dof_index> &new_dof_indices,
+ const DoFHandler<dim> &dof_handler,
+ const double tolerance)
+ {
+ Assert(dynamic_cast<const parallel::distributed::Triangulation<dim> *>(
+ &dof_handler.get_triangulation()) == nullptr,
+ ExcMessage(
+ "Lexicographic renumbering is not implemented for distributed "
+ "triangulations."));
+
+ std::map<types::global_dof_index, Point<dim>> dof_location_map;
+ DoFTools::map_dofs_to_support_points(MappingQ1<dim>(),
+ dof_handler,
+ dof_location_map);
+ std::vector<std::pair<types::global_dof_index, Point<dim>>>
+ 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<types::global_dof_index, Point<dim>> &p1,
+ const std::pair<types::global_dof_index, Point<dim>> &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 <int dim,
int spacedim,
typename Number,
const unsigned int,
const bool,
const std::vector<types::global_dof_index> &);
+
+ template void
+ lexicographic(DoFHandler<deal_II_dimension> &, const double);
\}
#endif
}
template std::vector<unsigned int>
lexicographic_to_hierarchic_numbering<deal_II_dimension>(unsigned int);
+ template std::pair<std::vector<unsigned int>, std::vector<unsigned int>>
+ cell_to_face_patch<deal_II_dimension>(const unsigned int &,
+ const unsigned int &,
+ const bool &,
+ const bool &);
+
#endif
\}
}
--- /dev/null
+// ------------------------------------------------------------------------
+//
+// 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 <deal.II/base/logstream.h>
+#include <deal.II/base/point.h>
+#include <deal.II/base/utilities.h>
+
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_renumbering.h>
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_system.h>
+#include <deal.II/fe/mapping_q1.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+
+#include <vector>
+
+#include "../tests.h"
+
+template <int dim>
+void
+test(const unsigned int fe_degree, unsigned int n_components)
+{
+ Triangulation<dim> 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<dim> fe(FE_Q<dim>(fe_degree), n_components);
+ DoFHandler<dim> 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<types::global_dof_index> local_dof_indices(fe.dofs_per_cell);
+
+ std::map<types::global_dof_index, Point<dim>> support_points;
+ MappingQ1<dim> 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;
+}
--- /dev/null
+
+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::
--- /dev/null
+// ------------------------------------------------------------------------
+//
+// 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 <deal.II/fe/fe_tools.h>
+
+#include "../tests.h"
+
+
+
+template <int dim>
+void
+test(unsigned int degree,
+ unsigned int direction,
+ bool cell_hierarchical_numbering,
+ bool is_continuous)
+{
+ std::pair<std::vector<unsigned int>, std::vector<unsigned int>> result =
+ FETools::cell_to_face_patch<dim>(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);
+}
--- /dev/null
+
+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