From c5d358ab7e97a61858ef5d9eb19c8a5428a5b0df Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Tue, 2 Feb 2021 21:08:41 -0700 Subject: [PATCH] Make get_gauss_type_quadrature() and get_nodal_type_quadrature() members of ReferenceCell::Type. --- include/deal.II/fe/fe_tools.templates.h | 11 +++---- include/deal.II/grid/reference_cell.h | 43 +++++++++++++------------ source/fe/mapping_fe_field.cc | 20 +++++++----- source/grid/reference_cell.cc | 43 ++++++++++++------------- source/grid/reference_cell.inst.in | 8 ++--- tests/simplex/fe_p_bubbles_02.cc | 2 +- 6 files changed, 64 insertions(+), 63 deletions(-) diff --git a/include/deal.II/fe/fe_tools.templates.h b/include/deal.II/fe/fe_tools.templates.h index b2142de629..efdea4bd2f 100644 --- a/include/deal.II/fe/fe_tools.templates.h +++ b/include/deal.II/fe/fe_tools.templates.h @@ -1598,7 +1598,7 @@ namespace FETools const unsigned int n1 = fe1.n_dofs_per_cell(); const unsigned int n2 = fe2.n_dofs_per_cell(); - const auto reference_cell_type = fe1.reference_cell_type(); + const ReferenceCell::Type reference_cell_type = fe1.reference_cell_type(); Assert(fe1.reference_cell_type() == fe2.reference_cell_type(), ExcNotImplemented()); @@ -1615,8 +1615,7 @@ namespace FETools std::max(fe1.tensor_degree(), fe2.tensor_degree()); Assert(degree != numbers::invalid_unsigned_int, ExcNotImplemented()); const auto quadrature = - ReferenceCell::get_gauss_type_quadrature(reference_cell_type, - degree + 1); + reference_cell_type.get_gauss_type_quadrature(degree + 1); // Set up FEValues. const UpdateFlags flags = @@ -1809,7 +1808,8 @@ namespace FETools ExcDimensionMismatch(matrices[i].m(), n)); } - const auto reference_cell_type = fe.reference_cell_type(); + const ReferenceCell::Type reference_cell_type = + fe.reference_cell_type(); // Set up meshes, one with a single // reference cell and refine it once @@ -1824,8 +1824,7 @@ namespace FETools reference_cell_type .template get_default_linear_mapping(); const auto &q_fine = - ReferenceCell::get_gauss_type_quadrature(reference_cell_type, - degree + 1); + reference_cell_type.get_gauss_type_quadrature(degree + 1); const unsigned int nq = q_fine.size(); FEValues fine(mapping, diff --git a/include/deal.II/grid/reference_cell.h b/include/deal.II/grid/reference_cell.h index 892fec8364..84b587b9b3 100644 --- a/include/deal.II/grid/reference_cell.h +++ b/include/deal.II/grid/reference_cell.h @@ -208,6 +208,28 @@ namespace ReferenceCell const Mapping & get_default_linear_mapping() const; + /** + * Return a Gauss-type quadrature matching the given reference cell (QGauss, + * Simplex::QGauss, Simplex::QGaussPyramid, Simplex::QGaussWedge). + * + * @param[in] n_points_1D The number of quadrature points in each direction + * (QGauss) or an indication of what polynomial degree needs to be + * integrated exactly for the other types. + */ + template + Quadrature + get_gauss_type_quadrature(const unsigned n_points_1D) const; + + /** + * Return a quadrature rule with the support points of the given reference + * cell. + * + * @note The weights of the quadrature object are left unfilled. + */ + template + const Quadrature & + get_nodal_type_quadrature() const; + /** * Return a text representation of the reference cell represented by the * current object. @@ -728,27 +750,6 @@ namespace ReferenceCell const Mapping & get_default_linear_mapping(const Triangulation &triangulation); - /** - * Return a Gauss-type quadrature matching the given reference cell(QGauss, - * Simplex::QGauss, Simplex::QGaussPyramid, Simplex::QGaussWedge) and - * @p n_points_1D the number of quadrature points in each direction (QGuass) - * or the indication of what polynomial degree to be integrated exactly. - */ - template - Quadrature - get_gauss_type_quadrature(const Type & reference_cell, - const unsigned n_points_1D); - - /** - * Return a quadrature rule with the support points of the given reference - * cell. - * - * @note The weights are not filled. - */ - template - const Quadrature & - get_nodal_type_quadrature(const Type &reference_cell); - namespace internal { /** diff --git a/source/fe/mapping_fe_field.cc b/source/fe/mapping_fe_field.cc index 3b9712365a..91e4c9ea80 100644 --- a/source/fe/mapping_fe_field.cc +++ b/source/fe/mapping_fe_field.cc @@ -207,8 +207,9 @@ MappingFEField::MappingFEField( true)) , fe_to_real(fe_mask.size(), numbers::invalid_unsigned_int) , fe_values(this->euler_dof_handler->get_fe(), - ReferenceCell::get_nodal_type_quadrature( - this->euler_dof_handler->get_fe().reference_cell_type()), + this->euler_dof_handler->get_fe() + .reference_cell_type() + .template get_nodal_type_quadrature(), update_values) { unsigned int size = 0; @@ -236,8 +237,9 @@ MappingFEField::MappingFEField( true)) , fe_to_real(fe_mask.size(), numbers::invalid_unsigned_int) , fe_values(this->euler_dof_handler->get_fe(), - ReferenceCell::get_nodal_type_quadrature( - this->euler_dof_handler->get_fe().reference_cell_type()), + this->euler_dof_handler->get_fe() + .reference_cell_type() + .template get_nodal_type_quadrature(), update_values) { unsigned int size = 0; @@ -276,8 +278,9 @@ MappingFEField::MappingFEField( true)) , fe_to_real(fe_mask.size(), numbers::invalid_unsigned_int) , fe_values(this->euler_dof_handler->get_fe(), - ReferenceCell::get_nodal_type_quadrature( - this->euler_dof_handler->get_fe().reference_cell_type()), + this->euler_dof_handler->get_fe() + .reference_cell_type() + .template get_nodal_type_quadrature(), update_values) { unsigned int size = 0; @@ -313,8 +316,9 @@ MappingFEField::MappingFEField( , fe_mask(mapping.fe_mask) , fe_to_real(mapping.fe_to_real) , fe_values(mapping.euler_dof_handler->get_fe(), - ReferenceCell::get_nodal_type_quadrature( - this->euler_dof_handler->get_fe().reference_cell_type()), + this->euler_dof_handler->get_fe() + .reference_cell_type() + .template get_nodal_type_quadrature(), update_values) {} diff --git a/source/grid/reference_cell.cc b/source/grid/reference_cell.cc index c187cdd4ea..ef0d1afcdc 100644 --- a/source/grid/reference_cell.cc +++ b/source/grid/reference_cell.cc @@ -179,18 +179,17 @@ namespace ReferenceCell template Quadrature - get_gauss_type_quadrature(const Type & reference_cell, - const unsigned n_points_1D) + Type::get_gauss_type_quadrature(const unsigned n_points_1D) const { - AssertDimension(dim, reference_cell.get_dimension()); + AssertDimension(dim, get_dimension()); - if (reference_cell == Type::get_hypercube()) + if (is_hyper_cube()) return QGauss(n_points_1D); - else if (reference_cell == Type::Tri || reference_cell == Type::Tet) + else if (is_simplex()) return Simplex::QGauss(n_points_1D); - else if (reference_cell == Type::Pyramid) + else if (*this == Type::Pyramid) return Simplex::QGaussPyramid(n_points_1D); - else if (reference_cell == Type::Wedge) + else if (*this == Type::Wedge) return Simplex::QGaussWedge(n_points_1D); else Assert(false, ExcNotImplemented()); @@ -202,10 +201,13 @@ namespace ReferenceCell template const Quadrature & - get_nodal_type_quadrature(const Type &reference_cell) + Type::get_nodal_type_quadrature() const { - AssertDimension(dim, reference_cell.get_dimension()); + AssertDimension(dim, get_dimension()); + // A function that is used to fill a quadrature object of the + // desired type the first time we encounter a particular + // reference cell const auto create_quadrature = [](const Type &reference_cell) { Triangulation tria; GridGenerator::reference_cell(reference_cell, tria); @@ -213,35 +215,30 @@ namespace ReferenceCell return Quadrature(tria.get_vertices()); }; - if (reference_cell == Type::get_hypercube()) + if (is_hyper_cube()) { - static const Quadrature quadrature = - create_quadrature(reference_cell); + static const Quadrature quadrature = create_quadrature(*this); return quadrature; } - else if (reference_cell == Type::Tri || reference_cell == Type::Tet) + else if (is_simplex()) { - static const Quadrature quadrature = - create_quadrature(reference_cell); + static const Quadrature quadrature = create_quadrature(*this); return quadrature; } - else if (reference_cell == Type::Pyramid) + else if (*this == Type::Pyramid) { - static const Quadrature quadrature = - create_quadrature(reference_cell); + static const Quadrature quadrature = create_quadrature(*this); return quadrature; } - else if (reference_cell == Type::Wedge) + else if (*this == Type::Wedge) { - static const Quadrature quadrature = - create_quadrature(reference_cell); + static const Quadrature quadrature = create_quadrature(*this); return quadrature; } else Assert(false, ExcNotImplemented()); - static Quadrature dummy; - + static const Quadrature dummy; return dummy; // never reached } diff --git a/source/grid/reference_cell.inst.in b/source/grid/reference_cell.inst.in index 418844ce49..be85dfd4f1 100644 --- a/source/grid/reference_cell.inst.in +++ b/source/grid/reference_cell.inst.in @@ -32,8 +32,8 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) for (deal_II_dimension : DIMENSIONS) { - template Quadrature get_gauss_type_quadrature( - const Type &reference_cell, const unsigned n_points_1D); - template const Quadrature &get_nodal_type_quadrature( - const Type &reference_cell); + template Quadrature Type::get_gauss_type_quadrature( + const unsigned n_points_1D) const; + template const Quadrature + &Type::get_nodal_type_quadrature() const; } diff --git a/tests/simplex/fe_p_bubbles_02.cc b/tests/simplex/fe_p_bubbles_02.cc index 9a45e95ed8..3cac6f262f 100644 --- a/tests/simplex/fe_p_bubbles_02.cc +++ b/tests/simplex/fe_p_bubbles_02.cc @@ -47,7 +47,7 @@ compute_nodal_quadrature(const FiniteElement &fe) const ReferenceCell::Type type = fe.reference_cell_type(); const Quadrature q_gauss = - ReferenceCell::get_gauss_type_quadrature(type, fe.tensor_degree() + 1); + type.get_gauss_type_quadrature(fe.tensor_degree() + 1); Triangulation tria; GridGenerator::reference_cell(type, tria); const Mapping &mapping = -- 2.39.5