From: prill Date: Tue, 7 Nov 2006 08:26:59 +0000 (+0000) Subject: Made FE_DGQ with arbitrary nodes a derivation from FE_DGQ with equidistant nodes. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=7fe5074c54b4ef67d79f0d89b2ad49d58de40978;p=dealii-svn.git Made FE_DGQ with arbitrary nodes a derivation from FE_DGQ with equidistant nodes. git-svn-id: https://svn.dealii.org/trunk@14169 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_dgq.h b/deal.II/deal.II/include/fe/fe_dgq.h index 7fd475fb9b..84e3c9a746 100644 --- a/deal.II/deal.II/include/fe/fe_dgq.h +++ b/deal.II/deal.II/include/fe/fe_dgq.h @@ -88,18 +88,6 @@ class FE_DGQ : public FE_Poly,dim> */ FE_DGQ (const unsigned int p); - /** - * Constructor for tensor product - * polynomials based on - * Polynomials::Lagrange - * interpolation of the support - * points in the quadrature rule - * points. The degree of - * these polynomials is - * points.size()-1. - */ - FE_DGQ (const Quadrature<1>& points); - /** * Return a string that uniquely * identifies a finite @@ -324,6 +312,21 @@ class FE_DGQ : public FE_Poly,dim> protected: + /** + * Constructor for tensor product + * polynomials based on + * Polynomials::Lagrange + * interpolation of the support + * points in the quadrature rule + * points. The degree of + * these polynomials is + * points.size()-1. + * + * Note: The FE_DGQ::clone function + * does not work properly for FE with + * arbitrary nodes! + */ + FE_DGQ (const Quadrature<1>& points); /** * @p clone function instead of @@ -389,6 +392,55 @@ class FE_DGQ : public FE_Poly,dim> template friend class MappingQ; }; + + +/** + * Implementation of scalar, discontinuous tensor product elements + * based on Lagrange polynomials with arbitrary nodes. + * + * See base class documentation for details. + * + * @author F. Prill 2006 + */ +template +class FE_DGQArbitraryNodes : public FE_DGQ +{ + public: + /** + * Constructor for tensor product + * polynomials based on + * Polynomials::Lagrange + * interpolation of the support + * points in the quadrature rule + * points. The degree of + * these polynomials is + * points.size()-1. + */ + FE_DGQArbitraryNodes (const Quadrature<1>& points); + + /** + * Return a string that uniquely + * identifies a finite + * element. This class returns + * FE_DGQArbitraryNodes(degree) , + * with dim and degree + * replaced by appropriate + * values. + */ + virtual std::string get_name () const; + + protected: + /** + * @p clone function instead of + * a copy constructor. + * + * This function is needed by the + * constructors of @p FESystem. + */ + virtual FiniteElement *clone() const; +}; + + /*@}*/ DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/deal.II/source/fe/fe_dgq.cc b/deal.II/deal.II/source/fe/fe_dgq.cc index 01bd899ccc..47d6365c09 100644 --- a/deal.II/deal.II/source/fe/fe_dgq.cc +++ b/deal.II/deal.II/source/fe/fe_dgq.cc @@ -242,18 +242,7 @@ template FiniteElement * FE_DGQ::clone() const { - // TODO[Prill] : There must be a better way - // to extract 1D quadrature points from the - // tensor product FE. - - // Construct a dummy quadrature formula - // containing the FE's nodes: - std::vector > qpoints(this->degree+1); - for (unsigned int i=0; i<=this->degree; ++i) - qpoints[i] = Point<1>(this->unit_support_points[i][0]); - Quadrature<1> pquadrature(qpoints); - - return new FE_DGQ(pquadrature); + return new FE_DGQ(this->degree); } @@ -661,6 +650,52 @@ FE_DGQ::memory_consumption () const +template +FE_DGQArbitraryNodes::FE_DGQArbitraryNodes (const Quadrature<1>& points) + : FE_DGQ(points) +{} + + + +template +std::string +FE_DGQArbitraryNodes::get_name () const +{ + // note that the + // FETools::get_fe_from_name + // function does not work for + // FE_DGQArbitraryNodes since + // there is no initialization by + // a degree value. + std::ostringstream namebuf; + namebuf << "FE_DGQArbitraryNodes<" << dim << ">(" << this->degree << ")"; + + return namebuf.str(); +} + + + +template +FiniteElement * +FE_DGQArbitraryNodes::clone() const +{ + // TODO[Prill] : There must be a better way + // to extract 1D quadrature points from the + // tensor product FE. + + // Construct a dummy quadrature formula + // containing the FE's nodes: + std::vector > qpoints(this->degree+1); + for (unsigned int i=0; i<=this->degree; ++i) + qpoints[i] = Point<1>(this->unit_support_points[i][0]); + Quadrature<1> pquadrature(qpoints); + + return new FE_DGQArbitraryNodes(pquadrature); +} + + + template class FE_DGQ; +template class FE_DGQArbitraryNodes; DEAL_II_NAMESPACE_CLOSE