From 059f3bbf2ef408226d5ff7d0d128e7fd0afb1f57 Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Fri, 9 Feb 2018 18:48:56 +0100 Subject: [PATCH] Switch FE_DGQHermite to the more appropriate polynomial basis. --- include/deal.II/fe/fe_dgq.h | 25 +++++++++++++++---------- source/fe/fe_dgq.cc | 16 +++++----------- tests/fe/element_constant_modes.cc | 1 - tests/fe/element_constant_modes.output | 6 ------ 4 files changed, 20 insertions(+), 28 deletions(-) diff --git a/include/deal.II/fe/fe_dgq.h b/include/deal.II/fe/fe_dgq.h index e8c832e183..cf7ff20b67 100644 --- a/include/deal.II/fe/fe_dgq.h +++ b/include/deal.II/fe/fe_dgq.h @@ -481,16 +481,20 @@ public: /** * Implementation of scalar, discontinuous tensor product elements based on - * Hermite polynomials, described by the polynomial space - * Polynomials::HermiteInterpolation. As opposed to the basic FE_DGQ element, + * Hermite-like polynomials, described by the polynomial space + * Polynomials::HermiteLikeInterpolation. As opposed to the basic FE_DGQ element, * these elements are not interpolatory and no support points are defined. * - * This element is only a Hermite polynomials for degrees larger or equal to - * three. For degrees zero to two, a usual Lagrange basis is selected. + * Note that Hermite polynomials are only available for degrees larger or + * equal to three, and thus the beneficial properties of + * Polynomials::HermiteLikeInterpolation with only two basis functions having + * a non-trivial value or derivative on a face per dimension is only present + * for higher degrees. To facilitate usage also for degrees zero to two, a + * usual Lagrange basis is constructed by this class. * * See the base class documentation in FE_DGQ for details. * - * @author Martin Kronbichler, 2017 + * @author Martin Kronbichler, 2017, 2018 */ template class FE_DGQHermite : public FE_DGQ @@ -498,15 +502,16 @@ class FE_DGQHermite : public FE_DGQ public: /** * Constructor for tensor product polynomials based on - * Polynomials::HermiteInterpolation. + * Polynomials::HermiteLikeInterpolation. */ FE_DGQHermite (const unsigned int degree); /** - * Return a list of constant modes of the element. For the Hermite basis of - * degree three and larger, it returns one row where the first two elements - * (corresponding to the value left and the value at the right) are set to - * true and all other elements are set to false. + * Return a list of constant modes of the element. Due to the special + * construction of the Hermite basis, there is no simple representation of + * the constant value 1 in terms of a 0/1 selection, so get_constant_modes() + * throws an exception in case it is called for degrees larger or equal to + * three. */ virtual std::pair, std::vector > get_constant_modes () const; diff --git a/source/fe/fe_dgq.cc b/source/fe/fe_dgq.cc index 31189df188..341f3f741c 100644 --- a/source/fe/fe_dgq.cc +++ b/source/fe/fe_dgq.cc @@ -882,7 +882,7 @@ FE_DGQHermite::FE_DGQHermite (const unsigned int degree) Polynomials::generate_complete_Lagrange_basis (internal::FE_DGQ::get_QGaussLobatto_points(degree)) : - Polynomials::HermiteInterpolation::generate_complete_basis(degree)) + Polynomials::HermiteLikeInterpolation::generate_complete_basis(degree)) {} @@ -895,17 +895,11 @@ FE_DGQHermite::get_constant_modes () const return this->FE_DGQ::get_constant_modes(); else { - // The first two basis functions in the Hermite polynomials represent - // the value 1 in the left and right end point of the element. Expand - // them into the tensor product. - AssertThrow(dim<=3, ExcNotImplemented()); - Table<2,bool> constant_modes(1, this->dofs_per_cell); - for (unsigned int i=0; i<(dim>2?2:1); ++i) - for (unsigned int j=0; j<(dim>1?2:1); ++j) - for (unsigned int k=0; k<2; ++k) - constant_modes(0,i*(this->degree+1)*(this->degree+1)+j*(this->degree+1)+k) = true; + AssertThrow(false, + ExcMessage("Constant mode cannot be represented by 0/1 vector")); return std::pair, std::vector > - (constant_modes, std::vector(1, 0)); + (Table<2,bool>(1, this->dofs_per_cell), + std::vector(1, 0)); } } diff --git a/tests/fe/element_constant_modes.cc b/tests/fe/element_constant_modes.cc index 14be190e6f..4fe22916e4 100644 --- a/tests/fe/element_constant_modes.cc +++ b/tests/fe/element_constant_modes.cc @@ -50,7 +50,6 @@ void test() print_constant_modes(FE_DGQ(1)); print_constant_modes(FE_DGQLegendre(2)); print_constant_modes(FE_DGQHermite(2)); - print_constant_modes(FE_DGQHermite(3)); print_constant_modes(FE_DGP(2)); print_constant_modes(FE_Q_Hierarchical(1)); print_constant_modes(FE_Q_Hierarchical(2)); diff --git a/tests/fe/element_constant_modes.output b/tests/fe/element_constant_modes.output index f25e1afb08..ee07412c5f 100644 --- a/tests/fe/element_constant_modes.output +++ b/tests/fe/element_constant_modes.output @@ -26,9 +26,6 @@ DEAL:: DEAL::Testing FE_DGQHermite<2>(2) DEAL::1 1 1 1 1 1 1 1 1 DEAL:: -DEAL::Testing FE_DGQHermite<2>(3) -DEAL::1 1 0 0 1 1 0 0 0 0 0 0 0 0 0 0 -DEAL:: DEAL::Testing FE_DGP<2>(2) DEAL::1 0 0 0 0 0 DEAL:: @@ -86,9 +83,6 @@ DEAL:: DEAL::Testing FE_DGQHermite<3>(2) DEAL::1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 DEAL:: -DEAL::Testing FE_DGQHermite<3>(3) -DEAL::1 1 0 0 1 1 0 0 0 0 0 0 0 0 0 0 1 1 0 0 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -DEAL:: DEAL::Testing FE_DGP<3>(2) DEAL::1 0 0 0 0 0 0 0 0 0 DEAL:: -- 2.39.5