From 947d2e4271ec902df9bf3182c7b59ff37460f087 Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Sun, 5 Apr 2020 16:08:47 +0200 Subject: [PATCH] Add additional lexicographic_to_hierarchic_numbering functions --- include/deal.II/fe/fe_tools.h | 21 +++++++++++++++++++- include/deal.II/fe/fe_tools.templates.h | 26 +++++++++++++++++++++++++ source/fe/fe_tools.cc | 14 +++++++++++++ source/fe/fe_tools.inst.in | 6 ++++++ 4 files changed, 66 insertions(+), 1 deletion(-) diff --git a/include/deal.II/fe/fe_tools.h b/include/deal.II/fe/fe_tools.h index 1d88e3e817..860f58817c 100644 --- a/include/deal.II/fe/fe_tools.h +++ b/include/deal.II/fe/fe_tools.h @@ -909,12 +909,23 @@ namespace FETools * correct size, which is equal to the number of degrees of freedom in the * finite element. */ - template void hierarchic_to_lexicographic_numbering(unsigned int degree, std::vector &h2l); + /** + * Like the previous function but instead of returning its result through + * the last argument return it as a value. + */ + template + std::vector + hierarchic_to_lexicographic_numbering(unsigned int degree); + + /** + * Like the previous functions but using a FiniteElementData instead of the + * polynomial degree. + */ template void hierarchic_to_lexicographic_numbering(const FiniteElementData &fe_data, @@ -934,6 +945,14 @@ namespace FETools * about the above function are also valid here. */ template + std::vector + lexicographic_to_hierarchic_numbering(unsigned int degree); + + /** + * Like the previous functions but using a FiniteElementData instead of the + * polynomial degree. + */ + template void lexicographic_to_hierarchic_numbering(const FiniteElementData &fe_data, std::vector & l2h); diff --git a/include/deal.II/fe/fe_tools.templates.h b/include/deal.II/fe/fe_tools.templates.h index d8f4c374b4..b4a594fa3a 100644 --- a/include/deal.II/fe/fe_tools.templates.h +++ b/include/deal.II/fe/fe_tools.templates.h @@ -3146,6 +3146,11 @@ namespace FETools // order switch (dim) { + case 0: + { + h2l[0] = 0; + break; + } case 1: { h2l[0] = 0; @@ -3288,6 +3293,17 @@ namespace FETools + template + std::vector + hierarchic_to_lexicographic_numbering(const unsigned int degree) + { + std::vector h2l(Utilities::pow(degree + 1, dim)); + hierarchic_to_lexicographic_numbering(degree, h2l); + return h2l; + } + + + template void hierarchic_to_lexicographic_numbering(const FiniteElementData &fe, @@ -3312,6 +3328,16 @@ namespace FETools + template + std::vector + lexicographic_to_hierarchic_numbering(const unsigned int degree) + { + return Utilities::invert_permutation( + hierarchic_to_lexicographic_numbering(degree)); + } + + + template void lexicographic_to_hierarchic_numbering(const FiniteElementData &fe, diff --git a/source/fe/fe_tools.cc b/source/fe/fe_tools.cc index 4b7989991e..c56e6c8725 100644 --- a/source/fe/fe_tools.cc +++ b/source/fe/fe_tools.cc @@ -21,4 +21,18 @@ DEAL_II_NAMESPACE_OPEN /*-------------- Explicit Instantiations -------------------------------*/ #include "fe_tools.inst" +// these do not fit into the templates of the dimension in the inst file +namespace FETools +{ + template void + hierarchic_to_lexicographic_numbering<0>(unsigned int, + std::vector &); + + template std::vector + hierarchic_to_lexicographic_numbering<0>(unsigned int); + + template std::vector + lexicographic_to_hierarchic_numbering<0>(unsigned int); +} // namespace FETools + DEAL_II_NAMESPACE_CLOSE diff --git a/source/fe/fe_tools.inst.in b/source/fe/fe_tools.inst.in index dd3a0fbeb7..80e2507720 100644 --- a/source/fe/fe_tools.inst.in +++ b/source/fe/fe_tools.inst.in @@ -271,6 +271,12 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) const FiniteElementData &, std::vector &); + template std::vector + hierarchic_to_lexicographic_numbering(unsigned int); + + template std::vector + lexicographic_to_hierarchic_numbering(unsigned int); + template std::vector hierarchic_to_lexicographic_numbering( const FiniteElementData &); -- 2.39.5