From 6006f779d2fff0a6fe72a5290548bae0622f39d3 Mon Sep 17 00:00:00 2001 From: grahambenharper Date: Sun, 1 Sep 2019 02:16:59 -0600 Subject: [PATCH] Change polynomial classes to derive from ScalarPolynomialsBase --- include/deal.II/base/polynomial_space.h | 69 +++++++------------ .../base/polynomials_rannacher_turek.h | 36 ++++++++-- .../deal.II/base/tensor_product_polynomials.h | 12 ++-- .../base/tensor_product_polynomials_bubbles.h | 12 ++-- .../base/tensor_product_polynomials_const.h | 12 ++-- include/deal.II/fe/fe_poly.h | 24 +++---- include/deal.II/fe/fe_poly_face.h | 12 ++-- source/base/polynomial_space.cc | 43 +++++++----- source/base/polynomials_bdm.cc | 12 ++-- source/base/polynomials_rannacher_turek.cc | 15 +++- source/base/tensor_product_polynomials.cc | 2 +- .../tensor_product_polynomials_bubbles.cc | 4 +- .../base/tensor_product_polynomials_const.cc | 4 +- source/fe/fe_bdm.cc | 12 ++-- source/fe/fe_dgp_nonparametric.cc | 36 +++++----- source/fe/mapping_q_generic.cc | 2 +- tests/base/anisotropic_1.cc | 2 +- tests/base/polynomial_test.cc | 2 +- 18 files changed, 168 insertions(+), 143 deletions(-) diff --git a/include/deal.II/base/polynomial_space.h b/include/deal.II/base/polynomial_space.h index 1e43d9a73a..368bda6d18 100644 --- a/include/deal.II/base/polynomial_space.h +++ b/include/deal.II/base/polynomial_space.h @@ -22,6 +22,7 @@ #include #include #include +#include #include #include @@ -95,7 +96,7 @@ DEAL_II_NAMESPACE_OPEN * 2005 */ template -class PolynomialSpace +class PolynomialSpace : public ScalarPolynomialsBase { public: /** @@ -142,12 +143,12 @@ public: * compute_grad_grad() functions, see below, in a loop over all polynomials. */ void - compute(const Point & unit_point, - std::vector & values, - std::vector> &grads, - std::vector> &grad_grads, - std::vector> &third_derivatives, - std::vector> &fourth_derivatives) const; + evaluate(const Point & unit_point, + std::vector & values, + std::vector> &grads, + std::vector> &grad_grads, + std::vector> &third_derivatives, + std::vector> &fourth_derivatives) const override; /** * Compute the value of the ith polynomial at unit point @@ -194,27 +195,20 @@ public: * given, then the result of this function is N in 1d, * N(N+1)/2 in 2d, and N(N+1)(N+2)/6 in 3d. */ - unsigned int - n() const; + static unsigned int + n_polynomials(const unsigned int n); /** - * Degree of the space. This is by definition the number of polynomials - * given to the constructor, NOT the maximal degree of a polynomial in this - * vector. The latter value is never checked and therefore left to the - * application. + * Return the name of the space, which is PolynomialSpace. */ - unsigned int - degree() const; + std::string + name() const override; /** - * Static function used in the constructor to compute the number of - * polynomials. - * - * @warning The argument `n` is not the maximal degree, but the number of - * onedimensional polynomials, thus the degree plus one. + * @copydoc ScalarPolynomialsBase::clone() */ - static unsigned int - n_polynomials(const unsigned int n); + virtual std::unique_ptr> + clone() const override; protected: /** @@ -234,11 +228,6 @@ private: */ const std::vector> polynomials; - /** - * Store the precomputed value which the n() function returns. - */ - const unsigned int n_pols; - /** * Index map for reordering the polynomials. */ @@ -270,16 +259,16 @@ PolynomialSpace<3>::compute_index(const unsigned int n) const; template template PolynomialSpace::PolynomialSpace(const std::vector &pols) - : polynomials(pols.begin(), pols.end()) - , n_pols(n_polynomials(polynomials.size())) - , index_map(n_pols) - , index_map_inverse(n_pols) + : ScalarPolynomialsBase(pols.size(), n_polynomials(pols.size())) + , polynomials(pols.begin(), pols.end()) + , index_map(n_polynomials(pols.size())) + , index_map_inverse(n_polynomials(pols.size())) { // per default set this index map // to identity. This map can be // changed by the user through the // set_numbering function - for (unsigned int i = 0; i < n_pols; ++i) + for (unsigned int i = 0; i < this->n(); ++i) { index_map[i] = i; index_map_inverse[i] = i; @@ -287,20 +276,12 @@ PolynomialSpace::PolynomialSpace(const std::vector &pols) } -template -inline unsigned int -PolynomialSpace::n() const -{ - return n_pols; -} - - template -inline unsigned int -PolynomialSpace::degree() const +inline std::string +PolynomialSpace::name() const { - return polynomials.size(); + return "PolynomialSpace"; } @@ -309,7 +290,7 @@ template void PolynomialSpace::output_indices(StreamType &out) const { - for (unsigned int i = 0; i < n_pols; ++i) + for (unsigned int i = 0; i < this->n(); ++i) { const std::array ix = compute_index(i); out << i << "\t"; diff --git a/include/deal.II/base/polynomials_rannacher_turek.h b/include/deal.II/base/polynomials_rannacher_turek.h index 37b8b3ac81..bf3fa41c46 100644 --- a/include/deal.II/base/polynomials_rannacher_turek.h +++ b/include/deal.II/base/polynomials_rannacher_turek.h @@ -18,6 +18,7 @@ #define dealii_polynomials_rannacher_turek_h #include +#include #include #include @@ -38,7 +39,7 @@ DEAL_II_NAMESPACE_OPEN * @date 2015 */ template -class PolynomialsRannacherTurek +class PolynomialsRannacherTurek : public ScalarPolynomialsBase { public: /** @@ -85,12 +86,24 @@ public: * zero. A size of zero means that we are not computing the vector entries. */ void - compute(const Point & unit_point, - std::vector & values, - std::vector> &grads, - std::vector> &grad_grads, - std::vector> &third_derivatives, - std::vector> &fourth_derivatives) const; + evaluate(const Point & unit_point, + std::vector & values, + std::vector> &grads, + std::vector> &grad_grads, + std::vector> &third_derivatives, + std::vector> &fourth_derivatives) const override; + + /** + * Return the name of the space, which is RannacherTurek. + */ + std::string + name() const override; + + /** + * @copydoc ScalarPolynomialsBase::clone() + */ + virtual std::unique_ptr> + clone() const override; }; @@ -204,6 +217,15 @@ PolynomialsRannacherTurek::compute_derivative(const unsigned int i, } + +template +inline std::string +PolynomialsRannacherTurek::name() const +{ + return "RannacherTurek"; +} + + DEAL_II_NAMESPACE_CLOSE #endif diff --git a/include/deal.II/base/tensor_product_polynomials.h b/include/deal.II/base/tensor_product_polynomials.h index 3a956a8b5e..9c39dec516 100644 --- a/include/deal.II/base/tensor_product_polynomials.h +++ b/include/deal.II/base/tensor_product_polynomials.h @@ -118,12 +118,12 @@ public: * over all tensor product polynomials. */ void - compute(const Point & unit_point, - std::vector & values, - std::vector> &grads, - std::vector> &grad_grads, - std::vector> &third_derivatives, - std::vector> &fourth_derivatives) const; + evaluate(const Point & unit_point, + std::vector & values, + std::vector> &grads, + std::vector> &grad_grads, + std::vector> &third_derivatives, + std::vector> &fourth_derivatives) const; /** * Compute the value of the ith tensor product polynomial at diff --git a/include/deal.II/base/tensor_product_polynomials_bubbles.h b/include/deal.II/base/tensor_product_polynomials_bubbles.h index 80559955c3..7a0f3dd018 100644 --- a/include/deal.II/base/tensor_product_polynomials_bubbles.h +++ b/include/deal.II/base/tensor_product_polynomials_bubbles.h @@ -69,12 +69,12 @@ public: * over all tensor product polynomials. */ void - compute(const Point & unit_point, - std::vector & values, - std::vector> &grads, - std::vector> &grad_grads, - std::vector> &third_derivatives, - std::vector> &fourth_derivatives) const; + evaluate(const Point & unit_point, + std::vector & values, + std::vector> &grads, + std::vector> &grad_grads, + std::vector> &third_derivatives, + std::vector> &fourth_derivatives) const; /** * Compute the value of the ith tensor product polynomial at diff --git a/include/deal.II/base/tensor_product_polynomials_const.h b/include/deal.II/base/tensor_product_polynomials_const.h index ff2680bd45..1671532506 100644 --- a/include/deal.II/base/tensor_product_polynomials_const.h +++ b/include/deal.II/base/tensor_product_polynomials_const.h @@ -75,12 +75,12 @@ public: * over all tensor product polynomials. */ void - compute(const Point & unit_point, - std::vector & values, - std::vector> &grads, - std::vector> &grad_grads, - std::vector> &third_derivatives, - std::vector> &fourth_derivatives) const; + evaluate(const Point & unit_point, + std::vector & values, + std::vector> &grads, + std::vector> &grad_grads, + std::vector> &third_derivatives, + std::vector> &fourth_derivatives) const; /** * Compute the value of the ith tensor product polynomial at diff --git a/include/deal.II/fe/fe_poly.h b/include/deal.II/fe/fe_poly.h index 57a4419522..d3917724d4 100644 --- a/include/deal.II/fe/fe_poly.h +++ b/include/deal.II/fe/fe_poly.h @@ -40,12 +40,12 @@ DEAL_II_NAMESPACE_OPEN * @code * static const unsigned int dimension; * - * void compute (const Point &unit_point, - * std::vector &values, - * std::vector > &grads, - * std::vector > &grad_grads, - * std::vector > &third_derivatives, - * std::vector > &fourth_derivatives) const; + * void evaluate (const Point &unit_point, + * std::vector &values, + * std::vector > &grads, + * std::vector > &grad_grads, + * std::vector > &third_derivatives, + * std::vector > &fourth_derivatives) const; * * double compute_value (const unsigned int i, * const Point &p) const; @@ -303,12 +303,12 @@ protected: update_3rd_derivatives)) for (unsigned int i = 0; i < n_q_points; ++i) { - poly_space.compute(quadrature.point(i), - values, - grads, - grad_grads, - third_derivatives, - fourth_derivatives); + poly_space.evaluate(quadrature.point(i), + values, + grads, + grad_grads, + third_derivatives, + fourth_derivatives); // the values of shape functions at quadrature points don't change. // consequently, write these values right into the output array if diff --git a/include/deal.II/fe/fe_poly_face.h b/include/deal.II/fe/fe_poly_face.h index 4a79199f40..1c5780fcf9 100644 --- a/include/deal.II/fe/fe_poly_face.h +++ b/include/deal.II/fe/fe_poly_face.h @@ -133,12 +133,12 @@ protected: std::vector(n_q_points)); for (unsigned int i = 0; i < n_q_points; ++i) { - poly_space.compute(quadrature.point(i), - values, - grads, - grad_grads, - empty_vector_of_3rd_order_tensors, - empty_vector_of_4th_order_tensors); + poly_space.evaluate(quadrature.point(i), + values, + grads, + grad_grads, + empty_vector_of_3rd_order_tensors, + empty_vector_of_4th_order_tensors); for (unsigned int k = 0; k < poly_space.n(); ++k) data.shape_values[k][i] = values[k]; diff --git a/source/base/polynomial_space.cc b/source/base/polynomial_space.cc index a87fcefa22..95f87b2be3 100644 --- a/source/base/polynomial_space.cc +++ b/source/base/polynomial_space.cc @@ -15,6 +15,7 @@ #include #include +#include #include DEAL_II_NAMESPACE_OPEN @@ -200,7 +201,7 @@ PolynomialSpace::compute_grad_grad(const unsigned int i, template void -PolynomialSpace::compute( +PolynomialSpace::evaluate( const Point & p, std::vector & values, std::vector> &grads, @@ -210,41 +211,42 @@ PolynomialSpace::compute( { const unsigned int n_1d = polynomials.size(); - Assert(values.size() == n_pols || values.size() == 0, - ExcDimensionMismatch2(values.size(), n_pols, 0)); - Assert(grads.size() == n_pols || grads.size() == 0, - ExcDimensionMismatch2(grads.size(), n_pols, 0)); - Assert(grad_grads.size() == n_pols || grad_grads.size() == 0, - ExcDimensionMismatch2(grad_grads.size(), n_pols, 0)); - Assert(third_derivatives.size() == n_pols || third_derivatives.size() == 0, - ExcDimensionMismatch2(third_derivatives.size(), n_pols, 0)); - Assert(fourth_derivatives.size() == n_pols || fourth_derivatives.size() == 0, - ExcDimensionMismatch2(fourth_derivatives.size(), n_pols, 0)); + Assert(values.size() == this->n() || values.size() == 0, + ExcDimensionMismatch2(values.size(), this->n(), 0)); + Assert(grads.size() == this->n() || grads.size() == 0, + ExcDimensionMismatch2(grads.size(), this->n(), 0)); + Assert(grad_grads.size() == this->n() || grad_grads.size() == 0, + ExcDimensionMismatch2(grad_grads.size(), this->n(), 0)); + Assert(third_derivatives.size() == this->n() || third_derivatives.size() == 0, + ExcDimensionMismatch2(third_derivatives.size(), this->n(), 0)); + Assert(fourth_derivatives.size() == this->n() || + fourth_derivatives.size() == 0, + ExcDimensionMismatch2(fourth_derivatives.size(), this->n(), 0)); unsigned int v_size = 0; bool update_values = false, update_grads = false, update_grad_grads = false; bool update_3rd_derivatives = false, update_4th_derivatives = false; - if (values.size() == n_pols) + if (values.size() == this->n()) { update_values = true; v_size = 1; } - if (grads.size() == n_pols) + if (grads.size() == this->n()) { update_grads = true; v_size = 2; } - if (grad_grads.size() == n_pols) + if (grad_grads.size() == this->n()) { update_grad_grads = true; v_size = 3; } - if (third_derivatives.size() == n_pols) + if (third_derivatives.size() == this->n()) { update_3rd_derivatives = true; v_size = 4; } - if (fourth_derivatives.size() == n_pols) + if (fourth_derivatives.size() == this->n()) { update_4th_derivatives = true; v_size = 5; @@ -396,6 +398,15 @@ PolynomialSpace::compute( } + +template +std::unique_ptr> +PolynomialSpace::clone() const +{ + return std_cxx14::make_unique>(*this); +} + + template class PolynomialSpace<1>; template class PolynomialSpace<2>; template class PolynomialSpace<3>; diff --git a/source/base/polynomials_bdm.cc b/source/base/polynomials_bdm.cc index 80a2abe0c1..34005fcede 100644 --- a/source/base/polynomials_bdm.cc +++ b/source/base/polynomials_bdm.cc @@ -98,12 +98,12 @@ PolynomialsBDM::evaluate( // will have first all polynomials // in the x-component, then y and // z. - polynomial_space.compute(unit_point, - p_values, - p_grads, - p_grad_grads, - p_third_derivatives, - p_fourth_derivatives); + polynomial_space.evaluate(unit_point, + p_values, + p_grads, + p_grad_grads, + p_third_derivatives, + p_fourth_derivatives); std::fill(values.begin(), values.end(), Tensor<1, dim>()); for (unsigned int i = 0; i < p_values.size(); ++i) diff --git a/source/base/polynomials_rannacher_turek.cc b/source/base/polynomials_rannacher_turek.cc index a896143dab..952c58e41a 100644 --- a/source/base/polynomials_rannacher_turek.cc +++ b/source/base/polynomials_rannacher_turek.cc @@ -16,12 +16,14 @@ #include #include +#include DEAL_II_NAMESPACE_OPEN template PolynomialsRannacherTurek::PolynomialsRannacherTurek() + : ScalarPolynomialsBase(2, dealii::GeometryInfo::faces_per_cell) { Assert(dim == 2, ExcNotImplemented()); } @@ -141,7 +143,7 @@ PolynomialsRannacherTurek::compute_grad_grad( template void -PolynomialsRannacherTurek::compute( +PolynomialsRannacherTurek::evaluate( const Point & unit_point, std::vector & values, std::vector> &grads, @@ -149,7 +151,7 @@ PolynomialsRannacherTurek::compute( std::vector> &third_derivatives, std::vector> &fourth_derivatives) const { - const unsigned int n_pols = dealii::GeometryInfo::faces_per_cell; + const unsigned int n_pols = this->n(); Assert(values.size() == n_pols || values.size() == 0, ExcDimensionMismatch(values.size(), n_pols)); Assert(grads.size() == n_pols || grads.size() == 0, @@ -187,6 +189,15 @@ PolynomialsRannacherTurek::compute( } + +template +std::unique_ptr> +PolynomialsRannacherTurek::clone() const +{ + return std_cxx14::make_unique>(*this); +} + + // explicit instantiations #include "polynomials_rannacher_turek.inst" diff --git a/source/base/tensor_product_polynomials.cc b/source/base/tensor_product_polynomials.cc index 71efbce20e..71d5f38186 100644 --- a/source/base/tensor_product_polynomials.cc +++ b/source/base/tensor_product_polynomials.cc @@ -241,7 +241,7 @@ TensorProductPolynomials::compute_grad_grad( template void -TensorProductPolynomials::compute( +TensorProductPolynomials::evaluate( const Point & p, std::vector & values, std::vector> &grads, diff --git a/source/base/tensor_product_polynomials_bubbles.cc b/source/base/tensor_product_polynomials_bubbles.cc index e927914d31..afd88bfb82 100644 --- a/source/base/tensor_product_polynomials_bubbles.cc +++ b/source/base/tensor_product_polynomials_bubbles.cc @@ -215,7 +215,7 @@ TensorProductPolynomialsBubbles::compute_grad_grad( template void -TensorProductPolynomialsBubbles::compute( +TensorProductPolynomialsBubbles::evaluate( const Point & p, std::vector & values, std::vector> &grads, @@ -273,7 +273,7 @@ TensorProductPolynomialsBubbles::compute( do_4th_derivatives = true; } - this->TensorProductPolynomials::compute( + this->TensorProductPolynomials::evaluate( p, values, grads, grad_grads, third_derivatives, fourth_derivatives); for (unsigned int i = this->n_tensor_pols; diff --git a/source/base/tensor_product_polynomials_const.cc b/source/base/tensor_product_polynomials_const.cc index bbeef38690..63bce45b98 100644 --- a/source/base/tensor_product_polynomials_const.cc +++ b/source/base/tensor_product_polynomials_const.cc @@ -86,7 +86,7 @@ TensorProductPolynomialsConst::compute_grad_grad(const unsigned int i, template void -TensorProductPolynomialsConst::compute( +TensorProductPolynomialsConst::evaluate( const Point & p, std::vector & values, std::vector> &grads, @@ -141,7 +141,7 @@ TensorProductPolynomialsConst::compute( do_4th_derivatives = true; } - this->TensorProductPolynomials::compute( + this->TensorProductPolynomials::evaluate( p, values, grads, grad_grads, third_derivatives, fourth_derivatives); // for dgq node: values =1, grads=0, grads_grads=0, third_derivatives=0, diff --git a/source/fe/fe_bdm.cc b/source/fe/fe_bdm.cc index 0a09f141f8..26b747877b 100644 --- a/source/fe/fe_bdm.cc +++ b/source/fe/fe_bdm.cc @@ -294,12 +294,12 @@ namespace internal for (unsigned int k = 0; k < quadrature.size(); ++k) { test_values[k].resize(poly.n()); - poly.compute(quadrature.point(k), - test_values[k], - dummy1, - dummy2, - dummy3, - dummy4); + poly.evaluate(quadrature.point(k), + test_values[k], + dummy1, + dummy2, + dummy3, + dummy4); for (unsigned int i = 0; i < poly.n(); ++i) { test_values[k][i] *= quadrature.weight(k); diff --git a/source/fe/fe_dgp_nonparametric.cc b/source/fe/fe_dgp_nonparametric.cc index 9c3bc1b760..f748a4ed03 100644 --- a/source/fe/fe_dgp_nonparametric.cc +++ b/source/fe/fe_dgp_nonparametric.cc @@ -332,12 +332,12 @@ FE_DGPNonparametric::fill_fe_values( if (fe_internal.update_each & (update_values | update_gradients)) for (unsigned int i = 0; i < n_q_points; ++i) { - polynomial_space.compute(mapping_data.quadrature_points[i], - values, - grads, - grad_grads, - empty_vector_of_3rd_order_tensors, - empty_vector_of_4th_order_tensors); + polynomial_space.evaluate(mapping_data.quadrature_points[i], + values, + grads, + grad_grads, + empty_vector_of_3rd_order_tensors, + empty_vector_of_4th_order_tensors); if (fe_internal.update_each & update_values) for (unsigned int k = 0; k < this->dofs_per_cell; ++k) @@ -388,12 +388,12 @@ FE_DGPNonparametric::fill_fe_face_values( if (fe_internal.update_each & (update_values | update_gradients)) for (unsigned int i = 0; i < n_q_points; ++i) { - polynomial_space.compute(mapping_data.quadrature_points[i], - values, - grads, - grad_grads, - empty_vector_of_3rd_order_tensors, - empty_vector_of_4th_order_tensors); + polynomial_space.evaluate(mapping_data.quadrature_points[i], + values, + grads, + grad_grads, + empty_vector_of_3rd_order_tensors, + empty_vector_of_4th_order_tensors); if (fe_internal.update_each & update_values) for (unsigned int k = 0; k < this->dofs_per_cell; ++k) @@ -445,12 +445,12 @@ FE_DGPNonparametric::fill_fe_subface_values( if (fe_internal.update_each & (update_values | update_gradients)) for (unsigned int i = 0; i < n_q_points; ++i) { - polynomial_space.compute(mapping_data.quadrature_points[i], - values, - grads, - grad_grads, - empty_vector_of_3rd_order_tensors, - empty_vector_of_4th_order_tensors); + polynomial_space.evaluate(mapping_data.quadrature_points[i], + values, + grads, + grad_grads, + empty_vector_of_3rd_order_tensors, + empty_vector_of_4th_order_tensors); if (fe_internal.update_each & update_values) for (unsigned int k = 0; k < this->dofs_per_cell; ++k) diff --git a/source/fe/mapping_q_generic.cc b/source/fe/mapping_q_generic.cc index 4c118116d3..b4e5d22f8b 100644 --- a/source/fe/mapping_q_generic.cc +++ b/source/fe/mapping_q_generic.cc @@ -296,7 +296,7 @@ namespace internal data.shape_fourth_derivatives.size() != 0) for (unsigned int point = 0; point < n_points; ++point) { - tensor_pols.compute( + tensor_pols.evaluate( unit_points[point], values, grads, grad2, grad3, grad4); if (data.shape_values.size() != 0) diff --git a/tests/base/anisotropic_1.cc b/tests/base/anisotropic_1.cc index 66cce181ef..e6f3327a91 100644 --- a/tests/base/anisotropic_1.cc +++ b/tests/base/anisotropic_1.cc @@ -42,7 +42,7 @@ check_poly(const Point & x, std::vector> third1(n), third2(n); std::vector> fourth1(n), fourth2(n); - p.compute(x, values1, gradients1, second1, third1, fourth1); + p.evaluate(x, values1, gradients1, second1, third1, fourth1); q.compute(x, values2, gradients2, second2, third2, fourth2); for (unsigned int k = 0; k < n; ++k) diff --git a/tests/base/polynomial_test.cc b/tests/base/polynomial_test.cc index 2e18d9fe37..7229f64a37 100644 --- a/tests/base/polynomial_test.cc +++ b/tests/base/polynomial_test.cc @@ -41,7 +41,7 @@ check_poly(const Point &x, const PolynomialType &p) std::vector> third(n); std::vector> fourth(n); - p.compute(x, values, gradients, second, third, fourth); + p.evaluate(x, values, gradients, second, third, fourth); for (unsigned int k = 0; k < n; ++k) { -- 2.39.5