From 3e72aa3436648de7d2ed6a1f44633fd0e517f121 Mon Sep 17 00:00:00 2001 From: grahambenharper Date: Sat, 17 Aug 2019 22:06:37 -0600 Subject: [PATCH] Remove template argument from FE_PolyTensor --- include/deal.II/fe/fe_abf.h | 2 +- include/deal.II/fe/fe_bdm.h | 2 +- include/deal.II/fe/fe_bernardi_raugel.h | 3 +- include/deal.II/fe/fe_dg_vector.h | 2 +- include/deal.II/fe/fe_dg_vector.templates.h | 6 +- include/deal.II/fe/fe_nedelec.h | 2 +- include/deal.II/fe/fe_poly_tensor.h | 28 ++++--- include/deal.II/fe/fe_raviart_thomas.h | 6 +- include/deal.II/fe/fe_rt_bubbles.h | 2 +- source/fe/fe_abf.cc | 4 +- source/fe/fe_bdm.cc | 4 +- source/fe/fe_bernardi_raugel.cc | 4 +- source/fe/fe_nedelec.cc | 4 +- source/fe/fe_poly_tensor.cc | 88 +++++++++++---------- source/fe/fe_poly_tensor.inst.in | 13 +-- source/fe/fe_raviart_thomas.cc | 4 +- source/fe/fe_raviart_thomas_nodal.cc | 19 +++-- source/fe/fe_rt_bubbles.cc | 19 +++-- 18 files changed, 103 insertions(+), 109 deletions(-) diff --git a/include/deal.II/fe/fe_abf.h b/include/deal.II/fe/fe_abf.h index 25c1fae511..3d26f0a9a3 100644 --- a/include/deal.II/fe/fe_abf.h +++ b/include/deal.II/fe/fe_abf.h @@ -100,7 +100,7 @@ DEAL_II_NAMESPACE_OPEN * Kanschat and Wolfgang Bangerth */ template -class FE_ABF : public FE_PolyTensor, dim> +class FE_ABF : public FE_PolyTensor { public: /** diff --git a/include/deal.II/fe/fe_bdm.h b/include/deal.II/fe/fe_bdm.h index 498e6a73be..ccec16b47f 100644 --- a/include/deal.II/fe/fe_bdm.h +++ b/include/deal.II/fe/fe_bdm.h @@ -56,7 +56,7 @@ DEAL_II_NAMESPACE_OPEN * @ingroup fe */ template -class FE_BDM : public FE_PolyTensor, dim> +class FE_BDM : public FE_PolyTensor { public: /** diff --git a/include/deal.II/fe/fe_bernardi_raugel.h b/include/deal.II/fe/fe_bernardi_raugel.h index ca3c94c8b9..6a4349f31e 100644 --- a/include/deal.II/fe/fe_bernardi_raugel.h +++ b/include/deal.II/fe/fe_bernardi_raugel.h @@ -84,8 +84,7 @@ DEAL_II_NAMESPACE_OPEN * */ template -class FE_BernardiRaugel - : public FE_PolyTensor, dim> +class FE_BernardiRaugel : public FE_PolyTensor { public: /** diff --git a/include/deal.II/fe/fe_dg_vector.h b/include/deal.II/fe/fe_dg_vector.h index 9e0ff6eee7..f2c871d664 100644 --- a/include/deal.II/fe/fe_dg_vector.h +++ b/include/deal.II/fe/fe_dg_vector.h @@ -52,7 +52,7 @@ DEAL_II_NAMESPACE_OPEN * @date 2010 */ template -class FE_DGVector : public FE_PolyTensor +class FE_DGVector : public FE_PolyTensor { public: /** diff --git a/include/deal.II/fe/fe_dg_vector.templates.h b/include/deal.II/fe/fe_dg_vector.templates.h index 05c0622a99..edba144969 100644 --- a/include/deal.II/fe/fe_dg_vector.templates.h +++ b/include/deal.II/fe/fe_dg_vector.templates.h @@ -33,8 +33,8 @@ DEAL_II_NAMESPACE_OPEN template FE_DGVector::FE_DGVector(const unsigned int deg, MappingKind map) - : FE_PolyTensor( - deg, + : FE_PolyTensor( + PolynomialType(deg), FiniteElementData(get_dpo_vector(deg), dim, deg + 1, @@ -69,7 +69,7 @@ std::string FE_DGVector::get_name() const { std::ostringstream namebuf; - namebuf << "FE_DGVector_" << this->poly_space.name() << "<" << dim << ">(" + namebuf << "FE_DGVector_" << this->poly_space->name() << "<" << dim << ">(" << this->degree - 1 << ")"; return namebuf.str(); } diff --git a/include/deal.II/fe/fe_nedelec.h b/include/deal.II/fe/fe_nedelec.h index 93c1541f5c..65f3f0adbe 100644 --- a/include/deal.II/fe/fe_nedelec.h +++ b/include/deal.II/fe/fe_nedelec.h @@ -110,7 +110,7 @@ DEAL_II_NAMESPACE_OPEN * @date 2009, 2010, 2011 */ template -class FE_Nedelec : public FE_PolyTensor, dim> +class FE_Nedelec : public FE_PolyTensor { public: /** diff --git a/include/deal.II/fe/fe_poly_tensor.h b/include/deal.II/fe/fe_poly_tensor.h index 4010161cb5..1fa69c4156 100644 --- a/include/deal.II/fe/fe_poly_tensor.h +++ b/include/deal.II/fe/fe_poly_tensor.h @@ -22,6 +22,7 @@ #include #include #include +#include #include #include @@ -142,21 +143,24 @@ DEAL_II_NAMESPACE_OPEN * @author Guido Kanschat * @date 2005 */ -template +template class FE_PolyTensor : public FiniteElement { public: /** * Constructor. - * - * @arg @c degree: constructor argument for poly. May be different from @p - * fe_data.degree. */ - FE_PolyTensor(const unsigned int degree, + FE_PolyTensor(const TensorPolynomialsBase &polynomials, const FiniteElementData & fe_data, const std::vector & restriction_is_additive_flags, const std::vector &nonzero_components); + + /** + * Copy constructor. + */ + FE_PolyTensor(const FE_PolyTensor &fe); + // for documentation, see the FiniteElement base class virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override; @@ -318,12 +322,12 @@ protected: if (update_flags & (update_values | update_gradients)) for (unsigned int k = 0; k < n_q_points; ++k) { - poly_space.compute(quadrature.point(k), - values, - grads, - grad_grads, - third_derivatives, - fourth_derivatives); + poly_space->compute(quadrature.point(k), + values, + grads, + grad_grads, + third_derivatives, + fourth_derivatives); if (update_flags & update_values) { @@ -471,7 +475,7 @@ protected: * The polynomial space. Its type is given by the template parameter * PolynomialType. */ - PolynomialType poly_space; + std::unique_ptr> poly_space; /** * The inverse of the matrix aij of node values diff --git a/include/deal.II/fe/fe_raviart_thomas.h b/include/deal.II/fe/fe_raviart_thomas.h index dd4c284a68..ea5371ad59 100644 --- a/include/deal.II/fe/fe_raviart_thomas.h +++ b/include/deal.II/fe/fe_raviart_thomas.h @@ -102,8 +102,7 @@ DEAL_II_NAMESPACE_OPEN * @author Guido Kanschat, 2005, based on previous work by Wolfgang Bangerth. */ template -class FE_RaviartThomas - : public FE_PolyTensor, dim> +class FE_RaviartThomas : public FE_PolyTensor { public: /** @@ -244,8 +243,7 @@ private: * @author Guido Kanschat, 2005, Zhu Liang, 2008 */ template -class FE_RaviartThomasNodal - : public FE_PolyTensor, dim> +class FE_RaviartThomasNodal : public FE_PolyTensor { public: /** diff --git a/include/deal.II/fe/fe_rt_bubbles.h b/include/deal.II/fe/fe_rt_bubbles.h index ef537cba05..f0cacf2238 100644 --- a/include/deal.II/fe/fe_rt_bubbles.h +++ b/include/deal.II/fe/fe_rt_bubbles.h @@ -87,7 +87,7 @@ DEAL_II_NAMESPACE_OPEN * @author Eldar Khattatov, 2018 */ template -class FE_RT_Bubbles : public FE_PolyTensor, dim> +class FE_RT_Bubbles : public FE_PolyTensor { public: /** diff --git a/source/fe/fe_abf.cc b/source/fe/fe_abf.cc index df9e67520f..6a62cd7788 100644 --- a/source/fe/fe_abf.cc +++ b/source/fe/fe_abf.cc @@ -46,8 +46,8 @@ DEAL_II_NAMESPACE_OPEN template FE_ABF::FE_ABF(const unsigned int deg) - : FE_PolyTensor, dim>( - deg, + : FE_PolyTensor( + PolynomialsABF(deg), FiniteElementData(get_dpo_vector(deg), dim, deg + 2, diff --git a/source/fe/fe_bdm.cc b/source/fe/fe_bdm.cc index 63bbd947e0..c80345675b 100644 --- a/source/fe/fe_bdm.cc +++ b/source/fe/fe_bdm.cc @@ -38,8 +38,8 @@ DEAL_II_NAMESPACE_OPEN template FE_BDM::FE_BDM(const unsigned int deg) - : FE_PolyTensor, dim>( - deg, + : FE_PolyTensor( + PolynomialsBDM(deg), FiniteElementData(get_dpo_vector(deg), dim, deg + 1, diff --git a/source/fe/fe_bernardi_raugel.cc b/source/fe/fe_bernardi_raugel.cc index 6e30098e84..f4fbd7874b 100644 --- a/source/fe/fe_bernardi_raugel.cc +++ b/source/fe/fe_bernardi_raugel.cc @@ -38,8 +38,8 @@ DEAL_II_NAMESPACE_OPEN template FE_BernardiRaugel::FE_BernardiRaugel(const unsigned int p) - : FE_PolyTensor, dim>( - p, + : FE_PolyTensor( + PolynomialsBernardiRaugel(p), FiniteElementData(get_dpo_vector(), dim, 2, diff --git a/source/fe/fe_nedelec.cc b/source/fe/fe_nedelec.cc index 81ce3afe43..c92aeed119 100644 --- a/source/fe/fe_nedelec.cc +++ b/source/fe/fe_nedelec.cc @@ -69,8 +69,8 @@ namespace internal template FE_Nedelec::FE_Nedelec(const unsigned int order) - : FE_PolyTensor, dim>( - order, + : FE_PolyTensor( + PolynomialsNedelec(order), FiniteElementData(get_dpo_vector(order), dim, order + 1, diff --git a/source/fe/fe_poly_tensor.cc b/source/fe/fe_poly_tensor.cc index 2cb944c386..c7222ae526 100644 --- a/source/fe/fe_poly_tensor.cc +++ b/source/fe/fe_poly_tensor.cc @@ -166,9 +166,9 @@ namespace internal -template -FE_PolyTensor::FE_PolyTensor( - const unsigned int degree, +template +FE_PolyTensor::FE_PolyTensor( + const TensorPolynomialsBase &polynomials, const FiniteElementData & fe_data, const std::vector & restriction_is_additive_flags, const std::vector &nonzero_components) @@ -176,7 +176,7 @@ FE_PolyTensor::FE_PolyTensor( restriction_is_additive_flags, nonzero_components) , mapping_kind({MappingKind::mapping_none}) - , poly_space(PolynomialType(degree)) + , poly_space(polynomials.clone()) { cached_point(0) = -1; // Set up the table converting @@ -192,19 +192,28 @@ FE_PolyTensor::FE_PolyTensor( -template +template +FE_PolyTensor::FE_PolyTensor(const FE_PolyTensor &fe) + : FiniteElement(fe) + , mapping_kind(fe.mapping_kind) + , poly_space(fe.poly_space->clone()) + , inverse_node_matrix(fe.inverse_node_matrix) +{} + + + +template bool -FE_PolyTensor::single_mapping_kind() const +FE_PolyTensor::single_mapping_kind() const { return mapping_kind.size() == 1; } -template +template MappingKind -FE_PolyTensor::get_mapping_kind( - const unsigned int i) const +FE_PolyTensor::get_mapping_kind(const unsigned int i) const { if (single_mapping_kind()) return mapping_kind[0]; @@ -215,11 +224,10 @@ FE_PolyTensor::get_mapping_kind( -template +template double -FE_PolyTensor::shape_value( - const unsigned int, - const Point &) const +FE_PolyTensor::shape_value(const unsigned int, + const Point &) const { Assert(false, (typename FiniteElement::ExcFENotPrimitive())); @@ -228,9 +236,9 @@ FE_PolyTensor::shape_value( -template +template double -FE_PolyTensor::shape_value_component( +FE_PolyTensor::shape_value_component( const unsigned int i, const Point & p, const unsigned int component) const @@ -243,11 +251,11 @@ FE_PolyTensor::shape_value_component( if (cached_point != p || cached_values.size() == 0) { cached_point = p; - cached_values.resize(poly_space.n()); + cached_values.resize(poly_space->n()); std::vector> dummy1; std::vector> dummy2; - poly_space.compute( + poly_space->compute( p, cached_values, cached_grads, cached_grad_grads, dummy1, dummy2); } @@ -262,11 +270,10 @@ FE_PolyTensor::shape_value_component( -template +template Tensor<1, dim> -FE_PolyTensor::shape_grad( - const unsigned int, - const Point &) const +FE_PolyTensor::shape_grad(const unsigned int, + const Point &) const { Assert(false, (typename FiniteElement::ExcFENotPrimitive())); return Tensor<1, dim>(); @@ -274,9 +281,9 @@ FE_PolyTensor::shape_grad( -template +template Tensor<1, dim> -FE_PolyTensor::shape_grad_component( +FE_PolyTensor::shape_grad_component( const unsigned int i, const Point & p, const unsigned int component) const @@ -289,11 +296,11 @@ FE_PolyTensor::shape_grad_component( if (cached_point != p || cached_grads.size() == 0) { cached_point = p; - cached_grads.resize(poly_space.n()); + cached_grads.resize(poly_space->n()); std::vector> dummy1; std::vector> dummy2; - poly_space.compute( + poly_space->compute( p, cached_values, cached_grads, cached_grad_grads, dummy1, dummy2); } @@ -309,11 +316,10 @@ FE_PolyTensor::shape_grad_component( -template +template Tensor<2, dim> -FE_PolyTensor::shape_grad_grad( - const unsigned int, - const Point &) const +FE_PolyTensor::shape_grad_grad(const unsigned int, + const Point &) const { Assert(false, (typename FiniteElement::ExcFENotPrimitive())); return Tensor<2, dim>(); @@ -321,9 +327,9 @@ FE_PolyTensor::shape_grad_grad( -template +template Tensor<2, dim> -FE_PolyTensor::shape_grad_grad_component( +FE_PolyTensor::shape_grad_grad_component( const unsigned int i, const Point & p, const unsigned int component) const @@ -336,11 +342,11 @@ FE_PolyTensor::shape_grad_grad_component( if (cached_point != p || cached_grad_grads.size() == 0) { cached_point = p; - cached_grad_grads.resize(poly_space.n()); + cached_grad_grads.resize(poly_space->n()); std::vector> dummy1; std::vector> dummy2; - poly_space.compute( + poly_space->compute( p, cached_values, cached_grads, cached_grad_grads, dummy1, dummy2); } @@ -359,9 +365,9 @@ FE_PolyTensor::shape_grad_grad_component( // Fill data of FEValues //--------------------------------------------------------------------------- -template +template void -FE_PolyTensor::fill_fe_values( +FE_PolyTensor::fill_fe_values( const typename Triangulation::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature & quadrature, @@ -934,9 +940,9 @@ FE_PolyTensor::fill_fe_values( -template +template void -FE_PolyTensor::fill_fe_face_values( +FE_PolyTensor::fill_fe_face_values( const typename Triangulation::cell_iterator &cell, const unsigned int face_no, const Quadrature & quadrature, @@ -1559,9 +1565,9 @@ FE_PolyTensor::fill_fe_face_values( -template +template void -FE_PolyTensor::fill_fe_subface_values( +FE_PolyTensor::fill_fe_subface_values( const typename Triangulation::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, @@ -2177,9 +2183,9 @@ FE_PolyTensor::fill_fe_subface_values( -template +template UpdateFlags -FE_PolyTensor::requires_update_flags( +FE_PolyTensor::requires_update_flags( const UpdateFlags flags) const { UpdateFlags out = update_default; diff --git a/source/fe/fe_poly_tensor.inst.in b/source/fe/fe_poly_tensor.inst.in index 7501ece4f9..13e1c27634 100644 --- a/source/fe/fe_poly_tensor.inst.in +++ b/source/fe/fe_poly_tensor.inst.in @@ -17,16 +17,5 @@ for (deal_II_dimension : DIMENSIONS) { - template class FE_PolyTensor, - deal_II_dimension>; - template class FE_PolyTensor, - deal_II_dimension>; - template class FE_PolyTensor, - deal_II_dimension>; - template class FE_PolyTensor, - deal_II_dimension>; - template class FE_PolyTensor, - deal_II_dimension>; - template class FE_PolyTensor, - deal_II_dimension>; + template class FE_PolyTensor; } diff --git a/source/fe/fe_raviart_thomas.cc b/source/fe/fe_raviart_thomas.cc index 022e33a9f4..7007d543b5 100644 --- a/source/fe/fe_raviart_thomas.cc +++ b/source/fe/fe_raviart_thomas.cc @@ -45,8 +45,8 @@ DEAL_II_NAMESPACE_OPEN template FE_RaviartThomas::FE_RaviartThomas(const unsigned int deg) - : FE_PolyTensor, dim>( - deg, + : FE_PolyTensor( + PolynomialsRaviartThomas(deg), FiniteElementData(get_dpo_vector(deg), dim, deg + 1, diff --git a/source/fe/fe_raviart_thomas_nodal.cc b/source/fe/fe_raviart_thomas_nodal.cc index 0e6dc86d7c..802979994b 100644 --- a/source/fe/fe_raviart_thomas_nodal.cc +++ b/source/fe/fe_raviart_thomas_nodal.cc @@ -38,16 +38,15 @@ DEAL_II_NAMESPACE_OPEN template FE_RaviartThomasNodal::FE_RaviartThomasNodal(const unsigned int deg) - : FE_PolyTensor, dim>( - deg, - FiniteElementData(get_dpo_vector(deg), - dim, - deg + 1, - FiniteElementData::Hdiv), - get_ria_vector(deg), - std::vector(PolynomialsRaviartThomas::compute_n_pols( - deg), - std::vector(dim, true))) + : FE_PolyTensor(PolynomialsRaviartThomas(deg), + FiniteElementData(get_dpo_vector(deg), + dim, + deg + 1, + FiniteElementData::Hdiv), + get_ria_vector(deg), + std::vector( + PolynomialsRaviartThomas::compute_n_pols(deg), + std::vector(dim, true))) { Assert(dim >= 2, ExcImpossibleInDim(dim)); const unsigned int n_dofs = this->dofs_per_cell; diff --git a/source/fe/fe_rt_bubbles.cc b/source/fe/fe_rt_bubbles.cc index 74ad3afba4..7319cf29e4 100644 --- a/source/fe/fe_rt_bubbles.cc +++ b/source/fe/fe_rt_bubbles.cc @@ -37,16 +37,15 @@ DEAL_II_NAMESPACE_OPEN template FE_RT_Bubbles::FE_RT_Bubbles(const unsigned int deg) - : FE_PolyTensor, dim>( - deg, - FiniteElementData(get_dpo_vector(deg), - dim, - deg + 1, - FiniteElementData::Hdiv), - get_ria_vector(deg), - std::vector(PolynomialsRT_Bubbles::compute_n_pols( - deg), - std::vector(dim, true))) + : FE_PolyTensor(PolynomialsRT_Bubbles(deg), + FiniteElementData(get_dpo_vector(deg), + dim, + deg + 1, + FiniteElementData::Hdiv), + get_ria_vector(deg), + std::vector( + PolynomialsRT_Bubbles::compute_n_pols(deg), + std::vector(dim, true))) { Assert(dim >= 2, ExcImpossibleInDim(dim)); Assert( -- 2.39.5