From 88c82c681421ebe0d2c8783e211cf6c56ccf26c8 Mon Sep 17 00:00:00 2001 From: Ralf Hartmann Date: Tue, 9 Mar 2004 14:45:44 +0000 Subject: [PATCH] Implementation of a complete polynomial space of order k based on monomials. Makes use of the PolynomialSpace class and imposes an ordering of the monomials which ensures identical ordering of Subspaces. git-svn-id: https://svn.dealii.org/trunk@8699 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/base/include/base/polynomials_p.h | 113 ++++++++++++++++++++++ deal.II/base/source/polynomials_p.cc | 98 +++++++++++++++++++ 2 files changed, 211 insertions(+) create mode 100644 deal.II/base/include/base/polynomials_p.h create mode 100644 deal.II/base/source/polynomials_p.cc diff --git a/deal.II/base/include/base/polynomials_p.h b/deal.II/base/include/base/polynomials_p.h new file mode 100644 index 0000000000..ee110547fc --- /dev/null +++ b/deal.II/base/include/base/polynomials_p.h @@ -0,0 +1,113 @@ +//------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2004 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//------------------------------------------------------------------- +#ifndef __deal2__polynomials_P_h +#define __deal2__polynomials_P_h + + +#include +#include +#include +#include +#include +#include +#include + +#include + + +/** + * @brief The complete polynomial space of order k based on + * the monomials. + * + * This class implements the polynomial space of order k + * based on the monomials ${1,x,x^2,...}$. I.e. in d + * dimensions it constructs all polynomials of the form $\prod_{i=1}^d + * x_i^{n_i}$, where $\sum_i n_i\leq k$. The base polynomials are + * given a specific ordering, e.g. in 2 dimensions: + * ${1,x,y,xy,x^2,y^2,x^2y,xy^2,x^3,y^3,...}$. The ordering of the + * monomials in $P_k1$ matches the ordering of the monomials in $P_k2$ + * for $k2>k1$. + * + * @author Ralf Hartmann, 2004 + */ +template +class PolynomialsP: public PolynomialSpace +{ + public: + /** + * Constructor. Creates all basis + * functions of $P_k$. + * @arg k: the degree of the + * polynomial space + */ + PolynomialsP (const unsigned int k); + + /** + * Returns the degree k + * of the polynomial space + * P_k. + * + * Note, that this number is + * PolynomialSpace::degree()-1, + * compare definition in + * @ref{PolynomialSpace}. + */ + unsigned int degree() const; + + /** + * For the nth + * polynomial $p_n(x,y,z)=x^i y^j + * z^k$ this function gives the + * degrees i,j,k in the x,y,z + * directions. + */ + void directional_degrees(unsigned int n, + unsigned int (°rees)[dim]) const; + + private: + + /** + * Fills the index_map. + */ + void create_polynomial_ordering(vector &index_map) const; + + /** + * Degree k of the + * polynomial space $P_k$, + * i.e. the number k + * which was given to the + * constructor. + */ + const unsigned int k; +}; + + + +template +inline unsigned int +PolynomialsP::degree() const +{ + return k; +} + + +template +inline void +PolynomialsP::directional_degrees(unsigned int n, + unsigned int (°rees)[dim]) const +{ + compute_index(n,degrees); +} + + +#endif diff --git a/deal.II/base/source/polynomials_p.cc b/deal.II/base/source/polynomials_p.cc new file mode 100644 index 0000000000..e2ce67889f --- /dev/null +++ b/deal.II/base/source/polynomials_p.cc @@ -0,0 +1,98 @@ +//------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2004 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//------------------------------------------------------------------- + +#include + + +template +PolynomialsP::PolynomialsP (const unsigned int k) + : + PolynomialSpace(Polynomials::Monomial::generate_complete_basis(k)), + k(k) +{ + vector index_map(n()); + create_polynomial_ordering(index_map); + set_polynomial_ordering(index_map); +} + + +template <> +void PolynomialsP<1>::create_polynomial_ordering( + vector &index_map) const +{ + Assert(index_map.size()==n(), ExcDimensionMismatch(index_map.size(), n())); + + // identity + for (unsigned int i=0; i +void PolynomialsP<2>::create_polynomial_ordering( + vector &index_map) const +{ + Assert(index_map.size()==n(), ExcDimensionMismatch(index_map.size(), n())); + Assert(k<=5, ExcNotImplemented()); + + // Given the number i of the + // polynomial in + // $1,x,y,xy,x2,y2,...$, + // index_map[i] gives the number of + // the polynomial in + // PolynomialSpace. + for (unsigned int i=0; i +void PolynomialsP<3>::create_polynomial_ordering( + vector &index_map) const +{ + Assert(index_map.size()==n(), ExcDimensionMismatch(index_map.size(), n())); + Assert(k<=3, ExcNotImplemented()); + + // Given the number i of the + // polynomial in + // $1,x,y,xy,x2,y2,...$, + // index_map[i] gives the number of + // the polynomial in + // PolynomialSpace. + for (unsigned int i=0; i; +template class PolynomialsP<2>; +template class PolynomialsP<3>; -- 2.39.5