]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added the header file for the implementation of the ABF elements.
authoroliver <oliver@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 25 Apr 2006 19:06:57 +0000 (19:06 +0000)
committeroliver <oliver@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 25 Apr 2006 19:06:57 +0000 (19:06 +0000)
git-svn-id: https://svn.dealii.org/trunk@12895 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_abf.h [new file with mode: 0644]

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 (file)
index 0000000..1372522
--- /dev/null
@@ -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 <base/config.h>
+#include <base/table.h>
+#include <base/polynomials_abf.h>
+#include <base/polynomial.h>
+#include <base/tensor_product_polynomials.h>
+#include <base/geometry_info.h>
+#include <fe/fe.h>
+#include <fe/fe_poly_tensor.h>
+
+#include <vector>
+
+template <int dim> class MappingQ;
+
+
+/*!@addtogroup fe */
+/*@{*/
+
+/**
+ * Implementation of Arnold-Boffi-Falk (ABF) elements, conforming with the
+ * space H<sup>div</sup>. 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.
+ *
+ * <h3>Interpolation</h3>
+ *
+ * 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:
+ *
+ * <h4>Node values on edges/faces</h4>
+ *
+ * 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 <i>k</i> on an edge/face is the space
+ * <i>Q<sub>k</sub></i>, the moments are taken with respect to this
+ * space.
+ *
+ * <h4>Interior node values</h4>
+ *
+ * Higher order RT spaces have interior nodes. These are moments taken
+ * with respect to the gradient of functions in <i>Q<sub>k</sub></i>
+ * on the cell (this space is the matching space for RT<sub>k</sub> in
+ * a mixed formulation).
+ *
+ * <h4>Generalized support points</h4>
+ *
+ * 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
+ * QGauss<sub>k+1</sub> on each face as well as QGauss<sub>k</sub> in
+ * the interior of the cell (or none for RT<sub>0</sub>).
+ *
+ *
+ * @author Oliver Kayser-Herold, 2006, based on previous work
+ * by Guido Kanschat and Wolfgang Bangerth
+ */
+template <int dim>
+class FE_ABF
+  :
+  public FE_PolyTensor<PolynomialsABF<dim>, 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
+                                     * <tt>FE_ABF<dim>(degree)</tt>, 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, <tt>base_element(0)</tt> is
+                                     * @p this, and all other
+                                     * indices throw an error.
+                                     */
+    virtual const FiniteElement<dim> &
+    base_element (const unsigned int index) const;
+
+                                     /**
+                                      * Multiplicity of base element
+                                      * @p index. Since this is an
+                                      * atomic element,
+                                      * <tt>element_multiplicity(0)</tt>
+                                      * 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<double>&                local_dofs,
+                            const std::vector<double>& values) const;
+    virtual void interpolate(std::vector<double>&                local_dofs,
+                            const std::vector<Vector<double> >& values,
+                            unsigned int offset = 0) const;
+    virtual void interpolate(
+      std::vector<double>& local_dofs,
+      const VectorSlice<const std::vector<std::vector<double> > >& values) const;
+    virtual unsigned int memory_consumption () const;
+    virtual FiniteElement<dim> * 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<unsigned int>
+    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<dim>::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<std::vector<Tensor<1,dim> > > 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<std::vector<Tensor<2,dim> > > 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 <int dim1> friend class FE_ABF;
+};
+
+
+
+/*@}*/
+
+/* -------------- declaration of explicit specializations ------------- */
+
+#ifndef DOXYGEN
+
+template <>
+std::vector<unsigned int> FE_ABF<1>::get_dpo_vector (const unsigned int);
+
+#endif // DOXYGEN
+
+#endif

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.