From 66721597cb932fdc44beb4efa356bd59646ecd9d Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Fri, 20 Jun 2025 11:29:56 +0200 Subject: [PATCH] Further cleanup in PolynomialsRaviartThomas --- .../deal.II/base/polynomials_raviart_thomas.h | 31 ++------ .../matrix_free/shape_info.templates.h | 2 +- source/base/polynomials_raviart_thomas.cc | 74 +++++-------------- source/fe/fe_raviart_thomas.cc | 11 ++- source/fe/fe_raviart_thomas_nodal.cc | 22 +++--- 5 files changed, 41 insertions(+), 99 deletions(-) diff --git a/include/deal.II/base/polynomials_raviart_thomas.h b/include/deal.II/base/polynomials_raviart_thomas.h index 85903dc3d6..a08908c219 100644 --- a/include/deal.II/base/polynomials_raviart_thomas.h +++ b/include/deal.II/base/polynomials_raviart_thomas.h @@ -36,7 +36,7 @@ DEAL_II_NAMESPACE_OPEN * This class implements the Hdiv-conforming, * Raviart-Thomas polynomials as described in the book by Brezzi and Fortin. * Most of the functionality comes from the vector-valued anisotropic - * polynomials class . + * polynomials class PolynomialsVectorAnisotropic. * * The Raviart-Thomas polynomials are constructed such that the divergence is * in the tensor product polynomial space Qk. Therefore, the @@ -51,15 +51,6 @@ template class PolynomialsRaviartThomas : public PolynomialsVectorAnisotropic { public: - /** - * Constructor. Creates all basis functions for Raviart-Thomas polynomials - * of given degree in normal and tangential directions. The usual - * FE_RaviartThomas and FE_RaviartThomasNodal classes will use `degree + 1` - * and `degree` in the two directions, respectively. - */ - PolynomialsRaviartThomas(const unsigned int degree_normal, - const unsigned int degree_tangential); - /** * Constructor, using the common Raviart-Thomas space of degree `k + 1` in * normal direction and `k` in the tangential directions. @@ -70,20 +61,6 @@ public: */ PolynomialsRaviartThomas(const unsigned int k); - /** - * Copy constructor. - */ - PolynomialsRaviartThomas(const PolynomialsRaviartThomas &other) = default; - - /** - * Return the number of polynomials in the space without requiring to - * build an object of PolynomialsRaviartThomas. This is required by the - * FiniteElement classes. - */ - static unsigned int - n_polynomials(const unsigned int normal_degree, - const unsigned int tangential_degree); - /** * Variant of the n_polynomials() function taking only a single argument * `degree`, assuming `degree + 1` in the normal direction and `degree` in @@ -92,13 +69,15 @@ public: static unsigned int n_polynomials(const unsigned int degree); + // Make respective two-argument method from base class available + using PolynomialsVectorAnisotropic::n_polynomials; + /** * Compute the lexicographic to hierarchic numbering underlying this class, * computed as a free function. */ static std::vector - get_lexicographic_numbering(const unsigned int normal_degree, - const unsigned int tangential_degree); + get_lexicographic_numbering(const unsigned int degree); /** * @copydoc TensorPolynomialsBase::clone() diff --git a/include/deal.II/matrix_free/shape_info.templates.h b/include/deal.II/matrix_free/shape_info.templates.h index b4a57b73db..f3fc16f505 100644 --- a/include/deal.II/matrix_free/shape_info.templates.h +++ b/include/deal.II/matrix_free/shape_info.templates.h @@ -270,7 +270,7 @@ namespace internal lexicographic_numbering = PolynomialsRaviartThomas::get_lexicographic_numbering( - fe_in.degree, fe_in.degree - 1); + fe_in.degree - 1); // To get the right shape_values of the RT element std::vector lex_normal, lex_tangent; diff --git a/source/base/polynomials_raviart_thomas.cc b/source/base/polynomials_raviart_thomas.cc index f8595a61c0..8ca30dcfa5 100644 --- a/source/base/polynomials_raviart_thomas.cc +++ b/source/base/polynomials_raviart_thomas.cc @@ -26,37 +26,13 @@ DEAL_II_NAMESPACE_OPEN -template -PolynomialsRaviartThomas::PolynomialsRaviartThomas( - const unsigned int normal_degree, - const unsigned int tangential_degree) - : PolynomialsVectorAnisotropic( - normal_degree, - tangential_degree, - get_lexicographic_numbering(normal_degree, tangential_degree)) -{} - - - template PolynomialsRaviartThomas::PolynomialsRaviartThomas(const unsigned int k) - : PolynomialsRaviartThomas(k + 1, k) + : PolynomialsVectorAnisotropic(k + 1, k, get_lexicographic_numbering(k)) {} -template -unsigned int -PolynomialsRaviartThomas::n_polynomials( - const unsigned int normal_degree, - const unsigned int tangential_degree) -{ - return PolynomialsVectorAnisotropic::n_polynomials(normal_degree, - tangential_degree); -} - - - template unsigned int PolynomialsRaviartThomas::n_polynomials(const unsigned int degree) @@ -69,44 +45,36 @@ PolynomialsRaviartThomas::n_polynomials(const unsigned int degree) template std::vector PolynomialsRaviartThomas::get_lexicographic_numbering( - const unsigned int normal_degree, - const unsigned int tangential_degree) + const unsigned int degree) { - const unsigned int n_dofs_face = - Utilities::pow(tangential_degree + 1, dim - 1); + const unsigned int n_dofs_face = Utilities::pow(degree + 1, dim - 1); std::vector lexicographic_numbering; + // component 1 for (unsigned int j = 0; j < n_dofs_face; ++j) { lexicographic_numbering.push_back(j); - if (normal_degree > 1) - for (unsigned int i = n_dofs_face * 2 * dim; - i < n_dofs_face * 2 * dim + normal_degree - 1; - ++i) - lexicographic_numbering.push_back(i + j * (normal_degree - 1)); + for (unsigned int i = n_dofs_face * 2 * dim; + i < n_dofs_face * 2 * dim + degree; + ++i) + lexicographic_numbering.push_back(i + j * degree); lexicographic_numbering.push_back(n_dofs_face + j); } // component 2 - unsigned int layers = (dim == 3) ? tangential_degree + 1 : 1; + unsigned int layers = (dim == 3) ? degree + 1 : 1; for (unsigned int k = 0; k < layers; ++k) { - unsigned int k_add = k * (tangential_degree + 1); - for (unsigned int j = n_dofs_face * 2; - j < n_dofs_face * 2 + tangential_degree + 1; + unsigned int k_add = k * (degree + 1); + for (unsigned int j = n_dofs_face * 2; j < n_dofs_face * 2 + degree + 1; ++j) lexicographic_numbering.push_back(j + k_add); - if (normal_degree > 1) - for (unsigned int i = n_dofs_face * (2 * dim + (normal_degree - 1)); - i < n_dofs_face * (2 * dim + (normal_degree - 1)) + - (normal_degree - 1) * (tangential_degree + 1); - ++i) - { - lexicographic_numbering.push_back(i + k_add * tangential_degree); - } - for (unsigned int j = n_dofs_face * 3; - j < n_dofs_face * 3 + tangential_degree + 1; + for (unsigned int i = n_dofs_face * (2 * dim + degree); + i < n_dofs_face * (2 * dim + degree) + degree * (degree + 1); + ++i) + lexicographic_numbering.push_back(i + k_add * degree); + for (unsigned int j = n_dofs_face * 3; j < n_dofs_face * 3 + degree + 1; ++j) lexicographic_numbering.push_back(j + k_add); } @@ -116,12 +84,10 @@ PolynomialsRaviartThomas::get_lexicographic_numbering( { for (unsigned int i = 4 * n_dofs_face; i < 5 * n_dofs_face; ++i) lexicographic_numbering.push_back(i); - if (normal_degree > 1) - for (unsigned int i = - 6 * n_dofs_face + n_dofs_face * 2 * (normal_degree - 1); - i < 6 * n_dofs_face + n_dofs_face * 3 * (normal_degree - 1); - ++i) - lexicographic_numbering.push_back(i); + for (unsigned int i = 6 * n_dofs_face + n_dofs_face * 2 * degree; + i < 6 * n_dofs_face + n_dofs_face * 3 * degree; + ++i) + lexicographic_numbering.push_back(i); for (unsigned int i = 5 * n_dofs_face; i < 6 * n_dofs_face; ++i) lexicographic_numbering.push_back(i); } diff --git a/source/fe/fe_raviart_thomas.cc b/source/fe/fe_raviart_thomas.cc index 4071971966..7308a8ea93 100644 --- a/source/fe/fe_raviart_thomas.cc +++ b/source/fe/fe_raviart_thomas.cc @@ -43,17 +43,16 @@ DEAL_II_NAMESPACE_OPEN template FE_RaviartThomas::FE_RaviartThomas(const unsigned int deg) : FE_PolyTensor( - PolynomialsRaviartThomas(deg + 1, deg), + PolynomialsRaviartThomas(deg), FiniteElementData(get_dpo_vector(deg), dim, deg + 1, FiniteElementData::Hdiv), - std::vector(PolynomialsRaviartThomas::n_polynomials(deg + 1, - deg), + std::vector(PolynomialsRaviartThomas::n_polynomials(deg), true), - std::vector( - PolynomialsRaviartThomas::n_polynomials(deg + 1, deg), - ComponentMask(std::vector(dim, true)))) + std::vector(PolynomialsRaviartThomas::n_polynomials( + deg), + ComponentMask(std::vector(dim, true)))) { Assert(dim >= 2, ExcImpossibleInDim(dim)); const unsigned int n_dofs = this->n_dofs_per_cell(); diff --git a/source/fe/fe_raviart_thomas_nodal.cc b/source/fe/fe_raviart_thomas_nodal.cc index e3e5948ea1..3b61ab611a 100644 --- a/source/fe/fe_raviart_thomas_nodal.cc +++ b/source/fe/fe_raviart_thomas_nodal.cc @@ -66,16 +66,15 @@ namespace template FE_RaviartThomasNodal::FE_RaviartThomasNodal(const unsigned int degree) - : FE_PolyTensor( - PolynomialsRaviartThomas(degree + 1, degree), - FiniteElementData(get_rt_dpo_vector(dim, degree), - dim, - degree + 1, - FiniteElementData::Hdiv), - std::vector(1, false), - std::vector( - PolynomialsRaviartThomas::n_polynomials(degree + 1, degree), - ComponentMask(std::vector(dim, true)))) + : FE_PolyTensor(PolynomialsRaviartThomas(degree), + FiniteElementData(get_rt_dpo_vector(dim, degree), + dim, + degree + 1, + FiniteElementData::Hdiv), + std::vector(1, false), + std::vector( + PolynomialsRaviartThomas::n_polynomials(degree), + ComponentMask(std::vector(dim, true)))) { Assert(dim >= 2, ExcImpossibleInDim(dim)); @@ -84,8 +83,7 @@ FE_RaviartThomasNodal::FE_RaviartThomasNodal(const unsigned int degree) // First, initialize the generalized support points and quadrature weights, // since they are required for interpolation. this->generalized_support_points = - PolynomialsRaviartThomas(degree + 1, degree) - .get_polynomial_support_points(); + PolynomialsRaviartThomas(degree).get_polynomial_support_points(); AssertDimension(this->generalized_support_points.size(), this->n_dofs_per_cell()); -- 2.39.5