From 6b4d2ee0c576f67cf45c59c3611266b88bd75896 Mon Sep 17 00:00:00 2001 From: oliver Date: Tue, 25 Apr 2006 19:06:57 +0000 Subject: [PATCH] Added the header file for the implementation of the ABF elements. git-svn-id: https://svn.dealii.org/trunk@12895 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_abf.h | 377 ++++++++++++++++++++++++++++ 1 file changed, 377 insertions(+) create mode 100644 deal.II/deal.II/include/fe/fe_abf.h diff --git a/deal.II/deal.II/include/fe/fe_abf.h b/deal.II/deal.II/include/fe/fe_abf.h new file mode 100644 index 0000000000..137252243d --- /dev/null +++ b/deal.II/deal.II/include/fe/fe_abf.h @@ -0,0 +1,377 @@ +//--------------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003, 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__fe_abf_h +#define __deal2__fe_abf_h + +#include +#include +#include +#include +#include +#include +#include +#include + +#include + +template class MappingQ; + + +/*!@addtogroup fe */ +/*@{*/ + +/** + * Implementation of Arnold-Boffi-Falk (ABF) elements, conforming with the + * space Hdiv. These elements generate vector fields with + * normal components continuous between mesh cells. + * + * These elements are based on an article from Arnold, Boffi and Falk: + * Quadrilateral H(div) finite elements, SIAM J. Numer. Anal. + * Vol.42, No.6, pp.2429-2451 + * + * In this article, the authors demonstrate that the usual RT elements + * and also BDM and other proposed finite dimensional subspaces of + * H(div) do not work properly on arbitrary FE grids. I.e. the + * convergence rates deteriorate on these meshes. As a solution the + * authors propose the ABF elements, which are implemented in this + * module. + * + * @todo Even if this element is implemented for two and three space + * dimensions, the definition of the node values relies on + * consistently oriented faces in 3D. Therefore, care should be taken + * on complicated meshes. + * + *

Interpolation

+ * + * The @ref GlossInterpolation "interpolation" operators associated + * with the RT element are constructed such that interpolation and + * computing the divergence are commuting operations. We require this + * from interpolating arbitrary functions as well as the #restriction + * matrices. It can be achieved by two interpolation schemes, the + * simplified one in FE_RaviartThomasNodal and the original one here: + * + *

Node values on edges/faces

+ * + * On edges or faces, the @ref GlossNodes "node values" are the moments of + * the normal component of the interpolated function with respect to + * the traces of the RT polynomials. Since the normal trace of the RT + * space of degree k on an edge/face is the space + * Qk, the moments are taken with respect to this + * space. + * + *

Interior node values

+ * + * Higher order RT spaces have interior nodes. These are moments taken + * with respect to the gradient of functions in Qk + * on the cell (this space is the matching space for RTk in + * a mixed formulation). + * + *

Generalized support points

+ * + * The node values above rely on integrals, which will be computed by + * quadrature rules themselves. The generalized support points are a + * set of points such that this quadrature can be performed with + * sufficient accuracy. The points needed are those of + * QGaussk+1 on each face as well as QGaussk in + * the interior of the cell (or none for RT0). + * + * + * @author Oliver Kayser-Herold, 2006, based on previous work + * by Guido Kanschat and Wolfgang Bangerth + */ +template +class FE_ABF + : + public FE_PolyTensor, dim> +{ + public: + /** + * Constructor for the ABF + * element of degree @p p. + */ + FE_ABF (const unsigned int p); + + /** + * Return a string that uniquely + * identifies a finite + * element. This class returns + * FE_ABF(degree), with + * @p dim and @p degree + * replaced by appropriate + * values. + */ + virtual std::string get_name () const; + + + /** + * Number of base elements in a + * mixed discretization. Here, + * this is of course equal to + * one. + */ + virtual unsigned int n_base_elements () const; + + /** + * Access to base element + * objects. Since this element is + * atomic, base_element(0) is + * @p this, and all other + * indices throw an error. + */ + virtual const FiniteElement & + base_element (const unsigned int index) const; + + /** + * Multiplicity of base element + * @p index. Since this is an + * atomic element, + * element_multiplicity(0) + * returns one, and all other + * indices will throw an error. + */ + virtual unsigned int element_multiplicity (const unsigned int index) const; + + /** + * Check whether a shape function + * may be non-zero on a face. + * + * Right now, this is only + * implemented for RT0 in + * 1D. Otherwise, returns always + * @p true. + */ + virtual bool has_support_on_face (const unsigned int shape_index, + const unsigned int face_index) const; + + virtual void interpolate(std::vector& local_dofs, + const std::vector& values) const; + virtual void interpolate(std::vector& local_dofs, + const std::vector >& values, + unsigned int offset = 0) const; + virtual void interpolate( + std::vector& local_dofs, + const VectorSlice > >& values) const; + virtual unsigned int memory_consumption () const; + virtual FiniteElement * clone() const; + + private: + /** + * The order of the + * ABF element. The + * lowest order elements are + * usually referred to as RT0, + * even though their shape + * functions are piecewise + * quadratics. + */ + const unsigned int rt_order; + + /** + * Only for internal use. Its + * full name is + * @p get_dofs_per_object_vector + * function and it creates the + * @p dofs_per_object vector that is + * needed within the constructor to + * be passed to the constructor of + * @p FiniteElementData. + */ + static std::vector + get_dpo_vector (const unsigned int degree); + + /** + * Initialize the @p + * generalized_support_points + * field of the FiniteElement + * class and fill the tables with + * interpolation weights + * (#boundary_weights and + * #interior_weights). Called + * from the constructor. + */ + void initialize_support_points (const unsigned int rt_degree); + + /** + * Given a set of flags indicating + * what quantities are requested + * from a @p FEValues object, + * return which of these can be + * precomputed once and for + * all. Often, the values of + * shape function at quadrature + * points can be precomputed, for + * example, in which case the + * return value of this function + * would be the logical and of + * the input @p flags and + * @p update_values. + * + * For the present kind of finite + * element, this is exactly the + * case. + */ + virtual UpdateFlags update_once (const UpdateFlags flags) const; + + /** + * This is the opposite to the + * above function: given a set of + * flags indicating what we want + * to know, return which of these + * need to be computed each time + * we visit a new cell. + * + * If for the computation of one + * quantity something else is + * also required (for example, we + * often need the covariant + * transformation when gradients + * need to be computed), include + * this in the result as well. + */ + virtual UpdateFlags update_each (const UpdateFlags flags) const; + + /** + * Fields of cell-independent data. + * + * For information about the + * general purpose of this class, + * see the documentation of the + * base class. + */ + class InternalData : public FiniteElement::InternalDataBase + { + public: + /** + * Array with shape function + * values in quadrature + * points. There is one row + * for each shape function, + * containing values for each + * quadrature point. Since + * the shape functions are + * vector-valued (with as + * many components as there + * are space dimensions), the + * value is a tensor. + * + * In this array, we store + * the values of the shape + * function in the quadrature + * points on the unit + * cell. The transformation + * to the real space cell is + * then simply done by + * multiplication with the + * Jacobian of the mapping. + */ + std::vector > > shape_values; + + /** + * Array with shape function + * gradients in quadrature + * points. There is one + * row for each shape + * function, containing + * values for each quadrature + * point. + * + * We store the gradients in + * the quadrature points on + * the unit cell. We then + * only have to apply the + * transformation (which is a + * matrix-vector + * multiplication) when + * visiting an actual cell. + */ + std::vector > > shape_gradients; + }; + + /** + * These are the factors + * multiplied to a function in + * the + * #generalized_face_support_points + * when computing the + * integration. They are + * organized such that there is + * one row for each generalized + * face support point and one + * column for each degree of + * freedom on the face. + */ + Table<2, double> boundary_weights; + /** + * Precomputed factors for + * interpolation of interior + * degrees of freedom. The + * rationale for this Table is + * the same as for + * #boundary_weights. Only, this + * table has a third coordinate + * for the space direction of the + * component evaluated. + */ + Table<3, double> interior_weights; + + + + /** + * These are the factors + * multiplied to a function in + * the + * #generalized_face_support_points + * when computing the + * integration. They are + * organized such that there is + * one row for each generalized + * face support point and one + * column for each degree of + * freedom on the face. + */ + Table<2, double> boundary_weights_abf; + /** + * Precomputed factors for + * interpolation of interior + * degrees of freedom. The + * rationale for this Table is + * the same as for + * #boundary_weights. Only, this + * table has a third coordinate + * for the space direction of the + * component evaluated. + */ + Table<3, double> interior_weights_abf; + + + /** + * Allow access from other + * dimensions. + */ + template friend class FE_ABF; +}; + + + +/*@}*/ + +/* -------------- declaration of explicit specializations ------------- */ + +#ifndef DOXYGEN + +template <> +std::vector FE_ABF<1>::get_dpo_vector (const unsigned int); + +#endif // DOXYGEN + +#endif -- 2.39.5