From 00021d3bf214028a3811dcf1fdfd4329b1cfb903 Mon Sep 17 00:00:00 2001 From: Ralf Hartmann Date: Wed, 7 Mar 2001 10:46:23 +0000 Subject: [PATCH] TensorProductPolynomials keep their own copy of Polynomials. git-svn-id: https://svn.dealii.org/trunk@4138 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_dgq.h | 10 ++-------- deal.II/deal.II/include/fe/fe_q.h | 10 ++-------- deal.II/deal.II/include/fe/mapping_q.h | 9 --------- deal.II/deal.II/source/fe/fe_dgq.cc | 26 ++++++++++---------------- deal.II/deal.II/source/fe/fe_q.cc | 8 ++------ deal.II/deal.II/source/fe/mapping_q.cc | 15 +++++---------- 6 files changed, 21 insertions(+), 57 deletions(-) diff --git a/deal.II/deal.II/include/fe/fe_dgq.h b/deal.II/deal.II/include/fe/fe_dgq.h index 2f4f473f12..ff2e01d6fa 100644 --- a/deal.II/deal.II/include/fe/fe_dgq.h +++ b/deal.II/deal.II/include/fe/fe_dgq.h @@ -188,14 +188,8 @@ class FE_DGQ : public FiniteElement const unsigned int degree; /** - * Vector of one-dimensional - * polynomials used. - */ - std::vector > polynomials; - - /** - * Implementation of the tensor - * product of polynomials. + * Pointer to the tensor + * product polynomials. */ TensorProductPolynomials* poly; diff --git a/deal.II/deal.II/include/fe/fe_q.h b/deal.II/deal.II/include/fe/fe_q.h index 1d65c0b486..4314c96e47 100644 --- a/deal.II/deal.II/include/fe/fe_q.h +++ b/deal.II/deal.II/include/fe/fe_q.h @@ -263,16 +263,10 @@ class FE_Q : public FiniteElement * shape function numbering on first face. */ std::vector face_renumber; - - /** - * Vector of one-dimensional - * polynomials used. - */ - std::vector polynomials; /** - * Implementation of the tensor - * product of polynomials. + * Pointer to the tensor + * product polynomials. */ TensorProductPolynomials* poly; diff --git a/deal.II/deal.II/include/fe/mapping_q.h b/deal.II/deal.II/include/fe/mapping_q.h index 2874fae674..fd0f76c5b1 100644 --- a/deal.II/deal.II/include/fe/mapping_q.h +++ b/deal.II/deal.II/include/fe/mapping_q.h @@ -426,15 +426,6 @@ class MappingQ : public MappingQ1 */ const unsigned int n_outer; - - /** - * Vector of one-dimensional - * polynomials used as shape - * functions for the Qp mapping - * of cell at the boundary. - */ - std::vector polynomials; - /** * Pointer to the * @p{dim}-dimensional tensor diff --git a/deal.II/deal.II/source/fe/fe_dgq.cc b/deal.II/deal.II/source/fe/fe_dgq.cc index 2250f1a059..9e9c75693b 100644 --- a/deal.II/deal.II/source/fe/fe_dgq.cc +++ b/deal.II/deal.II/source/fe/fe_dgq.cc @@ -52,31 +52,25 @@ FE_DGQ::InternalData::~InternalData () template -FE_DGQ::FE_DGQ (unsigned int degree) - : +FE_DGQ::FE_DGQ (unsigned int degree): FiniteElement (FiniteElementData(get_dpo_vector(degree),1), dummy), degree(degree), - polynomials(degree+1), poly(0) { - std::vector > > v(degree+1); if (degree==0) { - std::vector coeff(1); - coeff[0] = 1.; - polynomials[0] = Polynomial (coeff); - v[0] = &(polynomials[0]); - } else { + std::vector > v( + 1, Polynomial (std::vector (1,1.))); + poly = new TensorProductPolynomials (v); + } + else + { + std::vector v; for (unsigned int i=0;i<=degree;++i) - { - LagrangeEquidistant p(degree, i); - polynomials[i] = p; - v[i] = &(polynomials[i]); - } + v.push_back(LagrangeEquidistant(degree,i)); + poly = new TensorProductPolynomials (v); } - poly = new TensorProductPolynomials (v); - Assert (degree <= 10, ExcNotImplemented()); std::vector right; diff --git a/deal.II/deal.II/source/fe/fe_q.cc b/deal.II/deal.II/source/fe/fe_q.cc index 7b7fffb3fd..a1f585675a 100644 --- a/deal.II/deal.II/source/fe/fe_q.cc +++ b/deal.II/deal.II/source/fe/fe_q.cc @@ -47,15 +47,11 @@ FE_Q::FE_Q (unsigned int degree) degree(degree), renumber(dofs_per_cell, 0), face_renumber(dofs_per_face, 0), - polynomials(degree+1), poly(0) { - std::vector > > v(degree+1); + std::vector v; for (unsigned int i=0;i<=degree;++i) - { - polynomials[i] = LagrangeEquidistant(degree,i); - v[i] = &(polynomials[i]); - }; + v.push_back(LagrangeEquidistant(degree,i)); poly = new TensorProductPolynomials (v); build_renumbering (*this, degree, renumber); diff --git a/deal.II/deal.II/source/fe/mapping_q.cc b/deal.II/deal.II/source/fe/mapping_q.cc index 1ce7839043..04b73b163f 100644 --- a/deal.II/deal.II/source/fe/mapping_q.cc +++ b/deal.II/deal.II/source/fe/mapping_q.cc @@ -45,7 +45,6 @@ MappingQ<1>::MappingQ (unsigned int): degree(1), n_inner(0), n_outer(0), - polynomials(0), tensor_pols(0), n_shape_functions(2), renumber(0), @@ -80,7 +79,6 @@ MappingQ::MappingQ (unsigned int p): n_inner(power(degree-1, dim)), n_outer((dim==2) ? 4+4*(degree-1) :8+12*(degree-1)+6*(degree-1)*(degree-1)), - polynomials(p+1), tensor_pols(0), n_shape_functions(0), renumber(0), @@ -91,14 +89,11 @@ MappingQ::MappingQ (unsigned int p): // polynomials used as shape // functions for the Qp mapping of // cells at the boundary. - std::vector > > pol_pointers(p+1); - for (unsigned int i=0; i<=p; ++i) - { - LagrangeEquidistant lagrange_pol(p, i); - polynomials[i] = lagrange_pol; - pol_pointers[i] = &(polynomials[i]); - } - tensor_pols = new TensorProductPolynomials (pol_pointers); + std::vector v; + for (unsigned int i=0;i<=degree;++i) + v.push_back(LagrangeEquidistant(degree,i)); + + tensor_pols = new TensorProductPolynomials (v); n_shape_functions=tensor_pols->n_tensor_product_polynomials(); Assert(n_inner+n_outer==n_shape_functions, ExcInternalError()); -- 2.39.5