From f621f9b33239b00b340c2c3b51cc2a5cd2e7a574 Mon Sep 17 00:00:00 2001 From: oliver Date: Tue, 25 Apr 2006 19:04:04 +0000 Subject: [PATCH] Added the header file for the polynomial space, which builds the basis of the ABF elements. git-svn-id: https://svn.dealii.org/trunk@12892 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/base/include/base/polynomials_abf.h | 179 ++++++++++++++++++++ 1 file changed, 179 insertions(+) create mode 100644 deal.II/base/include/base/polynomials_abf.h diff --git a/deal.II/base/include/base/polynomials_abf.h b/deal.II/base/include/base/polynomials_abf.h new file mode 100644 index 0000000000..91944c0082 --- /dev/null +++ b/deal.II/base/include/base/polynomials_abf.h @@ -0,0 +1,179 @@ +//--------------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2004, 2005 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_abf_h +#define __deal2__polynomials_abf_h + + +#include +#include +#include +#include +#include +#include +#include +#include + +#include + +/** + * @addtogroup Polynomials + * @{ + */ + +/** + * This class implements the Hdiv-conforming, + * vector-valued Arnold-Boffi-Falk polynomials as described in the + * article by Arnold-Boffi-Falk: + * Quadrilateral H(div) finite elements, SIAM J. Numer. Anal. + * Vol.42, No.6, pp.2429-2451 + * + * + * The ABF polynomials are constructed such that the + * divergence is in the tensor product polynomial space + * Qk. Therefore, the polynomial order of each + * component must be two orders higher in the corresponding direction, + * yielding the polynomial spaces (Qk+2,k, + * Qk,k+2) and (Qk+2,k,k, + * Qk,k+2,k, Qk,k,k+2) in 2D and 3D, resp. + * + * @author Oliver Kayser-Herold, 2006 based on code from Guido Kanschat + */ +template +class PolynomialsABF +{ + public: + /** + * Constructor. Creates all basis + * functions for Raviart-Thomas polynomials + * of given degree. + * + * @arg k: the degree of the + * Raviart-Thomas-space, which is the degree + * of the largest tensor product + * polynomial space + * Qk contained. + */ + PolynomialsABF (const unsigned int k); + + /** + * Destructor deleting the polynomials. + */ + ~PolynomialsABF (); + + /** + * Computes the value and the + * first and second derivatives + * of each Raviart-Thomas + * polynomial at @p unit_point. + * + * The size of the vectors must + * either be zero or equal + * n(). In the + * first case, the function will + * not compute these values. + * + * If you need values or + * derivatives of all tensor + * product polynomials then use + * this function, rather than + * using any of the + * compute_value, + * compute_grad or + * compute_grad_grad + * functions, see below, in a + * loop over all tensor product + * polynomials. + */ + void compute (const Point &unit_point, + std::vector > &values, + std::vector > &grads, + std::vector > &grad_grads) const; + + /** + * Returns the number of ABF polynomials. + */ + unsigned int n () const; + + /** + * Returns the degree of the ABF + * space, which is two less than + * the highest polynomial degree. + */ + unsigned int degree () const; + + /** + * Return the number of + * polynomials in the space + * RT(degree) without + * requiring to build an object + * of PolynomialsABF. This is + * required by the FiniteElement + * classes. + */ + static unsigned int compute_n_pols(unsigned int degree); + + private: + /** + * The degree of this object as + * given to the constructor. + */ + const unsigned int my_degree; + + /** + * An object representing the + * polynomial space for a single + * component. We can re-use it by + * rotating the coordinates of + * the evaluation point. + */ + AnisotropicPolynomials* polynomial_space; + + /** + * Number of Raviart-Thomas + * polynomials. + */ + unsigned int n_pols; + + /** + * Auxiliary memory. + */ + mutable std::vector p_values; + + /** + * Auxiliary memory. + */ + mutable std::vector > p_grads; + + /** + * Auxiliary memory. + */ + mutable std::vector > p_grad_grads; +}; + +/** @} */ + +template +inline unsigned int +PolynomialsABF::n() const +{ + return n_pols; +} + +template +inline unsigned int +PolynomialsABF::degree() const +{ + return my_degree; +} + +#endif -- 2.39.5