]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
FiniteElementBase removed
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 25 Aug 2005 10:45:31 +0000 (10:45 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 25 Aug 2005 10:45:31 +0000 (10:45 +0000)
git-svn-id: https://svn.dealii.org/trunk@11326 0785d39b-7218-0410-832d-ea1e28bc413d

24 files changed:
deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_base.h
deal.II/deal.II/include/fe/fe_dgp_monomial.h
deal.II/deal.II/include/fe/fe_dgp_nonparametric.h
deal.II/deal.II/include/fe/fe_dgq.h
deal.II/deal.II/include/fe/fe_nedelec.h
deal.II/deal.II/include/fe/fe_poly.h
deal.II/deal.II/include/fe/fe_poly_tensor.h
deal.II/deal.II/include/fe/fe_q.h
deal.II/deal.II/include/fe/fe_q_hierarchical.h
deal.II/deal.II/include/fe/fe_raviart_thomas.h
deal.II/deal.II/include/fe/fe_system.h
deal.II/deal.II/include/fe/fe_tools.h
deal.II/deal.II/source/dofs/dof_tools.cc
deal.II/deal.II/source/fe/fe.cc
deal.II/deal.II/source/fe/fe_dgp_monomial.cc
deal.II/deal.II/source/fe/fe_dgq.cc
deal.II/deal.II/source/fe/fe_poly_tensor.cc
deal.II/deal.II/source/fe/fe_q.cc
deal.II/deal.II/source/fe/fe_raviart_thomas.cc
deal.II/deal.II/source/fe/fe_system.cc
deal.II/deal.II/source/fe/fe_tools.cc
deal.II/doc/news/changes.html

index 520b77e004cdfe9c2f9253826b0faf1e344b18d3..0e5245033c9fd2bfe706b96296dbb70530f4cab8 100644 (file)
@@ -130,7 +130,7 @@ template <int dim> class Mapping;
  * this is the point where it has the function value 1), is
  * located on the boundary. We do not check this directly, the
  * criterion is rather defined through the information the finite
- * element class gives: the FiniteElementBase class defines the
+ * element class gives: the FiniteElement class defines the
  * numbers of basis functions per vertex, per line, and so on and the
  * basis functions are numbered after this information; a basis
  * function is to be considered to be on the face of a cell (and thus
index 6c5e66119bdf06c23ea9847b29fe2fe84d15dcb8..817f1f9178b860178ee58296f6b81ae854738577 100644 (file)
@@ -30,9 +30,7 @@ template <int dim> class FECollection;
 
 /**
  * Common interface of all finite elements. Here, the functions to
- * fill the data fields of FEValues are declared. While
- * FiniteElementBase provides implementation of common
- * functionality, this class only serves as an abstract base class.
+ * fill the data fields of FEValues are declared.
  *
  * The interface of this class is very restrictive. The reason is that
  * finite element values should be accessed only by use of
@@ -54,8 +52,59 @@ template <int dim> class FECollection;
  * @author Wolfgang Bangerth, Guido Kanschat, Ralf Hartmann, 1998, 2000, 2001
  */
 template <int dim>
-class FiniteElement : public FiniteElementBase<dim>
+class FiniteElement : public Subscriptor,
+                     public FiniteElementData<dim>
 {
+  public:
+                                  /**
+                                   * Base class for internal data.
+                                   * Adds data for second derivatives to
+                                   * Mapping::InternalDataBase()
+                                   *
+                                   * For information about the
+                                   * general purpose of this class,
+                                   * see the documentation of the
+                                   * base class.
+                                   *
+                                   * @author Guido Kanschat, 2001
+                                   */
+  class InternalDataBase : public Mapping<dim>::InternalDataBase
+    {
+      public:      
+                                        /**
+                                         * Destructor. Needed to
+                                         * avoid memory leaks with
+                                         * difference quotients.
+                                         */
+       virtual ~InternalDataBase ();
+
+                                        /**
+                                         * Initialize some pointers
+                                         * used in the computation of
+                                         * second derivatives by
+                                         * finite differencing of
+                                         * gradients.
+                                         */
+       void initialize_2nd (const FiniteElement<dim> *element,
+                            const Mapping<dim>       &mapping,
+                            const Quadrature<dim>    &quadrature);
+       
+                                        /**
+                                         * Storage for @p FEValues
+                                         * objects needed to
+                                         * approximate second
+                                         * derivatives.
+                                         *
+                                         * The ordering is <tt>p+hx</tt>,
+                                         * <tt>p+hy</tt>, <tt>p+hz</tt>,
+                                         * @p p-hx, @p p-hy,
+                                         * @p p-hz, where unused
+                                         * entries in lower dimensions
+                                         * are missing.
+                                         */
+       std::vector<FEValues<dim>*> differences;
+    };
+
   public:
                                      /**
                                       * Copy constructor. This one is declared
@@ -87,6 +136,694 @@ class FiniteElement : public FiniteElementBase<dim>
                                      */
     virtual ~FiniteElement ();
     
+                                    /**
+                                     * Return a string that uniquely
+                                     * identifies a finite
+                                     * element. The general
+                                     * convention is that this is the
+                                     * class name, followed by the
+                                     * space dimension in angle
+                                     * brackets, and the polynomial
+                                     * degree and whatever else is
+                                     * necessary in parentheses. For
+                                     * example, <tt>FE_Q<2>(3)</tt> is the
+                                     * value returned for a cubic
+                                     * element in 2d.
+                                     *
+                                     * Systems of elements have their
+                                     * own naming convention, see the
+                                     * FESystem class.
+                                     */
+    virtual std::string get_name () const = 0;
+
+                                    /**
+                                     * @name Shape function access
+                                     * @{
+                                     */
+    
+                                    /**
+                                     * Return the value of the
+                                     * @p ith shape function at the
+                                     * point @p p. @p p is a point
+                                     * on the reference element. If
+                                     * the finite element is
+                                     * vector-valued, then return the
+                                     * value of the only non-zero
+                                     * component of the vector value
+                                     * of this shape function. If the
+                                     * shape function has more than
+                                     * one non-zero component (which
+                                     * we refer to with the term
+                                     * non-primitive), then derived
+                                     * classes implementing this
+                                     * function should throw an
+                                     * exception of type
+                                     * @p ExcShapeFunctionNotPrimitive. In
+                                     * that case, use the
+                                     * shape_value_component()
+                                     * function.
+                                     *
+                                     * An
+                                     * @p ExcUnitShapeValuesDoNotExist
+                                     * is thrown if the shape values
+                                     * of the @p FiniteElement under
+                                     * consideration depends on the
+                                     * shape of the cell in real
+                                     * space.
+                                     */
+    virtual double shape_value (const unsigned int  i,
+                               const Point<dim>   &p) const;
+
+                                    /**
+                                     * Just like for @p shape_value,
+                                     * but this function will be
+                                     * called when the shape function
+                                     * has more than one non-zero
+                                     * vector component. In that
+                                     * case, this function should
+                                     * return the value of the
+                                     * @p component-th vector
+                                     * component of the @p ith shape
+                                     * function at point @p p.
+                                     */
+    virtual double shape_value_component (const unsigned int i,
+                                         const Point<dim>   &p,
+                                         const unsigned int component) const;
+    
+                                    /**
+                                     * Return the gradient of the
+                                     * @p ith shape function at the
+                                     * point @p p. @p p is a point
+                                     * on the reference element, and
+                                     * likewise the gradient is the
+                                     * gradient on the unit cell with
+                                     * respect to unit cell
+                                     * coordinates. If
+                                     * the finite element is
+                                     * vector-valued, then return the
+                                     * value of the only non-zero
+                                     * component of the vector value
+                                     * of this shape function. If the
+                                     * shape function has more than
+                                     * one non-zero component (which
+                                     * we refer to with the term
+                                     * non-primitive), then derived
+                                     * classes implementing this
+                                     * function should throw an
+                                     * exception of type
+                                     * @p ExcShapeFunctionNotPrimitive. In
+                                     * that case, use the
+                                     * shape_grad_component()
+                                     * function.
+                                     *
+                                     * An
+                                     * @p ExcUnitShapeValuesDoNotExist
+                                     * is thrown if the shape values
+                                     * of the @p FiniteElement under
+                                     * consideration depends on the
+                                     * shape of the cell in real
+                                     * space.
+                                     */
+    virtual Tensor<1,dim> shape_grad (const unsigned int  i,
+                                     const Point<dim>   &p) const;
+
+                                    /**
+                                     * Just like for @p shape_grad,
+                                     * but this function will be
+                                     * called when the shape function
+                                     * has more than one non-zero
+                                     * vector component. In that
+                                     * case, this function should
+                                     * return the gradient of the
+                                     * @p component-th vector
+                                     * component of the @p ith shape
+                                     * function at point @p p.
+                                     */
+    virtual Tensor<1,dim> shape_grad_component (const unsigned int i,
+                                               const Point<dim>   &p,
+                                               const unsigned int component) const;
+
+                                    /**
+                                     * Return the tensor of second
+                                     * derivatives of the @p ith
+                                     * shape function at point @p p
+                                     * on the unit cell. The
+                                     * derivatives are derivatives on
+                                     * the unit cell with respect to
+                                     * unit cell coordinates. If
+                                     * the finite element is
+                                     * vector-valued, then return the
+                                     * value of the only non-zero
+                                     * component of the vector value
+                                     * of this shape function. If the
+                                     * shape function has more than
+                                     * one non-zero component (which
+                                     * we refer to with the term
+                                     * non-primitive), then derived
+                                     * classes implementing this
+                                     * function should throw an
+                                     * exception of type
+                                     * @p ExcShapeFunctionNotPrimitive. In
+                                     * that case, use the
+                                     * shape_grad_grad_component()
+                                     * function.
+                                     *
+                                     * An
+                                     * @p ExcUnitShapeValuesDoNotExist
+                                     * is thrown if the shape values
+                                     * of the @p FiniteElement under
+                                     * consideration depends on the
+                                     * shape of the cell in real
+                                     * space.
+                                     */
+    virtual Tensor<2,dim> shape_grad_grad (const unsigned int  i,
+                                          const Point<dim>   &p) const;
+
+                                    /**
+                                     * Just like for @p shape_grad_grad,
+                                     * but this function will be
+                                     * called when the shape function
+                                     * has more than one non-zero
+                                     * vector component. In that
+                                     * case, this function should
+                                     * return the gradient of the
+                                     * @p component-th vector
+                                     * component of the @p ith shape
+                                     * function at point @p p.
+                                     */
+    virtual Tensor<2,dim> shape_grad_grad_component (const unsigned int i,
+                                                    const Point<dim>   &p,
+                                                    const unsigned int component) const;
+                                    /**
+                                     * Check for non-zero values on a face.
+                                     *
+                                     * This function returns
+                                     * @p true, if the shape
+                                     * function @p shape_index has
+                                     * non-zero values on the face
+                                     * @p face_index.
+                                     */
+    virtual bool has_support_on_face (const unsigned int shape_index,
+                                     const unsigned int face_index) const = 0;
+    
+                                    //@}
+                                    /**
+                                     * @name Transfer and constraint matrices
+                                     * @{
+                                     */
+    
+                                    /**
+                                     * Projection from a fine grid
+                                     * space onto a coarse grid
+                                     * space. If this projection
+                                     * operator is associated with a
+                                     * matrix @p P, then the
+                                     * restriction of this matrix
+                                     * @p P_i to a single child cell
+                                     * is returned here.
+                                     *
+                                     * The matrix @p P is the
+                                     * concatenation or the sum of
+                                     * the cell matrices @p P_i,
+                                     * depending on the
+                                     * @p restriction_is_additive_flags
+                                     * given to the constructor. This
+                                     * distinguishes interpolation
+                                     * (concatenation) and projection
+                                     * with respect to scalar
+                                     * products (summation).
+                                     *
+                                     * Row and column indices are
+                                     * related to coarse grid and
+                                     * fine grid spaces,
+                                     * respectively, consistent with
+                                     * the definition of the
+                                     * associated operator.
+                                     *
+                                     * If projection matrices are not
+                                     * implemented in the derived
+                                     * finite element class, this
+                                     * function aborts with
+                                     * @p ExcProjectionVoid.
+                                     */
+    const FullMatrix<double> &
+    get_restriction_matrix (const unsigned int child) const;
+
+                                    /**
+                                     * Embedding matrix between grids.
+                                     * 
+                                     * The identity operator from a
+                                     * coarse grid space into a fine
+                                     * grid space is associated with
+                                     * a matrix @p P. The
+                                     * restriction of this matrix @p P_i to
+                                     * a single child cell is
+                                     * returned here.
+                                     *
+                                     * The matrix @p P is the
+                                     * concatenation, not the sum of
+                                     * the cell matrices
+                                     * @p P_i. That is, if the same
+                                     * non-zero entry <tt>j,k</tt> exists
+                                     * in in two different child
+                                     * matrices @p P_i, the value
+                                     * should be the same in both
+                                     * matrices and it is copied into
+                                     * the matrix @p P only once.
+                                     *
+                                     * Row and column indices are
+                                     * related to fine grid and
+                                     * coarse grid spaces,
+                                     * respectively, consistent with
+                                     * the definition of the
+                                     * associated operator.
+                                     *
+                                     * These matrices are used by
+                                     * routines assembling the
+                                     * prolongation matrix for
+                                     * multi-level methods.  Upon
+                                     * assembling the transfer matrix
+                                     * between cells using this
+                                     * matrix array, zero elements in
+                                     * the prolongation matrix are
+                                     * discarded and will not fill up
+                                     * the transfer matrix.
+                                     *
+                                     * If projection matrices are not
+                                     * implemented in the derived
+                                     * finite element class, this
+                                     * function aborts with
+                                     * @p ExcEmbeddingVoid. You can
+                                     * check whether this is the case
+                                     * by calling the
+                                     * prolongation_is_implemented().
+                                     */
+    const FullMatrix<double> &
+    get_prolongation_matrix (const unsigned int child) const;
+
+                                     /**
+                                      * Return whether this element implements
+                                      * its prolongation matrices. The return
+                                      * value also indicates whether a call to
+                                      * the @p get_prolongation_matrix
+                                      * function will generate an error or
+                                      * not.
+                                      *
+                                      * This function is mostly here in order
+                                      * to allow us to write more efficient
+                                      * test programs which we run on all
+                                      * kinds of weird elements, and for which
+                                      * we simply need to exclude certain
+                                      * tests in case something is not
+                                      * implemented. It will in general
+                                      * probably not be a great help in
+                                      * applications, since there is not much
+                                      * one can do if one needs these features
+                                      * and they are not implemented. This
+                                      * function could be used to check
+                                      * whether a call to
+                                      * <tt>get_prolongation_matrix()</tt> will
+                                      * succeed; however, one then still needs
+                                      * to cope with the lack of information
+                                      * this just expresses.
+                                      */
+    bool prolongation_is_implemented () const;
+
+                                     /**
+                                      * Return whether this element implements
+                                      * its restriction matrices. The return
+                                      * value also indicates whether a call to
+                                      * the @p get_restriction_matrix
+                                      * function will generate an error or
+                                      * not.
+                                      *
+                                      * This function is mostly here in order
+                                      * to allow us to write more efficient
+                                      * test programs which we run on all
+                                      * kinds of weird elements, and for which
+                                      * we simply need to exclude certain
+                                      * tests in case something is not
+                                      * implemented. It will in general
+                                      * probably not be a great help in
+                                      * applications, since there is not much
+                                      * one can do if one needs these features
+                                      * and they are not implemented. This
+                                      * function could be used to check
+                                      * whether a call to
+                                      * <tt>get_restriction_matrix()</tt> will
+                                      * succeed; however, one then still needs
+                                      * to cope with the lack of information
+                                      * this just expresses.
+                                      */
+    bool restriction_is_implemented () const;
+
+                                    /**
+                                     * Access the
+                                     * @p restriction_is_additive_flag
+                                     * field. See there for more
+                                     * information on its contents.
+                                     *
+                                     * The index must be between zero
+                                     * and the number of shape
+                                     * functions of this element.
+                                     */
+    bool restriction_is_additive (const unsigned int index) const;
+
+                                    /**
+                                     * Return a readonly reference to
+                                     * the matrix which describes the
+                                     * constraints at the interface
+                                     * between a refined and an
+                                     * unrefined cell.
+                                     * 
+                                     * The matrix is obviously empty
+                                     * in only one space dimension,
+                                     * since there are no constraints
+                                     * then.
+                                     *
+                                     * Note that some finite elements
+                                     * do not (yet) implement hanging
+                                     * node constraints. If this is
+                                     * the case, then this function
+                                     * will generate an exception,
+                                     * since no useful return value
+                                     * can be generated. If you
+                                     * should have a way to live with
+                                     * this, then you might want to
+                                     * use the
+                                     * @p constraints_are_implemented
+                                     * function to check up front
+                                     * whethehr this function will
+                                     * succeed or generate the
+                                     * exception.
+                                     */
+    const FullMatrix<double> & constraints () const;
+
+                                     /**
+                                      * Return whether this element
+                                      * implements its hanging node
+                                      * constraints. The return value
+                                      * also indicates whether a call
+                                      * to the @p constraint function
+                                      * will generate an error or not.
+                                      *
+                                      * This function is mostly here
+                                      * in order to allow us to write
+                                      * more efficient test programs
+                                      * which we run on all kinds of
+                                      * weird elements, and for which
+                                      * we simply need to exclude
+                                      * certain tests in case hanging
+                                      * node constraints are not
+                                      * implemented. It will in
+                                      * general probably not be a
+                                      * great help in applications,
+                                      * since there is not much one
+                                      * can do if one needs hanging
+                                      * node constraints and they are
+                                      * not implemented. This function
+                                      * could be used to check whether
+                                      * a call to <tt>constraints()</tt>
+                                      * will succeed; however, one
+                                      * then still needs to cope with
+                                      * the lack of information this
+                                      * just expresses.
+                                      */
+    bool constraints_are_implemented () const;
+
+                                    /**
+                                     * Return the matrix
+                                     * interpolating from the given
+                                     * finite element to the present
+                                     * one. The size of the matrix is
+                                     * then @p dofs_per_cell times
+                                     * <tt>source.dofs_per_cell</tt>.
+                                     *
+                                     * Derived elements will have to
+                                     * implement this function. They
+                                     * may only provide interpolation
+                                     * matrices for certain source
+                                     * finite elements, for example
+                                     * those from the same family. If
+                                     * they don't implement
+                                     * interpolation from a given
+                                     * element, then they must throw
+                                     * an exception of type
+                                     * FiniteElement<dim>::ExcInterpolationNotImplemented.
+                                     */
+    virtual void
+    get_interpolation_matrix (const FiniteElement<dim> &source,
+                             FullMatrix<double>           &matrix) const;
+                                    //@}
+    
+                                    /**
+                                     * Comparison operator. We also
+                                     * check for equality of the
+                                     * constraint matrix, which is
+                                     * quite an expensive operation.
+                                     * Do therefore use this function
+                                     * with care, if possible only
+                                     * for debugging purposes.
+                                     *
+                                     * Since this function is not
+                                     * that important, we avoid an
+                                     * implementational question
+                                     * about comparing arrays and do
+                                     * not compare the matrix arrays
+                                     * @p restriction and
+                                     * @p prolongation.
+                                     */
+    bool operator == (const FiniteElement<dim> &) const;
+
+                                    /**
+                                     * @name Index computations
+                                     * @{
+                                     */
+                                    /**
+                                     * Compute vector component and
+                                     * index of this shape function
+                                     * within the shape functions
+                                     * corresponding to this
+                                     * component from the index of a
+                                     * shape function within this
+                                     * finite element.
+                                     *
+                                     * If the element is scalar, then
+                                     * the component is always zero,
+                                     * and the index within this
+                                     * component is equal to the
+                                     * overall index.
+                                     *
+                                     * If the shape function
+                                     * referenced has more than one
+                                     * non-zero component, then it
+                                     * cannot be associated with one
+                                     * vector component, and an
+                                     * exception of type
+                                     * @p ExcShapeFunctionNotPrimitive
+                                     * will be raised.
+                                     *
+                                     * Note that if the element is
+                                     * composed of other (base)
+                                     * elements, and a base element
+                                     * has more than one component
+                                     * but all its shape functions
+                                     * are primitive (i.e. are
+                                     * non-zero in only one
+                                     * component), then this mapping
+                                     * contains valid
+                                     * information. However, the
+                                     * index of a shape function of
+                                     * this element within one
+                                     * component (i.e. the second
+                                     * number of the respective entry
+                                     * of this array) does not
+                                     * indicate the index of the
+                                     * respective shape function
+                                     * within the base element (since
+                                     * that has more than one
+                                     * vector-component). For this
+                                     * information, refer to the
+                                     * @p system_to_base_table field
+                                     * and the
+                                     * @p system_to_base_index
+                                     * function.
+                                     */
+    std::pair<unsigned int, unsigned int>
+    system_to_component_index (const unsigned int index) const;
+
+                                    /**
+                                     * Compute the shape function for
+                                     * the given vector component and
+                                     * index.
+                                     *
+                                     * If the element is scalar, then
+                                     * the component must be zero,
+                                     * and the index within this
+                                     * component is equal to the
+                                     * overall index.
+                                     *
+                                     * This is the opposite operation
+                                     * from the @p system_to_component_index
+                                     * function.
+                                     */
+   unsigned int component_to_system_index(const unsigned int component,
+                                          const unsigned int index) const;
+  
+                                    /**
+                                     * Same as above, but do it for
+                                     * shape functions and their
+                                     * indices on a face.
+                                     */
+    std::pair<unsigned int, unsigned int>
+    face_system_to_component_index (const unsigned int index) const;
+
+                                     /**
+                                      * Return for shape function
+                                      * @p index the base element it
+                                      * belongs to, the number of the
+                                      * copy of this base element
+                                      * (which is between zero and the
+                                      * multiplicity of this element),
+                                      * and the index of this shape
+                                      * function within this base
+                                      * element.
+                                      *
+                                      * If the element is not composed of
+                                     * others, then base and instance
+                                     * are always zero, and the index
+                                     * is equal to the number of the
+                                     * shape function. If the element
+                                     * is composed of single
+                                     * instances of other elements
+                                     * (i.e. all with multiplicity
+                                     * one) all of which are scalar,
+                                     * then base values and dof
+                                     * indices within this element
+                                     * are equal to the
+                                     * @p system_to_component_table. It
+                                     * differs only in case the
+                                     * element is composed of other
+                                     * elements and at least one of
+                                     * them is vector-valued itself.
+                                     *
+                                     * This function returns valid
+                                     * values also in the case of
+                                     * vector-valued
+                                     * (i.e. non-primitive) shape
+                                     * functions, in contrast to the
+                                     * @p system_to_component_index
+                                     * function.
+                                      */
+    std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
+    system_to_base_index (const unsigned int index) const;
+
+                                     /**
+                                      * Same as
+                                      * @p system_to_base_index, but
+                                      * for degrees of freedom located
+                                      * on a face.
+                                      */
+    std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
+    face_system_to_base_index (const unsigned int index) const;
+    
+                                    /**
+                                     * Return in which of the vector
+                                     * components of this finite
+                                     * element the @p ithe shape
+                                     * function is non-zero. The
+                                     * length of the returned array
+                                     * is equal to the number of
+                                     * vector components of this
+                                     * element.
+                                     *
+                                     * For most finite element
+                                     * spaces, the result of this
+                                     * function will be a vector with
+                                     * exactly one element being
+                                     * @p true, since for most
+                                     * spaces the individual vector
+                                     * components are independent. In
+                                     * that case, the component with
+                                     * the single zero is also the
+                                     * first element of what
+                                     * <tt>system_to_component_index(i)</tt>
+                                     * returns.
+                                     *
+                                     * Only for those
+                                     * spaces that couple the
+                                     * components, for example to
+                                     * make a shape function
+                                     * divergence free, will there be
+                                     * more than one @p true entry.
+                                     */
+    const std::vector<bool> &
+    get_nonzero_components (const unsigned int i) const;
+
+                                    /**
+                                     * Return in how many vector
+                                     * components the @p ith shape
+                                     * function is non-zero. This
+                                     * value equals the number of
+                                     * entries equal to @p true in
+                                     * the result of the
+                                     * @p get_nonzero_components
+                                     * function.
+                                     *
+                                     * For most finite element
+                                     * spaces, the result will be
+                                     * equal to one. It is not equal
+                                     * to one only for those ansatz
+                                     * spaces for which vector-valued
+                                     * shape functions couple the
+                                     * individual components, for
+                                     * example in order to make them
+                                     * divergence-free.
+                                     */
+    unsigned int
+    n_nonzero_components (const unsigned int i) const;
+
+                                    /**
+                                     * Return whether the @p ith
+                                     * shape function is primitive in
+                                     * the sense that the shape
+                                     * function is non-zero in only
+                                     * one vector
+                                     * component. Non-primitive shape
+                                     * functions would then, for
+                                     * example, be those of
+                                     * divergence free ansatz spaces,
+                                     * in which the individual vector
+                                     * components are coupled.
+                                     *
+                                     * The result of the function is
+                                     * @p true if and only if the
+                                     * result of
+                                     * <tt>n_nonzero_components(i)</tt> is
+                                     * equal to one.
+                                     */
+    bool
+    is_primitive (const unsigned int i) const;
+
+                                    /**
+                                     * Return whether the entire
+                                     * finite element is primitive,
+                                     * in the sense that all its
+                                     * shape functions are
+                                     * primitive. If the finite
+                                     * element is scalar, then this
+                                     * is always the case.
+                                     *
+                                     * Since this is an extremely
+                                     * common operation, the result
+                                     * is cached in the
+                                     * @p cached_primitivity
+                                     * variable which is computed in
+                                     * the constructor.
+                                     */
+    bool
+    is_primitive () const;
+    
                                     /**
                                      * Number of base elements in a
                                      * mixed discretization.
@@ -110,7 +847,7 @@ class FiniteElement : public FiniteElementBase<dim>
                                      * multiplicity.
                                      */
     virtual unsigned int n_base_elements () const = 0;
-    
+
                                     /**
                                      * Access to base element
                                      * objects. If the element is
@@ -136,17 +873,293 @@ class FiniteElement : public FiniteElementBase<dim>
     unsigned int
     element_multiplicity (const unsigned int index) const = 0;
     
+                                    /**
+                                     * Given a vector component,
+                                     * return an index which base
+                                     * element implements this
+                                     * component, and which vector
+                                     * component in this base element
+                                     * this is. This information is
+                                     * only of interest for
+                                     * vector-valued finite elements
+                                     * which are composed of several
+                                     * sub-elements. In that case,
+                                     * one may want to obtain
+                                     * information about the element
+                                     * implementing a certain vector
+                                     * component, which can be done
+                                     * using this function and the
+                                     * FESystem::@p base_element
+                                     * function.
+                                     *
+                                     * If this is a scalar finite
+                                     * element, then the return value
+                                     * is always equal to a pair of
+                                     * zeros.
+                                     */
+    std::pair<unsigned int,unsigned int>
+    component_to_base (const unsigned int component) const;
+                                    //@}
+    
                                     /**
-                                     * Check for non-zero values on a face.
+                                     * @name Support points and interpolation
+                                     * @{
+                                     */
+    
+                                    /**
+                                     * Return the support points of
+                                     * the trial functions on the
+                                     * unit cell, if the derived
+                                     * finite element defines some.
+                                     * Finite elements that allow
+                                     * some kind of interpolation
+                                     * operation usually have support
+                                     * points. On the other hand,
+                                     * elements that define their
+                                     * degrees of freedom by, for
+                                     * example, moments on faces, or
+                                     * as derivatives, don't have
+                                     * support points. In that case,
+                                     * the returned field is empty.
                                      *
-                                     * This function returns
-                                     * @p true, if the shape
-                                     * function @p shape_index has
-                                     * non-zero values on the face
-                                     * @p face_index.
+                                     * If the finite element defines
+                                     * support points, then their
+                                     * number equals the number of
+                                     * degrees of freedom of the
+                                     * element.  The order of points
+                                     * in the array matches that
+                                     * returned by the
+                                     * <tt>cell->get_dof_indices</tt>
+                                     * function.
+                                     *
+                                     * See the class documentation
+                                     * for details on support points.
                                      */
-    virtual bool has_support_on_face (const unsigned int shape_index,
-                                     const unsigned int face_index) const = 0;
+    const std::vector<Point<dim> > &
+    get_unit_support_points () const;    
+
+                                    /**
+                                     * Return whether a finite
+                                     * element has defined support
+                                     * points. If the result is true,
+                                     * then a call to the
+                                     * @p get_unit_support_points
+                                     * yields a non-empty array.
+                                     *
+                                     * The result may be false if an
+                                     * element is not defined by
+                                     * interpolating shape functions,
+                                     * for example by P-elements on
+                                     * quadrilaterals. It will
+                                     * usually only be true if the
+                                     * element constructs its shape
+                                     * functions by the requirement
+                                     * that they be one at a certain
+                                     * point and zero at all the
+                                     * points associated with the
+                                     * other shape functions.
+                                     *
+                                     * In composed elements (i.e. for
+                                     * the FESystem class, the
+                                     * result will be true if all all
+                                     * the base elements have defined
+                                     * support points.
+                                     */
+    bool has_support_points () const;
+
+                                     /**
+                                      * Return the position of the
+                                      * support point of the
+                                      * @p indexth shape function. If
+                                      * it does not exist, raise an
+                                      * exception.
+                                      *
+                                      * The default implementation
+                                      * simply returns the respective
+                                      * element from the array you get
+                                      * from
+                                      * get_unit_support_points(),
+                                      * but derived elements may
+                                      * overload this function. In
+                                      * particular, note that the
+                                      * FESystem class overloads
+                                      * it so that it can return the
+                                      * support points of individual
+                                      * base elements, of not all the
+                                      * base elements define support
+                                      * points. In this way, you can
+                                      * still ask for certain support
+                                      * points, even if
+                                      * @p get_unit_support_points
+                                      * only returns an empty array.
+                                      */
+    virtual
+    Point<dim>
+    unit_support_point (const unsigned int index) const;
+    
+                                    /**
+                                     * Return the support points of
+                                     * the trial functions on the
+                                     * unit face, if the derived
+                                     * finite element defines some.
+                                     * Finite elements that allow
+                                     * some kind of interpolation
+                                     * operation usually have support
+                                     * points. On the other hand,
+                                     * elements that define their
+                                     * degrees of freedom by, for
+                                     * example, moments on faces, or
+                                     * as derivatives, don't have
+                                     * support points. In that case,
+                                     * the returned field is empty
+                                     *
+                                     * Note that elements that have
+                                     * support points need not
+                                     * necessarily have some on the
+                                     * faces, even if the
+                                     * interpolation points are
+                                     * located physically on a
+                                     * face. For example, the
+                                     * discontinuous elements have
+                                     * interpolation points on the
+                                     * vertices, and for higher
+                                     * degree elements also on the
+                                     * faces, but they are not
+                                     * defined to be on faces since
+                                     * in that case degrees of
+                                     * freedom from both sides of a
+                                     * face (or from all adjacent
+                                     * elements to a vertex) would be
+                                     * identified with each other,
+                                     * which is not what we would
+                                     * like to have). Logically,
+                                     * these degrees of freedom are
+                                     * therefore defined to belong to
+                                     * the cell, rather than the face
+                                     * or vertex. In that case, the
+                                     * returned element would
+                                     * therefore have length zero.
+                                     *
+                                     * If the finite element defines
+                                     * support points, then their
+                                     * number equals the number of
+                                     * degrees of freedom on the face
+                                     * (@p dofs_per_face). The order
+                                     * of points in the array matches
+                                     * that returned by the
+                                     * <tt>cell->get_dof_indices</tt>
+                                     * function.
+                                     *
+                                     * See the class documentation
+                                     * for details on support points.
+                                     */
+    const std::vector<Point<dim-1> > &
+    get_unit_face_support_points () const;    
+
+                                    /**
+                                     * Return whether a finite
+                                     * element has defined support
+                                     * points on faces. If the result
+                                     * is true, then a call to the
+                                     * @p get_unit_support_points
+                                     * yields a non-empty array.
+                                     *
+                                     * For more information, see the
+                                     * documentation for the
+                                     * has_support_points()
+                                     * function.
+                                     */
+    bool has_face_support_points () const;
+
+                                     /**
+                                      * The function corresponding to
+                                      * the unit_support_point()
+                                      * function, but for faces. See
+                                      * there for more information.
+                                      */
+    virtual
+    Point<dim-1>
+    unit_face_support_point (const unsigned int index) const;
+    
+                                    /**
+                                     * Return a support point vector
+                                     * for generalized interpolation.
+                                     */
+    const std::vector<Point<dim> > &
+    get_generalized_support_points () const;    
+
+                                    /**
+                                     *
+                                     */
+    bool has_generalized_support_points () const;
+
+                                    /**
+                                     *
+                                     */
+    const std::vector<Point<dim-1> > &
+    get_generalized_face_support_points () const;
+
+                                    /**
+                                     * Return whether a finite
+                                     * element has defined support
+                                     * points on faces. If the result
+                                     * is true, then a call to the
+                                     * @p get_unit_support_points
+                                     * yields a non-empty array.
+                                     *
+                                     * For more information, see the
+                                     * documentation for the
+                                     * has_support_points()
+                                     * function.
+                                     */
+    bool has_generalized_face_support_points () const;
+
+                                    /**
+                                     * Interpolate a set of scalar
+                                     * values, computed in the
+                                     * generalized support points.
+                                     *
+                                     * @note This function is
+                                     * implemented in
+                                     * FiniteElement for the case
+                                     * that the element has support
+                                     * points. In this case, the
+                                     * resulting coefficients are
+                                     * just the values in the suport
+                                     * points. All other elements
+                                     * must reimplement it.
+                                     */
+    virtual void interpolate(std::vector<double>&       local_dofs,
+                            const std::vector<double>& values) const;
+      
+                                    /**
+                                     * Interpolate a set of vector
+                                     * values, computed in the
+                                     * generalized support points.
+                                     *
+                                     * Since a finite element often
+                                     * only interpolates part of a
+                                     * vector, <tt>offset</tt> is
+                                     * used to determine the first
+                                     * component of the vector to be
+                                     * interpolated. Maybe consider
+                                     * changing your data structures
+                                     * to use the next function.
+                                     */
+    virtual void interpolate(std::vector<double>&                local_dofs,
+                            const std::vector<Vector<double> >& values,
+                            unsigned int offset = 0) const;
+      
+                                    /**
+                                     * Interpolate a set of vector
+                                     * values, computed in the
+                                     * generalized support points.
+                                     */
+    virtual void interpolate(
+      std::vector<double>& local_dofs,
+      const VectorSlice<const std::vector<std::vector<double> > >& values) const;
+      
+                                    //@}
     
                                     /**
                                      * Determine an estimate for the
@@ -160,6 +1173,474 @@ class FiniteElement : public FiniteElementBase<dim>
                                      * rather than the class itself.
                                      */
     virtual unsigned int memory_consumption () const;
+                                    /**
+                                     * Exception
+                                     *
+                                     * @ingroup Exceptions
+                                     */
+    DeclException1 (ExcShapeFunctionNotPrimitive,
+                   int,
+                   << "The shape function with index " << arg1
+                   << " is not primitive, i.e. it is vector-valued and "
+                   << "has more than one non-zero vector component. This "
+                   << "function cannot be called for these shape functions. "
+                   << "Maybe you want to use the same function with the "
+                   << "_component suffix?");
+                                    /**
+                                     * Exception
+                                     *
+                                     * @ingroup Exceptions
+                                     */
+    DeclException0 (ExcFENotPrimitive);
+                                    /**
+                                     * Exception
+                                     *
+                                     * @ingroup Exceptions
+                                     */
+    DeclException0 (ExcUnitShapeValuesDoNotExist);
+
+                                    /**
+                                     * Attempt to access support
+                                     * points of a finite element
+                                     * which is not Lagrangian.
+                                     *
+                                     * @ingroup Exceptions
+                                     */
+    DeclException0 (ExcFEHasNoSupportPoints);
+
+                                    /**
+                                     * Attempt to access embedding
+                                     * matrices of a finite element
+                                     * which did not implement these
+                                     * matrices.
+                                     *
+                                     * @ingroup Exceptions
+                                     */
+    DeclException0 (ExcEmbeddingVoid);
+    
+                                    /**
+                                     * Attempt to access restriction
+                                     * matrices of a finite element
+                                     * which did not implement these
+                                     * matrices.
+                                     *
+                                     * Exception
+                                     * @ingroup Exceptions
+                                     */
+    DeclException0 (ExcProjectionVoid);
+    
+                                    /**
+                                     * Attempt to access constraint
+                                     * matrices of a finite element
+                                     * which did not implement these
+                                     * matrices.
+                                     *
+                                     * Exception
+                                     * @ingroup Exceptions
+                                     */
+    DeclException0 (ExcConstraintsVoid);
+    
+                                    /**
+                                     * Exception
+                                     * @ingroup Exceptions
+                                     */
+    DeclException2 (ExcWrongInterfaceMatrixSize,
+                   int, int,
+                   << "The interface matrix has a size of " << arg1
+                   << "x" << arg2
+                   << ", which is not reasonable in the present dimension.");
+                                    /**
+                                     * Exception
+                                     * @ingroup Exceptions
+                                     */
+    DeclException2 (ExcComponentIndexInvalid,
+                   int, int,
+                   << "The component-index pair (" << arg1 << ", " << arg2
+                   << ") is invalid, i.e. non-existent");
+                                     /**
+                                      * Exception
+                                     * @ingroup Exceptions
+                                      */
+    DeclException0 (ExcInterpolationNotImplemented);
+    
+  protected:  
+                                    /**
+                                     * Array of projection matrices. See
+                                     * get_restriction_matrix() above.
+                                     *
+                                     * Matrices in this array are
+                                     * automatically initialized to
+                                     * correct size. If the derived
+                                     * finite element class does not
+                                     * implement these matrices, they
+                                     * should be resized to zero
+                                     * size.
+                                     */
+    FullMatrix<double> restriction[GeometryInfo<dim>::children_per_cell];
+
+                                    /**
+                                     * Array of embedding matrices. See
+                                     * <tt>get_prolongation_matrix()</tt> above.
+                                     *
+                                     * Matrices in this array are
+                                     * automatically initialized to
+                                     * correct size. If the derived
+                                     * finite element class does not
+                                     * implement these matrices, they
+                                     * should be resized to zero
+                                     * size.
+                                     */
+    FullMatrix<double> prolongation[GeometryInfo<dim>::children_per_cell];
+
+                                    /**
+                                     * Specify the constraints which
+                                     * the dofs on the two sides of a
+                                     * cell interface underly if the
+                                     * line connects two cells of
+                                     * which one is refined once.
+                                     *
+                                     * For further details see the
+                                     * general description of the
+                                     * derived class.
+                                     *
+                                     * This field is obviously
+                                     * useless in one space dimension
+                                     * and has there a zero size.
+                                     */
+    FullMatrix<double> interface_constraints;
+
+                                     /**
+                                      * Return the size of interface
+                                      * constraint matrices. Since
+                                      * this is needed in every
+                                      * derived finite element class
+                                      * when initializing their size,
+                                      * it is placed into this
+                                      * function, to avoid having to
+                                      * recompute the
+                                      * dimension-dependent size of
+                                      * these matrices each time.
+                                      *
+                                      * Note that some elements do not
+                                      * implement the interface
+                                      * constraints for certain
+                                      * polynomial degrees. In this
+                                      * case, this function still
+                                      * returns the size these
+                                      * matrices should have when
+                                      * implemented, but the actual
+                                      * matrices are empty.
+                                      */
+    TableIndices<2>
+    interface_constraints_size () const;
+    
+                                    /**
+                                     * Store what
+                                     * @p system_to_component_index
+                                     * will return.
+                                     */
+    std::vector< std::pair<unsigned int, unsigned int> > system_to_component_table;
+
+                                     /**
+                                     * Map between linear dofs and
+                                     * component dofs on face. This
+                                     * is filled with default values
+                                     * in the constructor, but
+                                     * derived classes will have to
+                                     * overwrite the information if
+                                     * necessary.
+                                     *
+                                     * By component, we mean the
+                                     * vector component, not the base
+                                     * element. The information thus
+                                     * makes only sense if a shape
+                                     * function is non-zero in only
+                                     * one component.
+                                     */
+    std::vector< std::pair<unsigned int, unsigned int> > face_system_to_component_table;
+
+                                    /**
+                                     * For each shape function, store
+                                     * to which base element and
+                                     * which instance of this base
+                                     * element (in case its
+                                     * multiplicity is greater than
+                                     * one) it belongs, and its index
+                                     * within this base element. If
+                                     * the element is not composed of
+                                     * others, then base and instance
+                                     * are always zero, and the index
+                                     * is equal to the number of the
+                                     * shape function. If the element
+                                     * is composed of single
+                                     * instances of other elements
+                                     * (i.e. all with multiplicity
+                                     * one) all of which are scalar,
+                                     * then base values and dof
+                                     * indices within this element
+                                     * are equal to the
+                                     * @p system_to_component_table. It
+                                     * differs only in case the
+                                     * element is composed of other
+                                     * elements and at least one of
+                                     * them is vector-valued itself.
+                                     *
+                                     * This array has valid values
+                                     * also in the case of
+                                     * vector-valued
+                                     * (i.e. non-primitive) shape
+                                     * functions, in contrast to the
+                                     * @p system_to_component_table.
+                                     */
+    std::vector<std::pair<std::pair<unsigned int,unsigned int>,unsigned int> >
+    system_to_base_table;
+
+                                    /**
+                                     * Likewise for the indices on
+                                     * faces.
+                                     */
+    std::vector<std::pair<std::pair<unsigned int,unsigned int>,unsigned int> >
+    face_system_to_base_table;
+    
+                                    /**
+                                     * The base element establishing
+                                     * a component.
+                                     *
+                                     * This table converts a
+                                     * component number to a pair
+                                     * consisting of the
+                                     * @p base_element number, and
+                                     * the component within this base
+                                     * element. While component
+                                     * information contains
+                                     * multiplicity of base elements,
+                                     * the result allows access to
+                                     * shape functions of the base
+                                     * element.
+                                     *
+                                     * This variable is set to the
+                                     * correct size by the
+                                     * constructor of this class, but
+                                     * needs to be initialized by
+                                     * derived classes, unless its
+                                     * size is one and the only entry
+                                     * is a zero, which is the case
+                                     * for scalar elements. In that
+                                     * case, the initialization by
+                                     * the base class is sufficient.
+                                     */
+    std::vector<std::pair<unsigned int, unsigned int> > component_to_base_table;
+    
+                                    /**
+                                     * Projection matrices are
+                                     * concatenated or summed up.
+                                     *
+                                     * This flags decides on how the
+                                     * projection matrices of the
+                                     * children of the same father
+                                     * are put together to one
+                                     * operator. The possible modes
+                                     * are concatenation and
+                                     * summation.
+                                     *
+                                     * If the projection is defined
+                                     * by an interpolation operator,
+                                     * the child matrices are
+                                     * concatenated, i.e. values
+                                     * belonging to the same node
+                                     * functional are identified and
+                                     * enter the interpolated value
+                                     * only once. In this case, the
+                                     * flag must be @p false.
+                                     *
+                                     * For projections with respect
+                                     * to scalar products, the child
+                                     * matrices must be summed up to
+                                     * build the complete matrix. The
+                                     * flag should be @p true.
+                                     *
+                                     * For examples of use of these
+                                     * flags, see the places in the
+                                     * library where it is queried.
+                                     * 
+                                     * There is one flag per shape
+                                     * function, indicating whether
+                                     * it belongs to the class of
+                                     * shape functions that are
+                                     * additive in the restriction or
+                                     * not.
+                                     *
+                                     * Note that in previous versions
+                                     * of the library, there was one
+                                     * flag per vector component of
+                                     * the element. This is based on
+                                     * the fact that all the shape
+                                     * functions that belong to the
+                                     * same vector component must
+                                     * necessarily behave in the same
+                                     * way, to make things
+                                     * reasonable. However, the
+                                     * problem is that it is
+                                     * sometimes impossible to query
+                                     * this flag in the vector-valued
+                                     * case: this used to be done
+                                     * with the
+                                     * @p system_to_component_index
+                                     * function that returns which
+                                     * vector component a shape
+                                     * function is associated
+                                     * with. The point is that since
+                                     * we now support shape functions
+                                     * that are associated with more
+                                     * than one vector component (for
+                                     * example the shape functions of
+                                     * Raviart-Thomas, or Nedelec
+                                     * elements), that function can
+                                     * no more be used, so it can be
+                                     * difficult to find out which
+                                     * for vector component we would
+                                     * like to query the
+                                     * restriction-is-additive flags.
+                                     */
+    const std::vector<bool> restriction_is_additive_flags;
+
+                                    /**
+                                     * List of support points on the
+                                     * unit cell, in case the finite
+                                     * element has any. The
+                                     * constructor leaves this field
+                                     * empty, derived classes may
+                                     * write in some contents.
+                                     *
+                                     * Finite elements that allow
+                                     * some kind of interpolation
+                                     * operation usually have support
+                                     * points. On the other hand,
+                                     * elements that define their
+                                     * degrees of freedom by, for
+                                     * example, moments on faces, or
+                                     * as derivatives, don't have
+                                     * support points. In that case,
+                                     * this field remains empty.
+                                     */
+    std::vector<Point<dim> > unit_support_points;
+
+                                    /**
+                                     * Same for the faces. See the
+                                     * description of the
+                                     * @p get_unit_face_support_points
+                                     * function for a discussion of
+                                     * what contributes a face
+                                     * support point.
+                                     */
+    std::vector<Point<dim-1> > unit_face_support_points;
+    
+                                    /**
+                                     * Support points used for
+                                     * interpolation functions of
+                                     * non-Lagrangian elements.
+                                     */
+    std::vector<Point<dim> > generalized_support_points;
+    
+                                    /**
+                                     * Face support points used for
+                                     * interpolation functions of
+                                     * non-Lagrangian elements.
+                                     */    
+    std::vector<Point<dim-1> > generalized_face_support_points;
+    
+                                    /**
+                                     * For each shape function, give
+                                     * a vector of bools (with size
+                                     * equal to the number of vector
+                                     * components which this finite
+                                     * element has) indicating in
+                                     * which component each of these
+                                     * shape functions is non-zero.
+                                     *
+                                     * For primitive elements, there
+                                     * is only one non-zero
+                                     * component.
+                                     */
+    const std::vector<std::vector<bool> > nonzero_components;
+
+                                    /**
+                                     * This array holds how many
+                                     * values in the respective entry
+                                     * of the @p nonzero_components
+                                     * element are non-zero. The
+                                     * array is thus a short-cut to
+                                     * allow faster access to this
+                                     * information than if we had to
+                                     * count the non-zero entries
+                                     * upon each request for this
+                                     * information. The field is
+                                     * initialized in the constructor
+                                     * of this class.
+                                     */
+    const std::vector<unsigned int> n_nonzero_components_table;
+
+                                    /**
+                                     * Store whether all shape
+                                     * functions are primitive. Since
+                                     * finding this out is a very
+                                     * common operation, we cache the
+                                     * result, i.e. compute the value
+                                     * in the constructor for simpler
+                                     * access.
+                                     */
+    const bool cached_primitivity;
+
+                                     /**
+                                     * Compute second derivatives by
+                                     * finite differences of
+                                     * gradients.
+                                     */
+    void compute_2nd (const Mapping<dim>                      &mapping,
+                     const typename Triangulation<dim>::cell_iterator    &cell,
+                     const unsigned int                       offset,
+                     typename Mapping<dim>::InternalDataBase &mapping_internal,
+                     InternalDataBase                        &fe_internal,
+                     FEValuesData<dim>                       &data) const;
+
+                                    /**
+                                     * Given the pattern of nonzero
+                                     * components for each shape
+                                     * function, compute for each
+                                     * entry how many components are
+                                     * non-zero for each shape
+                                     * function. This function is
+                                     * used in the constructor of
+                                     * this class.
+                                     */
+    static
+    std::vector<unsigned int>
+    compute_n_nonzero_components (const std::vector<std::vector<bool> > &nonzero_components);
+    
+                                    /**
+                                     * Allow the FESystem class to access the
+                                     * restriction and prolongation matrices
+                                     * directly. Hence, FESystem has the
+                                     * possibility to see if these matrices
+                                     * are initialized without accessing
+                                     * these matrices through the
+                                     * @p get_restriction_matrix and
+                                     * @p get_prolongation_matrix
+                                     * functions. This is important as these
+                                     * functions include assertions that
+                                     * throw if the matrices are not already
+                                     * initialized.
+                                     */
+    template <int dim_> friend class FESystem;
+
+                                     /**
+                                      * Make the inner class a
+                                      * friend. This is not strictly
+                                      * necessary, but the Intel
+                                      * compiler seems to want this.
+                                      */
+    friend class InternalDataBase;
+    
 
                                     /**
                                      * Exception
@@ -244,6 +1725,20 @@ class FiniteElement : public FiniteElementBase<dim>
                                      */
     virtual FiniteElement<dim> *clone() const = 0;
     
+  private:
+                                    /**
+                                     * Second derivatives of shapes
+                                     * functions are not computed
+                                     * analytically, but by finite
+                                     * differences of the
+                                     * gradients. This static
+                                     * variable denotes the step
+                                     * length to be used for
+                                     * that. It's value is set to
+                                     * 1e-6.
+                                     */
+    static const double fd_step_length;
+
                                     /**
                                      * Prepare internal data
                                      * structures and fill in values
@@ -367,5 +1862,207 @@ class FiniteElement : public FiniteElementBase<dim>
 };
 
 /*@}*/
+//----------------------------------------------------------------------//
+
+template <int dim>  
+inline
+std::pair<unsigned int,unsigned int>
+FiniteElement<dim>::system_to_component_index (const unsigned int index) const
+{
+  Assert (index < system_to_component_table.size(),
+        ExcIndexRange(index, 0, system_to_component_table.size()));
+  Assert (is_primitive (index),
+         typename FiniteElement<dim>::ExcShapeFunctionNotPrimitive(index));
+  return system_to_component_table[index];
+}
+
+template <int dim>  
+inline
+unsigned int
+FiniteElement<dim>::component_to_system_index (const unsigned int component,
+                                                   const unsigned int index) const
+{
+   std::vector< std::pair<unsigned int, unsigned int> >::const_iterator
+      it = std::find(system_to_component_table.begin(), system_to_component_table.end(),
+                     std::pair<unsigned int, unsigned int>(component, index));
+
+   Assert(it != system_to_component_table.end(), ExcComponentIndexInvalid(component, index));
+   return std::distance(system_to_component_table.begin(), it);
+}
+
+
+
+template <int dim>  
+inline
+std::pair<unsigned int,unsigned int>
+FiniteElement<dim>::face_system_to_component_index (const unsigned int index) const
+{
+  Assert(index < face_system_to_component_table.size(),
+        ExcIndexRange(index, 0, face_system_to_component_table.size()));
+
+                                   // in debug mode, check whether the
+                                   // function is primitive, since
+                                   // otherwise the result may have no
+                                   // meaning
+                                   //
+                                   // since the primitivity tables are
+                                   // all geared towards cell dof
+                                   // indices, rather than face dof
+                                   // indices, we have to work a
+                                   // little bit...
+                                   //
+                                   // in 1d, the face index is equal
+                                   // to the cell index
+  Assert (((dim == 1) && is_primitive(index))
+          ||
+                                           // in 2d, construct it like
+                                           // this:
+          ((dim == 2) &&
+           is_primitive (index < (GeometryInfo<2>::vertices_per_face *
+                                  this->dofs_per_vertex)
+                         ?
+                         index
+                         :
+                         GeometryInfo<2>::vertices_per_cell *
+                         this->dofs_per_vertex +
+                         (index -
+                          GeometryInfo<2>::vertices_per_face *
+                          this->dofs_per_vertex)))
+          ||
+                                           // likewise in 3d, but more
+                                           // complicated
+          ((dim == 3) &&
+           is_primitive (index < (GeometryInfo<3>::vertices_per_face *
+                                  this->dofs_per_vertex)
+                         ?
+                         index
+                         :
+                         (index < (GeometryInfo<3>::vertices_per_face *
+                                   this->dofs_per_vertex
+                                   +
+                                   GeometryInfo<3>::lines_per_face *
+                                   this->dofs_per_line)
+                          ?
+                          GeometryInfo<3>::vertices_per_cell *
+                          this->dofs_per_vertex +
+                          (index -
+                           GeometryInfo<3>::vertices_per_face *
+                           this->dofs_per_vertex)
+                          :
+                          GeometryInfo<3>::vertices_per_cell *
+                          this->dofs_per_vertex +
+                          GeometryInfo<3>::lines_per_cell *
+                          this->dofs_per_line +
+                          (index -
+                           GeometryInfo<3>::vertices_per_face *
+                           this->dofs_per_vertex
+                           -
+                           GeometryInfo<3>::lines_per_face *
+                           this->dofs_per_line)))),
+          typename FiniteElement<dim>::ExcShapeFunctionNotPrimitive(index));
+
+  return face_system_to_component_table[index];
+}
+
+
+
+template <int dim>  
+inline
+std::pair<std::pair<unsigned int,unsigned int>,unsigned int>
+FiniteElement<dim>::system_to_base_index (const unsigned int index) const
+{
+  Assert (index < system_to_base_table.size(),
+        ExcIndexRange(index, 0, system_to_base_table.size()));
+  return system_to_base_table[index];
+}
+
+
+
+
+template <int dim>  
+inline
+std::pair<std::pair<unsigned int,unsigned int>,unsigned int>
+FiniteElement<dim>::face_system_to_base_index (const unsigned int index) const
+{
+  Assert(index < face_system_to_base_table.size(),
+        ExcIndexRange(index, 0, face_system_to_base_table.size()));
+  return face_system_to_base_table[index];
+}
+
+
+
+template <int dim>  
+inline
+std::pair<unsigned int,unsigned int>
+FiniteElement<dim>::component_to_base (const unsigned int index) const
+{
+  Assert(index < component_to_base_table.size(),
+        ExcIndexRange(index, 0, component_to_base_table.size()));
+
+  return component_to_base_table[index];
+}
+
+
+template <int dim>
+inline
+bool
+FiniteElement<dim>::restriction_is_additive (const unsigned int index) const
+{
+  Assert(index < this->dofs_per_cell,
+        ExcIndexRange(index, 0, this->dofs_per_cell));
+  return restriction_is_additive_flags[index];
+}
+
+
+template <int dim>
+inline
+const std::vector<bool> &
+FiniteElement<dim>::get_nonzero_components (const unsigned int i) const
+{
+  Assert (i < this->dofs_per_cell, ExcIndexRange (i, 0, this->dofs_per_cell));
+  return nonzero_components[i];
+}
+
+
+
+template <int dim>
+inline
+unsigned int
+FiniteElement<dim>::n_nonzero_components (const unsigned int i) const
+{
+  Assert (i < this->dofs_per_cell, ExcIndexRange (i, 0, this->dofs_per_cell));
+  return n_nonzero_components_table[i];
+}
+
+
+
+template <int dim>
+inline
+bool
+FiniteElement<dim>::is_primitive (const unsigned int i) const
+{
+  Assert (i < this->dofs_per_cell, ExcIndexRange (i, 0, this->dofs_per_cell));
+
+                                  // return primitivity of a shape
+                                  // function by checking whether it
+                                  // has more than one non-zero
+                                  // component or not. we could cache
+                                  // this value in an array of bools,
+                                  // but accessing a bit-vector (as
+                                  // std::vector<bool> is) is
+                                  // probably more expensive than
+                                  // just comparing against 1
+  return (n_nonzero_components_table[i] == 1);
+}
+
+
+template <int dim>
+inline
+bool
+FiniteElement<dim>::is_primitive () const
+{
+  return cached_primitivity;
+}
+
 
 #endif
index 1e5c1ac93b3af2d58da0c01d5d0b9283947f20da..33837eefe3b84e9bb2d3dbf68c34093689cdd914 100644 (file)
@@ -35,14 +35,14 @@ template<int dim> class FESystem;
 
 /**
  * Dimension independent data for finite elements. See the derived
- * class FiniteElementBase class for information on its use. All
+ * class FiniteElement class for information on its use. All
  * its data are available to the implementation in a concrete finite
  * element class.
  *
  * Remark on a change in implementation: it is now wrong to cast a
  * pointer to FiniteElement to a pointer to FiniteElementData and
  * delete it. The virtual destructor has been moved up. In a later
- * version, FiniteElementData and FiniteElementBase should be private
+ * version, FiniteElementData and FiniteElement should be private
  * base classes of FiniteElement.
  *
  * @author Wolfgang Bangerth, Guido Kanschat, 1998, 1999, 2000, 2001, 2003, 2005
@@ -560,1546 +560,6 @@ class FiniteElementData
  *
  * @author Wolfgang Bangerth, 1998, 2002; Ralf Hartmann, Guido Kanschat, 2001
  */
-template <int dim>
-class FiniteElementBase : public Subscriptor,
-                         public FiniteElementData<dim>
-{
-  public:
-                                  /**
-                                   * Base class for internal data.
-                                   * Adds data for second derivatives to
-                                   * Mapping::InternalDataBase()
-                                   *
-                                   * For information about the
-                                   * general purpose of this class,
-                                   * see the documentation of the
-                                   * base class.
-                                   *
-                                   * @author Guido Kanschat, 2001
-                                   */
-  class InternalDataBase : public Mapping<dim>::InternalDataBase
-    {
-      public:      
-                                        /**
-                                         * Destructor. Needed to
-                                         * avoid memory leaks with
-                                         * difference quotients.
-                                         */
-       virtual ~InternalDataBase ();
-
-                                        /**
-                                         * Initialize some pointers
-                                         * used in the computation of
-                                         * second derivatives by
-                                         * finite differencing of
-                                         * gradients.
-                                         */
-       void initialize_2nd (const FiniteElement<dim> *element,
-                            const Mapping<dim>       &mapping,
-                            const Quadrature<dim>    &quadrature);
-       
-                                        /**
-                                         * Storage for @p FEValues
-                                         * objects needed to
-                                         * approximate second
-                                         * derivatives.
-                                         *
-                                         * The ordering is <tt>p+hx</tt>,
-                                         * <tt>p+hy</tt>, <tt>p+hz</tt>,
-                                         * @p p-hx, @p p-hy,
-                                         * @p p-hz, where unused
-                                         * entries in lower dimensions
-                                         * are missing.
-                                         */
-       std::vector<FEValues<dim>*> differences;
-    };
-  
-                                    /**
-                                     * Construct an object of this
-                                     * type. You have to set some
-                                     * member variables, for example
-                                     * some matrices, explicitly
-                                     * after calling this base class'
-                                     * constructor. For this see the
-                                     * existing finite element
-                                     * classes. For the second and
-                                     * third parameter of this
-                                     * constructor, see the documentation
-                                     * of the respective member
-                                     * variables.
-                                     *
-                                     * @note Both vector parameters
-                                     * should have length
-                                     * <tt>dofs_per_cell</tt>. Nevertheless,
-                                     * it is allowed to use vectors
-                                     * of length one. In this case,
-                                     * the vector is resized to the
-                                     * correct length and filled with
-                                     * the entry value.
-                                     */
-    FiniteElementBase (const FiniteElementData<dim> &fe_data,
-                      const std::vector<bool> &restriction_is_additive_flags,
-                      const std::vector<std::vector<bool> > &nonzero_components);
-
-                                    /**
-                                     * Return a string that uniquely
-                                     * identifies a finite
-                                     * element. The general
-                                     * convention is that this is the
-                                     * class name, followed by the
-                                     * space dimension in angle
-                                     * brackets, and the polynomial
-                                     * degree and whatever else is
-                                     * necessary in parentheses. For
-                                     * example, <tt>FE_Q<2>(3)</tt> is the
-                                     * value returned for a cubic
-                                     * element in 2d.
-                                     *
-                                     * Systems of elements have their
-                                     * own naming convention, see the
-                                     * FESystem class.
-                                     */
-    virtual std::string get_name () const = 0;
-
-                                    /**
-                                     * @name Shape function access
-                                     * @{
-                                     */
-    
-                                    /**
-                                     * Return the value of the
-                                     * @p ith shape function at the
-                                     * point @p p. @p p is a point
-                                     * on the reference element. If
-                                     * the finite element is
-                                     * vector-valued, then return the
-                                     * value of the only non-zero
-                                     * component of the vector value
-                                     * of this shape function. If the
-                                     * shape function has more than
-                                     * one non-zero component (which
-                                     * we refer to with the term
-                                     * non-primitive), then derived
-                                     * classes implementing this
-                                     * function should throw an
-                                     * exception of type
-                                     * @p ExcShapeFunctionNotPrimitive. In
-                                     * that case, use the
-                                     * shape_value_component()
-                                     * function.
-                                     *
-                                     * An
-                                     * @p ExcUnitShapeValuesDoNotExist
-                                     * is thrown if the shape values
-                                     * of the @p FiniteElement under
-                                     * consideration depends on the
-                                     * shape of the cell in real
-                                     * space.
-                                     */
-    virtual double shape_value (const unsigned int  i,
-                               const Point<dim>   &p) const;
-
-                                    /**
-                                     * Just like for @p shape_value,
-                                     * but this function will be
-                                     * called when the shape function
-                                     * has more than one non-zero
-                                     * vector component. In that
-                                     * case, this function should
-                                     * return the value of the
-                                     * @p component-th vector
-                                     * component of the @p ith shape
-                                     * function at point @p p.
-                                     */
-    virtual double shape_value_component (const unsigned int i,
-                                         const Point<dim>   &p,
-                                         const unsigned int component) const;
-    
-                                    /**
-                                     * Return the gradient of the
-                                     * @p ith shape function at the
-                                     * point @p p. @p p is a point
-                                     * on the reference element, and
-                                     * likewise the gradient is the
-                                     * gradient on the unit cell with
-                                     * respect to unit cell
-                                     * coordinates. If
-                                     * the finite element is
-                                     * vector-valued, then return the
-                                     * value of the only non-zero
-                                     * component of the vector value
-                                     * of this shape function. If the
-                                     * shape function has more than
-                                     * one non-zero component (which
-                                     * we refer to with the term
-                                     * non-primitive), then derived
-                                     * classes implementing this
-                                     * function should throw an
-                                     * exception of type
-                                     * @p ExcShapeFunctionNotPrimitive. In
-                                     * that case, use the
-                                     * shape_grad_component()
-                                     * function.
-                                     *
-                                     * An
-                                     * @p ExcUnitShapeValuesDoNotExist
-                                     * is thrown if the shape values
-                                     * of the @p FiniteElement under
-                                     * consideration depends on the
-                                     * shape of the cell in real
-                                     * space.
-                                     */
-    virtual Tensor<1,dim> shape_grad (const unsigned int  i,
-                                     const Point<dim>   &p) const;
-
-                                    /**
-                                     * Just like for @p shape_grad,
-                                     * but this function will be
-                                     * called when the shape function
-                                     * has more than one non-zero
-                                     * vector component. In that
-                                     * case, this function should
-                                     * return the gradient of the
-                                     * @p component-th vector
-                                     * component of the @p ith shape
-                                     * function at point @p p.
-                                     */
-    virtual Tensor<1,dim> shape_grad_component (const unsigned int i,
-                                               const Point<dim>   &p,
-                                               const unsigned int component) const;
-
-                                    /**
-                                     * Return the tensor of second
-                                     * derivatives of the @p ith
-                                     * shape function at point @p p
-                                     * on the unit cell. The
-                                     * derivatives are derivatives on
-                                     * the unit cell with respect to
-                                     * unit cell coordinates. If
-                                     * the finite element is
-                                     * vector-valued, then return the
-                                     * value of the only non-zero
-                                     * component of the vector value
-                                     * of this shape function. If the
-                                     * shape function has more than
-                                     * one non-zero component (which
-                                     * we refer to with the term
-                                     * non-primitive), then derived
-                                     * classes implementing this
-                                     * function should throw an
-                                     * exception of type
-                                     * @p ExcShapeFunctionNotPrimitive. In
-                                     * that case, use the
-                                     * shape_grad_grad_component()
-                                     * function.
-                                     *
-                                     * An
-                                     * @p ExcUnitShapeValuesDoNotExist
-                                     * is thrown if the shape values
-                                     * of the @p FiniteElement under
-                                     * consideration depends on the
-                                     * shape of the cell in real
-                                     * space.
-                                     */
-    virtual Tensor<2,dim> shape_grad_grad (const unsigned int  i,
-                                          const Point<dim>   &p) const;
-
-                                    /**
-                                     * Just like for @p shape_grad_grad,
-                                     * but this function will be
-                                     * called when the shape function
-                                     * has more than one non-zero
-                                     * vector component. In that
-                                     * case, this function should
-                                     * return the gradient of the
-                                     * @p component-th vector
-                                     * component of the @p ith shape
-                                     * function at point @p p.
-                                     */
-    virtual Tensor<2,dim> shape_grad_grad_component (const unsigned int i,
-                                                    const Point<dim>   &p,
-                                                    const unsigned int component) const;
-                                    //@}
-                                    /**
-                                     * @name Transfer and constraint matrices
-                                     * @{
-                                     */
-    
-                                    /**
-                                     * Projection from a fine grid
-                                     * space onto a coarse grid
-                                     * space. If this projection
-                                     * operator is associated with a
-                                     * matrix @p P, then the
-                                     * restriction of this matrix
-                                     * @p P_i to a single child cell
-                                     * is returned here.
-                                     *
-                                     * The matrix @p P is the
-                                     * concatenation or the sum of
-                                     * the cell matrices @p P_i,
-                                     * depending on the
-                                     * @p restriction_is_additive_flags
-                                     * given to the constructor. This
-                                     * distinguishes interpolation
-                                     * (concatenation) and projection
-                                     * with respect to scalar
-                                     * products (summation).
-                                     *
-                                     * Row and column indices are
-                                     * related to coarse grid and
-                                     * fine grid spaces,
-                                     * respectively, consistent with
-                                     * the definition of the
-                                     * associated operator.
-                                     *
-                                     * If projection matrices are not
-                                     * implemented in the derived
-                                     * finite element class, this
-                                     * function aborts with
-                                     * @p ExcProjectionVoid.
-                                     */
-    const FullMatrix<double> &
-    get_restriction_matrix (const unsigned int child) const;
-
-                                    /**
-                                     * Embedding matrix between grids.
-                                     * 
-                                     * The identity operator from a
-                                     * coarse grid space into a fine
-                                     * grid space is associated with
-                                     * a matrix @p P. The
-                                     * restriction of this matrix @p P_i to
-                                     * a single child cell is
-                                     * returned here.
-                                     *
-                                     * The matrix @p P is the
-                                     * concatenation, not the sum of
-                                     * the cell matrices
-                                     * @p P_i. That is, if the same
-                                     * non-zero entry <tt>j,k</tt> exists
-                                     * in in two different child
-                                     * matrices @p P_i, the value
-                                     * should be the same in both
-                                     * matrices and it is copied into
-                                     * the matrix @p P only once.
-                                     *
-                                     * Row and column indices are
-                                     * related to fine grid and
-                                     * coarse grid spaces,
-                                     * respectively, consistent with
-                                     * the definition of the
-                                     * associated operator.
-                                     *
-                                     * These matrices are used by
-                                     * routines assembling the
-                                     * prolongation matrix for
-                                     * multi-level methods.  Upon
-                                     * assembling the transfer matrix
-                                     * between cells using this
-                                     * matrix array, zero elements in
-                                     * the prolongation matrix are
-                                     * discarded and will not fill up
-                                     * the transfer matrix.
-                                     *
-                                     * If projection matrices are not
-                                     * implemented in the derived
-                                     * finite element class, this
-                                     * function aborts with
-                                     * @p ExcEmbeddingVoid. You can
-                                     * check whether this is the case
-                                     * by calling the
-                                     * prolongation_is_implemented().
-                                     */
-    const FullMatrix<double> &
-    get_prolongation_matrix (const unsigned int child) const;
-
-                                     /**
-                                      * Return whether this element implements
-                                      * its prolongation matrices. The return
-                                      * value also indicates whether a call to
-                                      * the @p get_prolongation_matrix
-                                      * function will generate an error or
-                                      * not.
-                                      *
-                                      * This function is mostly here in order
-                                      * to allow us to write more efficient
-                                      * test programs which we run on all
-                                      * kinds of weird elements, and for which
-                                      * we simply need to exclude certain
-                                      * tests in case something is not
-                                      * implemented. It will in general
-                                      * probably not be a great help in
-                                      * applications, since there is not much
-                                      * one can do if one needs these features
-                                      * and they are not implemented. This
-                                      * function could be used to check
-                                      * whether a call to
-                                      * <tt>get_prolongation_matrix()</tt> will
-                                      * succeed; however, one then still needs
-                                      * to cope with the lack of information
-                                      * this just expresses.
-                                      */
-    bool prolongation_is_implemented () const;
-
-                                     /**
-                                      * Return whether this element implements
-                                      * its restriction matrices. The return
-                                      * value also indicates whether a call to
-                                      * the @p get_restriction_matrix
-                                      * function will generate an error or
-                                      * not.
-                                      *
-                                      * This function is mostly here in order
-                                      * to allow us to write more efficient
-                                      * test programs which we run on all
-                                      * kinds of weird elements, and for which
-                                      * we simply need to exclude certain
-                                      * tests in case something is not
-                                      * implemented. It will in general
-                                      * probably not be a great help in
-                                      * applications, since there is not much
-                                      * one can do if one needs these features
-                                      * and they are not implemented. This
-                                      * function could be used to check
-                                      * whether a call to
-                                      * <tt>get_restriction_matrix()</tt> will
-                                      * succeed; however, one then still needs
-                                      * to cope with the lack of information
-                                      * this just expresses.
-                                      */
-    bool restriction_is_implemented () const;
-
-                                    /**
-                                     * Access the
-                                     * @p restriction_is_additive_flag
-                                     * field. See there for more
-                                     * information on its contents.
-                                     *
-                                     * The index must be between zero
-                                     * and the number of shape
-                                     * functions of this element.
-                                     */
-    bool restriction_is_additive (const unsigned int index) const;
-
-                                    /**
-                                     * Return a readonly reference to
-                                     * the matrix which describes the
-                                     * constraints at the interface
-                                     * between a refined and an
-                                     * unrefined cell.
-                                     * 
-                                     * The matrix is obviously empty
-                                     * in only one space dimension,
-                                     * since there are no constraints
-                                     * then.
-                                     *
-                                     * Note that some finite elements
-                                     * do not (yet) implement hanging
-                                     * node constraints. If this is
-                                     * the case, then this function
-                                     * will generate an exception,
-                                     * since no useful return value
-                                     * can be generated. If you
-                                     * should have a way to live with
-                                     * this, then you might want to
-                                     * use the
-                                     * @p constraints_are_implemented
-                                     * function to check up front
-                                     * whethehr this function will
-                                     * succeed or generate the
-                                     * exception.
-                                     */
-    const FullMatrix<double> & constraints () const;
-
-                                     /**
-                                      * Return whether this element
-                                      * implements its hanging node
-                                      * constraints. The return value
-                                      * also indicates whether a call
-                                      * to the @p constraint function
-                                      * will generate an error or not.
-                                      *
-                                      * This function is mostly here
-                                      * in order to allow us to write
-                                      * more efficient test programs
-                                      * which we run on all kinds of
-                                      * weird elements, and for which
-                                      * we simply need to exclude
-                                      * certain tests in case hanging
-                                      * node constraints are not
-                                      * implemented. It will in
-                                      * general probably not be a
-                                      * great help in applications,
-                                      * since there is not much one
-                                      * can do if one needs hanging
-                                      * node constraints and they are
-                                      * not implemented. This function
-                                      * could be used to check whether
-                                      * a call to <tt>constraints()</tt>
-                                      * will succeed; however, one
-                                      * then still needs to cope with
-                                      * the lack of information this
-                                      * just expresses.
-                                      */
-    bool constraints_are_implemented () const;
-
-                                    /**
-                                     * Return the matrix
-                                     * interpolating from the given
-                                     * finite element to the present
-                                     * one. The size of the matrix is
-                                     * then @p dofs_per_cell times
-                                     * <tt>source.dofs_per_cell</tt>.
-                                     *
-                                     * Derived elements will have to
-                                     * implement this function. They
-                                     * may only provide interpolation
-                                     * matrices for certain source
-                                     * finite elements, for example
-                                     * those from the same family. If
-                                     * they don't implement
-                                     * interpolation from a given
-                                     * element, then they must throw
-                                     * an exception of type
-                                     * FiniteElementBase<dim>::ExcInterpolationNotImplemented.
-                                     */
-    virtual void
-    get_interpolation_matrix (const FiniteElementBase<dim> &source,
-                             FullMatrix<double>           &matrix) const;
-                                    //@}
-    
-                                    /**
-                                     * Comparison operator. We also
-                                     * check for equality of the
-                                     * constraint matrix, which is
-                                     * quite an expensive operation.
-                                     * Do therefore use this function
-                                     * with care, if possible only
-                                     * for debugging purposes.
-                                     *
-                                     * Since this function is not
-                                     * that important, we avoid an
-                                     * implementational question
-                                     * about comparing arrays and do
-                                     * not compare the matrix arrays
-                                     * @p restriction and
-                                     * @p prolongation.
-                                     */
-    bool operator == (const FiniteElementBase<dim> &) const;
-
-                                    /**
-                                     * @name Index computations
-                                     * @{
-                                     */
-                                    /**
-                                     * Compute vector component and
-                                     * index of this shape function
-                                     * within the shape functions
-                                     * corresponding to this
-                                     * component from the index of a
-                                     * shape function within this
-                                     * finite element.
-                                     *
-                                     * If the element is scalar, then
-                                     * the component is always zero,
-                                     * and the index within this
-                                     * component is equal to the
-                                     * overall index.
-                                     *
-                                     * If the shape function
-                                     * referenced has more than one
-                                     * non-zero component, then it
-                                     * cannot be associated with one
-                                     * vector component, and an
-                                     * exception of type
-                                     * @p ExcShapeFunctionNotPrimitive
-                                     * will be raised.
-                                     *
-                                     * Note that if the element is
-                                     * composed of other (base)
-                                     * elements, and a base element
-                                     * has more than one component
-                                     * but all its shape functions
-                                     * are primitive (i.e. are
-                                     * non-zero in only one
-                                     * component), then this mapping
-                                     * contains valid
-                                     * information. However, the
-                                     * index of a shape function of
-                                     * this element within one
-                                     * component (i.e. the second
-                                     * number of the respective entry
-                                     * of this array) does not
-                                     * indicate the index of the
-                                     * respective shape function
-                                     * within the base element (since
-                                     * that has more than one
-                                     * vector-component). For this
-                                     * information, refer to the
-                                     * @p system_to_base_table field
-                                     * and the
-                                     * @p system_to_base_index
-                                     * function.
-                                     */
-    std::pair<unsigned int, unsigned int>
-    system_to_component_index (const unsigned int index) const;
-
-                                    /**
-                                     * Compute the shape function for
-                                     * the given vector component and
-                                     * index.
-                                     *
-                                     * If the element is scalar, then
-                                     * the component must be zero,
-                                     * and the index within this
-                                     * component is equal to the
-                                     * overall index.
-                                     *
-                                     * This is the opposite operation
-                                     * from the @p system_to_component_index
-                                     * function.
-                                     */
-   unsigned int component_to_system_index(const unsigned int component,
-                                          const unsigned int index) const;
-  
-                                    /**
-                                     * Same as above, but do it for
-                                     * shape functions and their
-                                     * indices on a face.
-                                     */
-    std::pair<unsigned int, unsigned int>
-    face_system_to_component_index (const unsigned int index) const;
-
-                                     /**
-                                      * Return for shape function
-                                      * @p index the base element it
-                                      * belongs to, the number of the
-                                      * copy of this base element
-                                      * (which is between zero and the
-                                      * multiplicity of this element),
-                                      * and the index of this shape
-                                      * function within this base
-                                      * element.
-                                      *
-                                      * If the element is not composed of
-                                     * others, then base and instance
-                                     * are always zero, and the index
-                                     * is equal to the number of the
-                                     * shape function. If the element
-                                     * is composed of single
-                                     * instances of other elements
-                                     * (i.e. all with multiplicity
-                                     * one) all of which are scalar,
-                                     * then base values and dof
-                                     * indices within this element
-                                     * are equal to the
-                                     * @p system_to_component_table. It
-                                     * differs only in case the
-                                     * element is composed of other
-                                     * elements and at least one of
-                                     * them is vector-valued itself.
-                                     *
-                                     * This function returns valid
-                                     * values also in the case of
-                                     * vector-valued
-                                     * (i.e. non-primitive) shape
-                                     * functions, in contrast to the
-                                     * @p system_to_component_index
-                                     * function.
-                                      */
-    std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
-    system_to_base_index (const unsigned int index) const;
-
-                                     /**
-                                      * Same as
-                                      * @p system_to_base_index, but
-                                      * for degrees of freedom located
-                                      * on a face.
-                                      */
-    std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
-    face_system_to_base_index (const unsigned int index) const;
-    
-                                    /**
-                                     * Given a vector component,
-                                     * return an index which base
-                                     * element implements this
-                                     * component, and which vector
-                                     * component in this base element
-                                     * this is. This information is
-                                     * only of interest for
-                                     * vector-valued finite elements
-                                     * which are composed of several
-                                     * sub-elements. In that case,
-                                     * one may want to obtain
-                                     * information about the element
-                                     * implementing a certain vector
-                                     * component, which can be done
-                                     * using this function and the
-                                     * FESystem::@p base_element
-                                     * function.
-                                     *
-                                     * If this is a scalar finite
-                                     * element, then the return value
-                                     * is always equal to a pair of
-                                     * zeros.
-                                     */
-    std::pair<unsigned int,unsigned int>
-    component_to_base (const unsigned int component) const;
-                                    //@}
-    
-                                    /**
-                                     * @name Support points and interpolation
-                                     * @{
-                                     */
-    
-                                    /**
-                                     * Return the support points of
-                                     * the trial functions on the
-                                     * unit cell, if the derived
-                                     * finite element defines some.
-                                     * Finite elements that allow
-                                     * some kind of interpolation
-                                     * operation usually have support
-                                     * points. On the other hand,
-                                     * elements that define their
-                                     * degrees of freedom by, for
-                                     * example, moments on faces, or
-                                     * as derivatives, don't have
-                                     * support points. In that case,
-                                     * the returned field is empty.
-                                     *
-                                     * If the finite element defines
-                                     * support points, then their
-                                     * number equals the number of
-                                     * degrees of freedom of the
-                                     * element.  The order of points
-                                     * in the array matches that
-                                     * returned by the
-                                     * <tt>cell->get_dof_indices</tt>
-                                     * function.
-                                     *
-                                     * See the class documentation
-                                     * for details on support points.
-                                     */
-    const std::vector<Point<dim> > &
-    get_unit_support_points () const;    
-
-                                    /**
-                                     * Return whether a finite
-                                     * element has defined support
-                                     * points. If the result is true,
-                                     * then a call to the
-                                     * @p get_unit_support_points
-                                     * yields a non-empty array.
-                                     *
-                                     * The result may be false if an
-                                     * element is not defined by
-                                     * interpolating shape functions,
-                                     * for example by P-elements on
-                                     * quadrilaterals. It will
-                                     * usually only be true if the
-                                     * element constructs its shape
-                                     * functions by the requirement
-                                     * that they be one at a certain
-                                     * point and zero at all the
-                                     * points associated with the
-                                     * other shape functions.
-                                     *
-                                     * In composed elements (i.e. for
-                                     * the FESystem class, the
-                                     * result will be true if all all
-                                     * the base elements have defined
-                                     * support points.
-                                     */
-    bool has_support_points () const;
-
-                                     /**
-                                      * Return the position of the
-                                      * support point of the
-                                      * @p indexth shape function. If
-                                      * it does not exist, raise an
-                                      * exception.
-                                      *
-                                      * The default implementation
-                                      * simply returns the respective
-                                      * element from the array you get
-                                      * from
-                                      * get_unit_support_points(),
-                                      * but derived elements may
-                                      * overload this function. In
-                                      * particular, note that the
-                                      * FESystem class overloads
-                                      * it so that it can return the
-                                      * support points of individual
-                                      * base elements, of not all the
-                                      * base elements define support
-                                      * points. In this way, you can
-                                      * still ask for certain support
-                                      * points, even if
-                                      * @p get_unit_support_points
-                                      * only returns an empty array.
-                                      */
-    virtual
-    Point<dim>
-    unit_support_point (const unsigned int index) const;
-    
-                                    /**
-                                     * Return the support points of
-                                     * the trial functions on the
-                                     * unit face, if the derived
-                                     * finite element defines some.
-                                     * Finite elements that allow
-                                     * some kind of interpolation
-                                     * operation usually have support
-                                     * points. On the other hand,
-                                     * elements that define their
-                                     * degrees of freedom by, for
-                                     * example, moments on faces, or
-                                     * as derivatives, don't have
-                                     * support points. In that case,
-                                     * the returned field is empty
-                                     *
-                                     * Note that elements that have
-                                     * support points need not
-                                     * necessarily have some on the
-                                     * faces, even if the
-                                     * interpolation points are
-                                     * located physically on a
-                                     * face. For example, the
-                                     * discontinuous elements have
-                                     * interpolation points on the
-                                     * vertices, and for higher
-                                     * degree elements also on the
-                                     * faces, but they are not
-                                     * defined to be on faces since
-                                     * in that case degrees of
-                                     * freedom from both sides of a
-                                     * face (or from all adjacent
-                                     * elements to a vertex) would be
-                                     * identified with each other,
-                                     * which is not what we would
-                                     * like to have). Logically,
-                                     * these degrees of freedom are
-                                     * therefore defined to belong to
-                                     * the cell, rather than the face
-                                     * or vertex. In that case, the
-                                     * returned element would
-                                     * therefore have length zero.
-                                     *
-                                     * If the finite element defines
-                                     * support points, then their
-                                     * number equals the number of
-                                     * degrees of freedom on the face
-                                     * (@p dofs_per_face). The order
-                                     * of points in the array matches
-                                     * that returned by the
-                                     * <tt>cell->get_dof_indices</tt>
-                                     * function.
-                                     *
-                                     * See the class documentation
-                                     * for details on support points.
-                                     */
-    const std::vector<Point<dim-1> > &
-    get_unit_face_support_points () const;    
-
-                                    /**
-                                     * Return whether a finite
-                                     * element has defined support
-                                     * points on faces. If the result
-                                     * is true, then a call to the
-                                     * @p get_unit_support_points
-                                     * yields a non-empty array.
-                                     *
-                                     * For more information, see the
-                                     * documentation for the
-                                     * has_support_points()
-                                     * function.
-                                     */
-    bool has_face_support_points () const;
-
-                                     /**
-                                      * The function corresponding to
-                                      * the unit_support_point()
-                                      * function, but for faces. See
-                                      * there for more information.
-                                      */
-    virtual
-    Point<dim-1>
-    unit_face_support_point (const unsigned int index) const;
-    
-                                    /**
-                                     * Return a support point vector
-                                     * for generalized interpolation.
-                                     */
-    const std::vector<Point<dim> > &
-    get_generalized_support_points () const;    
-
-                                    /**
-                                     *
-                                     */
-    bool has_generalized_support_points () const;
-
-                                    /**
-                                     *
-                                     */
-    const std::vector<Point<dim-1> > &
-    get_generalized_face_support_points () const;
-
-                                    /**
-                                     * Return whether a finite
-                                     * element has defined support
-                                     * points on faces. If the result
-                                     * is true, then a call to the
-                                     * @p get_unit_support_points
-                                     * yields a non-empty array.
-                                     *
-                                     * For more information, see the
-                                     * documentation for the
-                                     * has_support_points()
-                                     * function.
-                                     */
-    bool has_generalized_face_support_points () const;
-
-                                    /**
-                                     * Interpolate a set of scalar
-                                     * values, computed in the
-                                     * generalized support points.
-                                     *
-                                     * @note This function is
-                                     * implemented in
-                                     * FiniteElementBase for the case
-                                     * that the element has support
-                                     * points. In this case, the
-                                     * resulting coefficients are
-                                     * just the values in the suport
-                                     * points. All other elements
-                                     * must reimplement it.
-                                     */
-    virtual void interpolate(std::vector<double>&       local_dofs,
-                            const std::vector<double>& values) const;
-      
-                                    /**
-                                     * Interpolate a set of vector
-                                     * values, computed in the
-                                     * generalized support points.
-                                     *
-                                     * Since a finite element often
-                                     * only interpolates part of a
-                                     * vector, <tt>offset</tt> is
-                                     * used to determine the first
-                                     * component of the vector to be
-                                     * interpolated. Maybe consider
-                                     * changing your data structures
-                                     * to use the next function.
-                                     */
-    virtual void interpolate(std::vector<double>&                local_dofs,
-                            const std::vector<Vector<double> >& values,
-                            unsigned int offset = 0) const;
-      
-                                    /**
-                                     * Interpolate a set of vector
-                                     * values, computed in the
-                                     * generalized support points.
-                                     */
-    virtual void interpolate(
-      std::vector<double>& local_dofs,
-      const VectorSlice<const std::vector<std::vector<double> > >& values) const;
-      
-                                    //@}
-    
-                                    /**
-                                     * Return in which of the vector
-                                     * components of this finite
-                                     * element the @p ithe shape
-                                     * function is non-zero. The
-                                     * length of the returned array
-                                     * is equal to the number of
-                                     * vector components of this
-                                     * element.
-                                     *
-                                     * For most finite element
-                                     * spaces, the result of this
-                                     * function will be a vector with
-                                     * exactly one element being
-                                     * @p true, since for most
-                                     * spaces the individual vector
-                                     * components are independent. In
-                                     * that case, the component with
-                                     * the single zero is also the
-                                     * first element of what
-                                     * <tt>system_to_component_index(i)</tt>
-                                     * returns.
-                                     *
-                                     * Only for those
-                                     * spaces that couple the
-                                     * components, for example to
-                                     * make a shape function
-                                     * divergence free, will there be
-                                     * more than one @p true entry.
-                                     */
-    const std::vector<bool> &
-    get_nonzero_components (const unsigned int i) const;
-
-                                    /**
-                                     * Return in how many vector
-                                     * components the @p ith shape
-                                     * function is non-zero. This
-                                     * value equals the number of
-                                     * entries equal to @p true in
-                                     * the result of the
-                                     * @p get_nonzero_components
-                                     * function.
-                                     *
-                                     * For most finite element
-                                     * spaces, the result will be
-                                     * equal to one. It is not equal
-                                     * to one only for those ansatz
-                                     * spaces for which vector-valued
-                                     * shape functions couple the
-                                     * individual components, for
-                                     * example in order to make them
-                                     * divergence-free.
-                                     */
-    unsigned int
-    n_nonzero_components (const unsigned int i) const;
-
-                                    /**
-                                     * Return whether the @p ith
-                                     * shape function is primitive in
-                                     * the sense that the shape
-                                     * function is non-zero in only
-                                     * one vector
-                                     * component. Non-primitive shape
-                                     * functions would then, for
-                                     * example, be those of
-                                     * divergence free ansatz spaces,
-                                     * in which the individual vector
-                                     * components are coupled.
-                                     *
-                                     * The result of the function is
-                                     * @p true if and only if the
-                                     * result of
-                                     * <tt>n_nonzero_components(i)</tt> is
-                                     * equal to one.
-                                     */
-    bool
-    is_primitive (const unsigned int i) const;
-
-                                    /**
-                                     * Return whether the entire
-                                     * finite element is primitive,
-                                     * in the sense that all its
-                                     * shape functions are
-                                     * primitive. If the finite
-                                     * element is scalar, then this
-                                     * is always the case.
-                                     *
-                                     * Since this is an extremely
-                                     * common operation, the result
-                                     * is cached in the
-                                     * @p cached_primitivity
-                                     * variable which is computed in
-                                     * the constructor.
-                                     */
-    bool
-    is_primitive () const;
-    
-                                    /**
-                                     * Determine an estimate for the
-                                     * memory consumption (in bytes)
-                                     * of this object.
-                                     *
-                                     * This function is not
-                                     * virtual. Use a
-                                     * FiniteElement object to
-                                     * get the actual size of a
-                                     * concrete element.
-                                     */
-    unsigned int memory_consumption () const;
-
-                                    /**
-                                     * Exception
-                                     *
-                                     * @ingroup Exceptions
-                                     */
-    DeclException1 (ExcShapeFunctionNotPrimitive,
-                   int,
-                   << "The shape function with index " << arg1
-                   << " is not primitive, i.e. it is vector-valued and "
-                   << "has more than one non-zero vector component. This "
-                   << "function cannot be called for these shape functions. "
-                   << "Maybe you want to use the same function with the "
-                   << "_component suffix?");
-                                    /**
-                                     * Exception
-                                     *
-                                     * @ingroup Exceptions
-                                     */
-    DeclException0 (ExcFENotPrimitive);
-                                    /**
-                                     * Exception
-                                     *
-                                     * @ingroup Exceptions
-                                     */
-    DeclException0 (ExcUnitShapeValuesDoNotExist);
-
-                                    /**
-                                     * Attempt to access support
-                                     * points of a finite element
-                                     * which is not Lagrangian.
-                                     *
-                                     * @ingroup Exceptions
-                                     */
-    DeclException0 (ExcFEHasNoSupportPoints);
-
-                                    /**
-                                     * Attempt to access embedding
-                                     * matrices of a finite element
-                                     * which did not implement these
-                                     * matrices.
-                                     *
-                                     * @ingroup Exceptions
-                                     */
-    DeclException0 (ExcEmbeddingVoid);
-    
-                                    /**
-                                     * Attempt to access restriction
-                                     * matrices of a finite element
-                                     * which did not implement these
-                                     * matrices.
-                                     *
-                                     * Exception
-                                     * @ingroup Exceptions
-                                     */
-    DeclException0 (ExcProjectionVoid);
-    
-                                    /**
-                                     * Attempt to access constraint
-                                     * matrices of a finite element
-                                     * which did not implement these
-                                     * matrices.
-                                     *
-                                     * Exception
-                                     * @ingroup Exceptions
-                                     */
-    DeclException0 (ExcConstraintsVoid);
-    
-                                    /**
-                                     * Exception
-                                     * @ingroup Exceptions
-                                     */
-    DeclException2 (ExcWrongInterfaceMatrixSize,
-                   int, int,
-                   << "The interface matrix has a size of " << arg1
-                   << "x" << arg2
-                   << ", which is not reasonable in the present dimension.");
-                                    /**
-                                     * Exception
-                                     * @ingroup Exceptions
-                                     */
-    DeclException2 (ExcComponentIndexInvalid,
-                   int, int,
-                   << "The component-index pair (" << arg1 << ", " << arg2
-                   << ") is invalid, i.e. non-existent");
-                                     /**
-                                      * Exception
-                                     * @ingroup Exceptions
-                                      */
-    DeclException0 (ExcInterpolationNotImplemented);
-    
-  protected:  
-                                    /**
-                                     * Array of projection matrices. See
-                                     * get_restriction_matrix() above.
-                                     *
-                                     * Matrices in this array are
-                                     * automatically initialized to
-                                     * correct size. If the derived
-                                     * finite element class does not
-                                     * implement these matrices, they
-                                     * should be resized to zero
-                                     * size.
-                                     */
-    FullMatrix<double> restriction[GeometryInfo<dim>::children_per_cell];
-
-                                    /**
-                                     * Array of embedding matrices. See
-                                     * <tt>get_prolongation_matrix()</tt> above.
-                                     *
-                                     * Matrices in this array are
-                                     * automatically initialized to
-                                     * correct size. If the derived
-                                     * finite element class does not
-                                     * implement these matrices, they
-                                     * should be resized to zero
-                                     * size.
-                                     */
-    FullMatrix<double> prolongation[GeometryInfo<dim>::children_per_cell];
-
-                                    /**
-                                     * Specify the constraints which
-                                     * the dofs on the two sides of a
-                                     * cell interface underly if the
-                                     * line connects two cells of
-                                     * which one is refined once.
-                                     *
-                                     * For further details see the
-                                     * general description of the
-                                     * derived class.
-                                     *
-                                     * This field is obviously
-                                     * useless in one space dimension
-                                     * and has there a zero size.
-                                     */
-    FullMatrix<double> interface_constraints;
-
-                                     /**
-                                      * Return the size of interface
-                                      * constraint matrices. Since
-                                      * this is needed in every
-                                      * derived finite element class
-                                      * when initializing their size,
-                                      * it is placed into this
-                                      * function, to avoid having to
-                                      * recompute the
-                                      * dimension-dependent size of
-                                      * these matrices each time.
-                                      *
-                                      * Note that some elements do not
-                                      * implement the interface
-                                      * constraints for certain
-                                      * polynomial degrees. In this
-                                      * case, this function still
-                                      * returns the size these
-                                      * matrices should have when
-                                      * implemented, but the actual
-                                      * matrices are empty.
-                                      */
-    TableIndices<2>
-    interface_constraints_size () const;
-    
-                                    /**
-                                     * Store what
-                                     * @p system_to_component_index
-                                     * will return.
-                                     */
-    std::vector< std::pair<unsigned int, unsigned int> > system_to_component_table;
-
-                                     /**
-                                     * Map between linear dofs and
-                                     * component dofs on face. This
-                                     * is filled with default values
-                                     * in the constructor, but
-                                     * derived classes will have to
-                                     * overwrite the information if
-                                     * necessary.
-                                     *
-                                     * By component, we mean the
-                                     * vector component, not the base
-                                     * element. The information thus
-                                     * makes only sense if a shape
-                                     * function is non-zero in only
-                                     * one component.
-                                     */
-    std::vector< std::pair<unsigned int, unsigned int> > face_system_to_component_table;
-
-                                    /**
-                                     * For each shape function, store
-                                     * to which base element and
-                                     * which instance of this base
-                                     * element (in case its
-                                     * multiplicity is greater than
-                                     * one) it belongs, and its index
-                                     * within this base element. If
-                                     * the element is not composed of
-                                     * others, then base and instance
-                                     * are always zero, and the index
-                                     * is equal to the number of the
-                                     * shape function. If the element
-                                     * is composed of single
-                                     * instances of other elements
-                                     * (i.e. all with multiplicity
-                                     * one) all of which are scalar,
-                                     * then base values and dof
-                                     * indices within this element
-                                     * are equal to the
-                                     * @p system_to_component_table. It
-                                     * differs only in case the
-                                     * element is composed of other
-                                     * elements and at least one of
-                                     * them is vector-valued itself.
-                                     *
-                                     * This array has valid values
-                                     * also in the case of
-                                     * vector-valued
-                                     * (i.e. non-primitive) shape
-                                     * functions, in contrast to the
-                                     * @p system_to_component_table.
-                                     */
-    std::vector<std::pair<std::pair<unsigned int,unsigned int>,unsigned int> >
-    system_to_base_table;
-
-                                    /**
-                                     * Likewise for the indices on
-                                     * faces.
-                                     */
-    std::vector<std::pair<std::pair<unsigned int,unsigned int>,unsigned int> >
-    face_system_to_base_table;
-    
-                                    /**
-                                     * The base element establishing
-                                     * a component.
-                                     *
-                                     * This table converts a
-                                     * component number to a pair
-                                     * consisting of the
-                                     * @p base_element number, and
-                                     * the component within this base
-                                     * element. While component
-                                     * information contains
-                                     * multiplicity of base elements,
-                                     * the result allows access to
-                                     * shape functions of the base
-                                     * element.
-                                     *
-                                     * This variable is set to the
-                                     * correct size by the
-                                     * constructor of this class, but
-                                     * needs to be initialized by
-                                     * derived classes, unless its
-                                     * size is one and the only entry
-                                     * is a zero, which is the case
-                                     * for scalar elements. In that
-                                     * case, the initialization by
-                                     * the base class is sufficient.
-                                     */
-    std::vector<std::pair<unsigned int, unsigned int> > component_to_base_table;
-    
-                                    /**
-                                     * Projection matrices are
-                                     * concatenated or summed up.
-                                     *
-                                     * This flags decides on how the
-                                     * projection matrices of the
-                                     * children of the same father
-                                     * are put together to one
-                                     * operator. The possible modes
-                                     * are concatenation and
-                                     * summation.
-                                     *
-                                     * If the projection is defined
-                                     * by an interpolation operator,
-                                     * the child matrices are
-                                     * concatenated, i.e. values
-                                     * belonging to the same node
-                                     * functional are identified and
-                                     * enter the interpolated value
-                                     * only once. In this case, the
-                                     * flag must be @p false.
-                                     *
-                                     * For projections with respect
-                                     * to scalar products, the child
-                                     * matrices must be summed up to
-                                     * build the complete matrix. The
-                                     * flag should be @p true.
-                                     *
-                                     * For examples of use of these
-                                     * flags, see the places in the
-                                     * library where it is queried.
-                                     * 
-                                     * There is one flag per shape
-                                     * function, indicating whether
-                                     * it belongs to the class of
-                                     * shape functions that are
-                                     * additive in the restriction or
-                                     * not.
-                                     *
-                                     * Note that in previous versions
-                                     * of the library, there was one
-                                     * flag per vector component of
-                                     * the element. This is based on
-                                     * the fact that all the shape
-                                     * functions that belong to the
-                                     * same vector component must
-                                     * necessarily behave in the same
-                                     * way, to make things
-                                     * reasonable. However, the
-                                     * problem is that it is
-                                     * sometimes impossible to query
-                                     * this flag in the vector-valued
-                                     * case: this used to be done
-                                     * with the
-                                     * @p system_to_component_index
-                                     * function that returns which
-                                     * vector component a shape
-                                     * function is associated
-                                     * with. The point is that since
-                                     * we now support shape functions
-                                     * that are associated with more
-                                     * than one vector component (for
-                                     * example the shape functions of
-                                     * Raviart-Thomas, or Nedelec
-                                     * elements), that function can
-                                     * no more be used, so it can be
-                                     * difficult to find out which
-                                     * for vector component we would
-                                     * like to query the
-                                     * restriction-is-additive flags.
-                                     */
-    const std::vector<bool> restriction_is_additive_flags;
-
-                                    /**
-                                     * List of support points on the
-                                     * unit cell, in case the finite
-                                     * element has any. The
-                                     * constructor leaves this field
-                                     * empty, derived classes may
-                                     * write in some contents.
-                                     *
-                                     * Finite elements that allow
-                                     * some kind of interpolation
-                                     * operation usually have support
-                                     * points. On the other hand,
-                                     * elements that define their
-                                     * degrees of freedom by, for
-                                     * example, moments on faces, or
-                                     * as derivatives, don't have
-                                     * support points. In that case,
-                                     * this field remains empty.
-                                     */
-    std::vector<Point<dim> > unit_support_points;
-
-                                    /**
-                                     * Same for the faces. See the
-                                     * description of the
-                                     * @p get_unit_face_support_points
-                                     * function for a discussion of
-                                     * what contributes a face
-                                     * support point.
-                                     */
-    std::vector<Point<dim-1> > unit_face_support_points;
-    
-                                    /**
-                                     * Support points used for
-                                     * interpolation functions of
-                                     * non-Lagrangian elements.
-                                     */
-    std::vector<Point<dim> > generalized_support_points;
-    
-                                    /**
-                                     * Face support points used for
-                                     * interpolation functions of
-                                     * non-Lagrangian elements.
-                                     */    
-    std::vector<Point<dim-1> > generalized_face_support_points;
-    
-                                    /**
-                                     * For each shape function, give
-                                     * a vector of bools (with size
-                                     * equal to the number of vector
-                                     * components which this finite
-                                     * element has) indicating in
-                                     * which component each of these
-                                     * shape functions is non-zero.
-                                     *
-                                     * For primitive elements, there
-                                     * is only one non-zero
-                                     * component.
-                                     */
-    const std::vector<std::vector<bool> > nonzero_components;
-
-                                    /**
-                                     * This array holds how many
-                                     * values in the respective entry
-                                     * of the @p nonzero_components
-                                     * element are non-zero. The
-                                     * array is thus a short-cut to
-                                     * allow faster access to this
-                                     * information than if we had to
-                                     * count the non-zero entries
-                                     * upon each request for this
-                                     * information. The field is
-                                     * initialized in the constructor
-                                     * of this class.
-                                     */
-    const std::vector<unsigned int> n_nonzero_components_table;
-
-                                    /**
-                                     * Store whether all shape
-                                     * functions are primitive. Since
-                                     * finding this out is a very
-                                     * common operation, we cache the
-                                     * result, i.e. compute the value
-                                     * in the constructor for simpler
-                                     * access.
-                                     */
-    const bool cached_primitivity;
-
-                                     /**
-                                     * Compute second derivatives by
-                                     * finite differences of
-                                     * gradients.
-                                     */
-    void compute_2nd (const Mapping<dim>                      &mapping,
-                     const typename Triangulation<dim>::cell_iterator    &cell,
-                     const unsigned int                       offset,
-                     typename Mapping<dim>::InternalDataBase &mapping_internal,
-                     InternalDataBase                        &fe_internal,
-                     FEValuesData<dim>                       &data) const;
-
-  private:
-                                    /**
-                                     * Second derivatives of shapes
-                                     * functions are not computed
-                                     * analytically, but by finite
-                                     * differences of the
-                                     * gradients. This static
-                                     * variable denotes the step
-                                     * length to be used for
-                                     * that. It's value is set to
-                                     * 1e-6.
-                                     */
-    static const double fd_step_length;
-
-                                    /**
-                                     * Given the pattern of nonzero
-                                     * components for each shape
-                                     * function, compute for each
-                                     * entry how many components are
-                                     * non-zero for each shape
-                                     * function. This function is
-                                     * used in the constructor of
-                                     * this class.
-                                     */
-    static
-    std::vector<unsigned int>
-    compute_n_nonzero_components (const std::vector<std::vector<bool> > &nonzero_components);
-    
-                                    /**
-                                     * Allow the FESystem class to access the
-                                     * restriction and prolongation matrices
-                                     * directly. Hence, FESystem has the
-                                     * possibility to see if these matrices
-                                     * are initialized without accessing
-                                     * these matrices through the
-                                     * @p get_restriction_matrix and
-                                     * @p get_prolongation_matrix
-                                     * functions. This is important as these
-                                     * functions include assertions that
-                                     * throw if the matrices are not already
-                                     * initialized.
-                                     */
-    template <int dim_> friend class FESystem;
-
-                                     /**
-                                      * Make the inner class a
-                                      * friend. This is not strictly
-                                      * necessary, but the Intel
-                                      * compiler seems to want this.
-                                      */
-    friend class InternalDataBase;
-};
 /*@}*/
 
 
@@ -2186,208 +646,6 @@ FiniteElementData<dim>::conforms (Conformity space) const
   return ((space & conforming_space) != 0);
 }
 
-//----------------------------------------------------------------------//
-
-template <int dim>  
-inline
-std::pair<unsigned int,unsigned int>
-FiniteElementBase<dim>::system_to_component_index (const unsigned int index) const
-{
-  Assert (index < system_to_component_table.size(),
-        ExcIndexRange(index, 0, system_to_component_table.size()));
-  Assert (is_primitive (index),
-         typename FiniteElementBase<dim>::ExcShapeFunctionNotPrimitive(index));
-  return system_to_component_table[index];
-}
-
-template <int dim>  
-inline
-unsigned int
-FiniteElementBase<dim>::component_to_system_index (const unsigned int component,
-                                                   const unsigned int index) const
-{
-   std::vector< std::pair<unsigned int, unsigned int> >::const_iterator
-      it = std::find(system_to_component_table.begin(), system_to_component_table.end(),
-                     std::pair<unsigned int, unsigned int>(component, index));
-
-   Assert(it != system_to_component_table.end(), ExcComponentIndexInvalid(component, index));
-   return std::distance(system_to_component_table.begin(), it);
-}
-
-
-
-template <int dim>  
-inline
-std::pair<unsigned int,unsigned int>
-FiniteElementBase<dim>::face_system_to_component_index (const unsigned int index) const
-{
-  Assert(index < face_system_to_component_table.size(),
-        ExcIndexRange(index, 0, face_system_to_component_table.size()));
-
-                                   // in debug mode, check whether the
-                                   // function is primitive, since
-                                   // otherwise the result may have no
-                                   // meaning
-                                   //
-                                   // since the primitivity tables are
-                                   // all geared towards cell dof
-                                   // indices, rather than face dof
-                                   // indices, we have to work a
-                                   // little bit...
-                                   //
-                                   // in 1d, the face index is equal
-                                   // to the cell index
-  Assert (((dim == 1) && is_primitive(index))
-          ||
-                                           // in 2d, construct it like
-                                           // this:
-          ((dim == 2) &&
-           is_primitive (index < (GeometryInfo<2>::vertices_per_face *
-                                  this->dofs_per_vertex)
-                         ?
-                         index
-                         :
-                         GeometryInfo<2>::vertices_per_cell *
-                         this->dofs_per_vertex +
-                         (index -
-                          GeometryInfo<2>::vertices_per_face *
-                          this->dofs_per_vertex)))
-          ||
-                                           // likewise in 3d, but more
-                                           // complicated
-          ((dim == 3) &&
-           is_primitive (index < (GeometryInfo<3>::vertices_per_face *
-                                  this->dofs_per_vertex)
-                         ?
-                         index
-                         :
-                         (index < (GeometryInfo<3>::vertices_per_face *
-                                   this->dofs_per_vertex
-                                   +
-                                   GeometryInfo<3>::lines_per_face *
-                                   this->dofs_per_line)
-                          ?
-                          GeometryInfo<3>::vertices_per_cell *
-                          this->dofs_per_vertex +
-                          (index -
-                           GeometryInfo<3>::vertices_per_face *
-                           this->dofs_per_vertex)
-                          :
-                          GeometryInfo<3>::vertices_per_cell *
-                          this->dofs_per_vertex +
-                          GeometryInfo<3>::lines_per_cell *
-                          this->dofs_per_line +
-                          (index -
-                           GeometryInfo<3>::vertices_per_face *
-                           this->dofs_per_vertex
-                           -
-                           GeometryInfo<3>::lines_per_face *
-                           this->dofs_per_line)))),
-          typename FiniteElementBase<dim>::ExcShapeFunctionNotPrimitive(index));
-
-  return face_system_to_component_table[index];
-}
-
-
-
-template <int dim>  
-inline
-std::pair<std::pair<unsigned int,unsigned int>,unsigned int>
-FiniteElementBase<dim>::system_to_base_index (const unsigned int index) const
-{
-  Assert (index < system_to_base_table.size(),
-        ExcIndexRange(index, 0, system_to_base_table.size()));
-  return system_to_base_table[index];
-}
-
-
-
-
-template <int dim>  
-inline
-std::pair<std::pair<unsigned int,unsigned int>,unsigned int>
-FiniteElementBase<dim>::face_system_to_base_index (const unsigned int index) const
-{
-  Assert(index < face_system_to_base_table.size(),
-        ExcIndexRange(index, 0, face_system_to_base_table.size()));
-  return face_system_to_base_table[index];
-}
-
-
-
-template <int dim>  
-inline
-std::pair<unsigned int,unsigned int>
-FiniteElementBase<dim>::component_to_base (const unsigned int index) const
-{
-  Assert(index < component_to_base_table.size(),
-        ExcIndexRange(index, 0, component_to_base_table.size()));
-
-  return component_to_base_table[index];
-}
-
-
-template <int dim>
-inline
-bool
-FiniteElementBase<dim>::restriction_is_additive (const unsigned int index) const
-{
-  Assert(index < this->dofs_per_cell,
-        ExcIndexRange(index, 0, this->dofs_per_cell));
-  return restriction_is_additive_flags[index];
-}
-
-
-template <int dim>
-inline
-const std::vector<bool> &
-FiniteElementBase<dim>::get_nonzero_components (const unsigned int i) const
-{
-  Assert (i < this->dofs_per_cell, ExcIndexRange (i, 0, this->dofs_per_cell));
-  return nonzero_components[i];
-}
-
-
-
-template <int dim>
-inline
-unsigned int
-FiniteElementBase<dim>::n_nonzero_components (const unsigned int i) const
-{
-  Assert (i < this->dofs_per_cell, ExcIndexRange (i, 0, this->dofs_per_cell));
-  return n_nonzero_components_table[i];
-}
-
-
-
-template <int dim>
-inline
-bool
-FiniteElementBase<dim>::is_primitive (const unsigned int i) const
-{
-  Assert (i < this->dofs_per_cell, ExcIndexRange (i, 0, this->dofs_per_cell));
-
-                                  // return primitivity of a shape
-                                  // function by checking whether it
-                                  // has more than one non-zero
-                                  // component or not. we could cache
-                                  // this value in an array of bools,
-                                  // but accessing a bit-vector (as
-                                  // std::vector<bool> is) is
-                                  // probably more expensive than
-                                  // just comparing against 1
-  return (n_nonzero_components_table[i] == 1);
-}
-
-
-template <int dim>
-inline
-bool
-FiniteElementBase<dim>::is_primitive () const
-{
-  return cached_primitivity;
-}
-
 
 
 
index 2409122a6e9174a9854e2fe132727b43d0b6deee..5698ce39fb6cc679f76632499beb0aee55ae9a0e 100644 (file)
@@ -68,11 +68,11 @@ class FE_DGPMonomial : public FE_Poly<PolynomialsP<dim>,dim>
                                      * element is also a @p FE_Q
                                      * element. Otherwise, an
                                      * exception of type
-                                     * FiniteElementBase<dim>::ExcInterpolationNotImplemented
+                                     * FiniteElement<dim>::ExcInterpolationNotImplemented
                                      * is thrown.
                                      */
     virtual void
-    get_interpolation_matrix (const FiniteElementBase<dim> &source,
+    get_interpolation_matrix (const FiniteElement<dim> &source,
                              FullMatrix<double>           &matrix) const;
     
                                     /**
index 6f6274308a86a05a04d2a55b5c3e3bbc1e916ee8..0b7fab840051b4cc6dfd710b8912bc11ddbc7c85 100644 (file)
@@ -70,7 +70,7 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      * Return the value of the
                                      * @p ith shape function at the
                                      * point @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -84,7 +84,7 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      * component of the @p ith shape
                                      * function at the point
                                      * @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -104,7 +104,7 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      * Return the gradient of the
                                      * @p ith shape function at the
                                      * point @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -118,7 +118,7 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      * component of the @p ith shape
                                      * function at the point
                                      * @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -139,7 +139,7 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      * derivatives of the @p ith
                                      * shape function at point @p p
                                      * on the unit cell.  See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -153,7 +153,7 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      * component of the @p ith shape
                                      * function at the point
                                      * @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -418,7 +418,7 @@ class FE_DGPNonparametric : public FiniteElement<dim>
                                      * see the documentation of the
                                      * base class.
                                      */
-    class InternalData : public FiniteElementBase<dim>::InternalDataBase
+    class InternalData : public FiniteElement<dim>::InternalDataBase
     {
       public:
                                       // have some scratch arrays
index 79804159159036bcbf4b0efd5870abdde6ebfe71..27fff613f09177c36b91b7125566ebe87340397e 100644 (file)
@@ -115,11 +115,11 @@ class FE_DGQ : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                      * element is also a @p FE_DGQ
                                      * element. Otherwise, an
                                      * exception of type
-                                     * FiniteElementBase<dim>::ExcInterpolationNotImplemented
+                                     * FiniteElement<dim>::ExcInterpolationNotImplemented
                                      * is thrown.
                                      */
     virtual void
-    get_interpolation_matrix (const FiniteElementBase<dim> &source,
+    get_interpolation_matrix (const FiniteElement<dim> &source,
                              FullMatrix<double>           &matrix) const;
     
                                     /**
index f5b06b293527b8b198ed0b3975deab3e3fd5415c..e4f8b12f20eb44bd011c361a5b3cd8d924106760 100644 (file)
@@ -205,7 +205,7 @@ class FE_Nedelec : public FiniteElement<dim>
                                      * component of the @p ith shape
                                      * function at the point
                                      * @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -220,7 +220,7 @@ class FE_Nedelec : public FiniteElement<dim>
                                      * component of the @p ith shape
                                      * function at the point
                                      * @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -235,7 +235,7 @@ class FE_Nedelec : public FiniteElement<dim>
                                      * component of the @p ith shape
                                      * function at the point
                                      * @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -482,7 +482,7 @@ class FE_Nedelec : public FiniteElement<dim>
                                     /**
                                      * Initialize the
                                      * @p unit_support_points field
-                                     * of the FiniteElementBase
+                                     * of the FiniteElement
                                      * class. Called from the
                                      * constructor.
                                      */
@@ -491,7 +491,7 @@ class FE_Nedelec : public FiniteElement<dim>
                                     /**
                                      * Initialize the
                                      * @p unit_face_support_points field
-                                     * of the FiniteElementBase
+                                     * of the FiniteElement
                                      * class. Called from the
                                      * constructor.
                                      */
@@ -549,7 +549,7 @@ class FE_Nedelec : public FiniteElement<dim>
                                      * see the documentation of the
                                      * base class.
                                      */
-    class InternalData : public FiniteElementBase<dim>::InternalDataBase
+    class InternalData : public FiniteElement<dim>::InternalDataBase
     {
       public:
                                         /**
index 4b15ec7b2fa98d4d6e99b5a8a0bf0bf3256e34b9..bffcdc241dbbeafa71477b1925905dbc8b700539 100644 (file)
@@ -42,7 +42,7 @@
  *
  * This class is not a fully implemented FiniteElement class. Instead
  * there are several pure virtual functions declared in the
- * FiniteElement and FiniteElementBase classes which cannot
+ * FiniteElement and FiniteElement classes which cannot
  * implemented by this class but are left for implementation in
  * derived classes.
  *
@@ -83,7 +83,7 @@ class FE_Poly : public FiniteElement<dim>
                                      * Return the value of the
                                      * <tt>i</tt>th shape function at
                                      * the point <tt>p</tt>. See the
-                                     * FiniteElementBase base class
+                                     * FiniteElement base class
                                      * for more information about the
                                      * semantics of this function.
                                      */
@@ -96,7 +96,7 @@ class FE_Poly : public FiniteElement<dim>
                                      * component of the <tt>i</tt>th
                                      * shape function at the point
                                      * <tt>p</tt>. See the
-                                     * FiniteElementBase base class
+                                     * FiniteElement base class
                                      * for more information about the
                                      * semantics of this function.
                                      *
@@ -115,7 +115,7 @@ class FE_Poly : public FiniteElement<dim>
                                      * Return the gradient of the
                                      * <tt>i</tt>th shape function at
                                      * the point <tt>p</tt>. See the
-                                     * FiniteElementBase base class
+                                     * FiniteElement base class
                                      * for more information about the
                                      * semantics of this function.
                                      */
@@ -128,7 +128,7 @@ class FE_Poly : public FiniteElement<dim>
                                      * component of the <tt>i</tt>th
                                      * shape function at the point
                                      * <tt>p</tt>. See the
-                                     * FiniteElementBase base class
+                                     * FiniteElement base class
                                      * for more information about the
                                      * semantics of this function.
                                      *
@@ -149,7 +149,7 @@ class FE_Poly : public FiniteElement<dim>
                                      * <tt>i</tt>th shape function at
                                      * point <tt>p</tt> on the unit
                                      * cell. See the
-                                     * FiniteElementBase base class
+                                     * FiniteElement base class
                                      * for more information about the
                                      * semantics of this function.
                                      */
@@ -162,7 +162,7 @@ class FE_Poly : public FiniteElement<dim>
                                      * vector component of the
                                      * <tt>i</tt>th shape function at
                                      * the point <tt>p</tt>. See the
-                                     * FiniteElementBase base class
+                                     * FiniteElement base class
                                      * for more information about the
                                      * semantics of this function.
                                      *
@@ -336,7 +336,7 @@ class FE_Poly : public FiniteElement<dim>
                                      * see the documentation of the
                                      * base class.
                                      */
-    class InternalData : public FiniteElementBase<dim>::InternalDataBase
+    class InternalData : public FiniteElement<dim>::InternalDataBase
     {
       public:
                                         /**
index 7378a58d638076b3fc01635bebc39c1c58c15f2b..c0cf51b345562b4754a07eedcca579c282d7fe9a 100644 (file)
@@ -318,7 +318,7 @@ class FE_PolyTensor : public FiniteElement<dim>
                                      * accessed by indices
                                      * <i>(i,k)</i>.
                                      */
-    class InternalData : public FiniteElementBase<dim>::InternalDataBase
+    class InternalData : public FiniteElement<dim>::InternalDataBase
     {
       public:
                                         /**
index 8ca83af7810c0cad38ee6641a1a655eecef92003..21588ebe374f3fbf4c1c04a240380eeeec5baf88 100644 (file)
@@ -269,11 +269,11 @@ class FE_Q : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                      * element is also a @p FE_Q
                                      * element. Otherwise, an
                                      * exception of type
-                                     * FiniteElementBase<dim>::ExcInterpolationNotImplemented
+                                     * FiniteElement<dim>::ExcInterpolationNotImplemented
                                      * is thrown.
                                      */
     virtual void
-    get_interpolation_matrix (const FiniteElementBase<dim> &source,
+    get_interpolation_matrix (const FiniteElement<dim> &source,
                              FullMatrix<double>           &matrix) const;
     
                                     /**
@@ -364,7 +364,7 @@ class FE_Q : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                     /**
                                      * Initialize the
                                      * @p unit_support_points field
-                                     * of the FiniteElementBase
+                                     * of the FiniteElement
                                      * class. Called from the
                                      * constructor.
                                      */
@@ -373,7 +373,7 @@ class FE_Q : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                     /**
                                      * Initialize the
                                      * @p unit_face_support_points field
-                                     * of the FiniteElementBase
+                                     * of the FiniteElement
                                      * class. Called from the
                                      * constructor.
                                      */
index 3ba126f3ad1f8757d277a23a4dfa5a243c197d36..1c8aa5e7689821e302380b7ed1980233bd7f2678 100644 (file)
@@ -409,7 +409,7 @@ class FE_Q_Hierarchical : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                     /**
                                      * Initialize the
                                      * @p unit_support_points field
-                                     * of the FiniteElementBase
+                                     * of the FiniteElement
                                      * class. Called from the
                                      * constructor.
                                      */
@@ -418,7 +418,7 @@ class FE_Q_Hierarchical : public FE_Poly<TensorProductPolynomials<dim>,dim>
                                     /**
                                      * Initialize the
                                      * @p unit_face_support_points field
-                                     * of the FiniteElementBase
+                                     * of the FiniteElement
                                      * class. Called from the
                                      * constructor.
                                      */
index 0cea7b8eb62539778451765f8a1d4907a1a82bea..2f7f3f07a01a131aa767537dcbafdef3aaf03c62 100644 (file)
@@ -163,7 +163,7 @@ class FE_RaviartThomas : public FiniteElement<dim>
                                      * component of the @p ith shape
                                      * function at the point
                                      * @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -178,7 +178,7 @@ class FE_RaviartThomas : public FiniteElement<dim>
                                      * component of the @p ith shape
                                      * function at the point
                                      * @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -193,7 +193,7 @@ class FE_RaviartThomas : public FiniteElement<dim>
                                      * component of the @p ith shape
                                      * function at the point
                                      * @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -237,11 +237,11 @@ class FE_RaviartThomas : public FiniteElement<dim>
                                      * element is also a Raviart
                                      * Thomas element. Otherwise, an
                                      * exception of type
-                                     * FiniteElementBase<dim>::ExcInterpolationNotImplemented
+                                     * FiniteElement<dim>::ExcInterpolationNotImplemented
                                      * is thrown.
                                      */
     virtual void
-    get_interpolation_matrix (const FiniteElementBase<dim> &source,
+    get_interpolation_matrix (const FiniteElement<dim> &source,
                              FullMatrix<double>           &matrix) const;
 
                                     /**
@@ -439,7 +439,7 @@ class FE_RaviartThomas : public FiniteElement<dim>
                                     /**
                                      * Initialize the
                                      * @p unit_support_points field
-                                     * of the FiniteElementBase
+                                     * of the FiniteElement
                                      * class. Called from the
                                      * constructor.
                                      */
@@ -448,7 +448,7 @@ class FE_RaviartThomas : public FiniteElement<dim>
                                     /**
                                      * Initialize the
                                      * @p unit_face_support_points field
-                                     * of the FiniteElementBase
+                                     * of the FiniteElement
                                      * class. Called from the
                                      * constructor.
                                      */
@@ -501,7 +501,7 @@ class FE_RaviartThomas : public FiniteElement<dim>
                                      * see the documentation of the
                                      * base class.
                                      */
-    class InternalData : public FiniteElementBase<dim>::InternalDataBase
+    class InternalData : public FiniteElement<dim>::InternalDataBase
     {
       public:
                                         /**
@@ -655,8 +655,8 @@ class FE_RaviartThomasNodal
     get_ria_vector (const unsigned int degree);
                                     /**
                                      * Initialize the
-                                     * FiniteElementBase<dim>::unit_support_points
-                                     * and FiniteElementBase<dim>::unit_face_support_points
+                                     * FiniteElement<dim>::unit_support_points
+                                     * and FiniteElement<dim>::unit_face_support_points
                                      * fields. Called from the
                                      * constructor.
                                      */
@@ -739,7 +739,7 @@ FE_RaviartThomas<3>::initialize_restriction ();
 template <>
 void
 FE_RaviartThomas<1>::
-get_interpolation_matrix (const FiniteElementBase<1> &,
+get_interpolation_matrix (const FiniteElement<1> &,
                          FullMatrix<double>         &) const;
 
 /// @endif
index 1928683a2a7214c8d23ee3ca80eb571a16490102..64f34b053991bac6385a935afa0cfe6088e04e4f 100644 (file)
@@ -174,7 +174,7 @@ class FESystem : public FiniteElement<dim>
                                      * component of the @p ith shape
                                      * function at the point
                                      * @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -227,7 +227,7 @@ class FESystem : public FiniteElement<dim>
                                      * component of the @p ith shape
                                      * function at the point
                                      * @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -279,7 +279,7 @@ class FESystem : public FiniteElement<dim>
                                      * component of the @p ith shape
                                      * function at the point
                                      * @p p. See the
-                                     * FiniteElementBase base
+                                     * FiniteElement base
                                      * class for more information
                                      * about the semantics of this
                                      * function.
@@ -314,11 +314,11 @@ class FESystem : public FiniteElement<dim>
                                      * @p get_interpolation_matrix
                                      * functions. Otherwise, an
                                      * exception of type
-                                     * FiniteElementBase<dim>::ExcInterpolationNotImplemented
+                                     * FiniteElement<dim>::ExcInterpolationNotImplemented
                                      * is thrown.
                                      */
     virtual void
-    get_interpolation_matrix (const FiniteElementBase<dim> &source,
+    get_interpolation_matrix (const FiniteElement<dim> &source,
                              FullMatrix<double>           &matrix) const;
 
                                      /** 
@@ -558,7 +558,7 @@ class FESystem : public FiniteElement<dim>
                                     /**
                                      * Initialize the
                                      * @p unit_support_points field
-                                     * of the FiniteElementBase
+                                     * of the FiniteElement
                                      * class. Called from the
                                      * constructor.
                                      */
@@ -567,7 +567,7 @@ class FESystem : public FiniteElement<dim>
                                     /**
                                      * Initialize the
                                      * @p unit_face_support_points field
-                                     * of the FiniteElementBase
+                                     * of the FiniteElement
                                      * class. Called from the
                                      * constructor.
                                      */
@@ -743,7 +743,7 @@ class FESystem : public FiniteElement<dim>
                                      * @p InternalData objects for
                                      * each of the base elements.
                                      */
-    class InternalData : public FiniteElementBase<dim>::InternalDataBase
+    class InternalData : public FiniteElement<dim>::InternalDataBase
     {
       public:
                                         /**
@@ -780,7 +780,7 @@ class FESystem : public FiniteElement<dim>
                                          * element.
                                          */
        void set_fe_data(const unsigned int                        base_no,
-                        typename FiniteElementBase<dim>::InternalDataBase *);
+                        typename FiniteElement<dim>::InternalDataBase *);
 
                                         /**
                                          * Gives read-access to the
@@ -788,7 +788,7 @@ class FESystem : public FiniteElement<dim>
                                          * @p InternalData of the
                                          * @p base_noth base element.
                                          */    
-       typename FiniteElementBase<dim>::InternalDataBase &
+       typename FiniteElement<dim>::InternalDataBase &
        get_fe_data (const unsigned int base_no) const;
 
 
@@ -872,7 +872,7 @@ class FESystem : public FiniteElement<dim>
                                          * elements, irrespective of
                                          * their multiplicity.
                                          */
-       typename std::vector<typename FiniteElementBase<dim>::InternalDataBase *> base_fe_datas;
+       typename std::vector<typename FiniteElement<dim>::InternalDataBase *> base_fe_datas;
 
                                         /**
                                          * Pointers to the
index 5bdb67003e47d300fb29fbf25bd400791fdad03e..5c88bc49fc136ae1cc1007e740c4cf04cc7158d3 100644 (file)
@@ -195,7 +195,7 @@ class FETools
                                      * @param matrices A pointer to
                                      * <i>2<sup>dim</sup></i> FullMatrix
                                      * objects. This is the format
-                                     * used in FiniteElementBase,
+                                     * used in FiniteElement,
                                      * where we want to use ths
                                      * function mostly.
                                      */
@@ -241,7 +241,7 @@ class FETools
                                      * matrices.  @arg matrices A pointer to
                                      * <tt>GeometryInfo::children_per_cell</tt>=2<sup>dim</sup>
                                      * FullMatrix objects. This is the format
-                                     * used in FiniteElementBase, where we
+                                     * used in FiniteElement, where we
                                      * want to use this function mostly.
                                      */
     template <int dim, typename number>
index 3e8a53a2f60ea937461b1faf3002ffc7b86d43b8..bc82fc01c4ad8e9315a74a44753df8458f37770f 100644 (file)
@@ -92,7 +92,7 @@ DoFTools::compute_row_length_vector(
   const Table<2,Coupling>& flux_couplings)
 {
   const FiniteElement<dim>& fe = dofs.get_fe();
-  Assert (fe.is_primitive(), typename FiniteElementBase<dim>::ExcFENotPrimitive());
+  Assert (fe.is_primitive(), typename FiniteElement<dim>::ExcFENotPrimitive());
   Assert (row_lengths.size() == dofs.n_dofs(),
          ExcDimensionMismatch(row_lengths.size(), dofs.n_dofs()));
   Assert (couplings.n_rows()==fe.n_components(),
@@ -329,7 +329,7 @@ DoFTools::compute_row_length_vector(
   const Table<2,Coupling>& flux_couplings)
 {
   const FiniteElement<dim>& fe = dofs.get_fe();
-  Assert (fe.is_primitive(), typename FiniteElementBase<dim>::ExcFENotPrimitive());
+  Assert (fe.is_primitive(), typename FiniteElement<dim>::ExcFENotPrimitive());
   Assert (row_lengths.size() == dofs.n_dofs(),
          ExcDimensionMismatch(row_lengths.size(), dofs.n_dofs()));
   Assert (couplings.n_rows()==fe.n_components(),
index 1463402d32e22107ae02899c12f2820f6c1204ce..6747e018324e55ba612af4379659648e02e17be7 100644 (file)
 
 
 
-/*------------------------------- FiniteElementBase ----------------------*/
+/*------------------------------- FiniteElement ----------------------*/
 
 
 template <int dim>
-const double FiniteElementBase<dim>::fd_step_length = 1.0e-6;
+const double FiniteElement<dim>::fd_step_length = 1.0e-6;
 
 
 template <int dim>
 void
-FiniteElementBase<dim>::
+FiniteElement<dim>::
 InternalDataBase::initialize_2nd (const FiniteElement<dim> *element,
                                  const Mapping<dim>       &mapping,
                                  const Quadrature<dim>    &quadrature)
@@ -89,7 +89,7 @@ InternalDataBase::initialize_2nd (const FiniteElement<dim> *element,
 
 
 template <int dim>
-FiniteElementBase<dim>::InternalDataBase::~InternalDataBase ()
+FiniteElement<dim>::InternalDataBase::~InternalDataBase ()
 {
   for (unsigned int i=0; i<differences.size (); ++i)
     if (differences[i] != 0)
@@ -106,7 +106,7 @@ FiniteElementBase<dim>::InternalDataBase::~InternalDataBase ()
 
 
 template <int dim>
-FiniteElementBase<dim>::FiniteElementBase (
+FiniteElement<dim>::FiniteElement (
   const FiniteElementData<dim> &fe_data,
   const std::vector<bool> &r_i_a_f,
   const std::vector<std::vector<bool> > &nonzero_c)
@@ -199,9 +199,31 @@ FiniteElementBase<dim>::FiniteElementBase (
 }
 
 
+template <int dim>
+FiniteElement<dim>::FiniteElement (const FiniteElement<dim> &)
+               :
+               Subscriptor(),
+               FiniteElementData<dim>()
+{
+  Assert (false,
+          ExcMessage ("Finite element objects don't support copying "
+                      "semantics through the copy constructor. If "
+                      "you want to copy a finite element, use the "
+                      "clone() function."));
+}
+
+
+
+template <int dim>
+FiniteElement<dim>::~FiniteElement ()
+{}
+
+
+
+
 template <int dim>
 double
-FiniteElementBase<dim>::shape_value (const unsigned int,
+FiniteElement<dim>::shape_value (const unsigned int,
                                     const Point<dim> &) const
 {
   AssertThrow(false, ExcUnitShapeValuesDoNotExist());
@@ -212,7 +234,7 @@ FiniteElementBase<dim>::shape_value (const unsigned int,
 
 template <int dim>
 double
-FiniteElementBase<dim>::shape_value_component (const unsigned int,
+FiniteElement<dim>::shape_value_component (const unsigned int,
                                               const Point<dim> &,
                                               const unsigned int) const
 {
@@ -224,7 +246,7 @@ FiniteElementBase<dim>::shape_value_component (const unsigned int,
 
 template <int dim>
 Tensor<1,dim>
-FiniteElementBase<dim>::shape_grad (const unsigned int,
+FiniteElement<dim>::shape_grad (const unsigned int,
                                    const Point<dim> &) const
 {
   AssertThrow(false, ExcUnitShapeValuesDoNotExist());
@@ -235,7 +257,7 @@ FiniteElementBase<dim>::shape_grad (const unsigned int,
 
 template <int dim>
 Tensor<1,dim>
-FiniteElementBase<dim>::shape_grad_component (const unsigned int,
+FiniteElement<dim>::shape_grad_component (const unsigned int,
                                              const Point<dim> &,
                                              const unsigned int) const
 {
@@ -247,7 +269,7 @@ FiniteElementBase<dim>::shape_grad_component (const unsigned int,
 
 template <int dim>
 Tensor<2,dim>
-FiniteElementBase<dim>::shape_grad_grad (const unsigned int,
+FiniteElement<dim>::shape_grad_grad (const unsigned int,
                                         const Point<dim> &) const
 {
   AssertThrow(false, ExcUnitShapeValuesDoNotExist());
@@ -258,7 +280,7 @@ FiniteElementBase<dim>::shape_grad_grad (const unsigned int,
 
 template <int dim>
 Tensor<2,dim>
-FiniteElementBase<dim>::shape_grad_grad_component (const unsigned int,
+FiniteElement<dim>::shape_grad_grad_component (const unsigned int,
                                                   const Point<dim> &,
                                                   const unsigned int) const
 {
@@ -269,7 +291,7 @@ FiniteElementBase<dim>::shape_grad_grad_component (const unsigned int,
 
 template <int dim>
 const FullMatrix<double> &
-FiniteElementBase<dim>::get_restriction_matrix (const unsigned int child) const
+FiniteElement<dim>::get_restriction_matrix (const unsigned int child) const
 {
   Assert (child<GeometryInfo<dim>::children_per_cell,
          ExcIndexRange(child, 0, GeometryInfo<dim>::children_per_cell));
@@ -281,7 +303,7 @@ FiniteElementBase<dim>::get_restriction_matrix (const unsigned int child) const
 
 template <int dim>
 const FullMatrix<double> &
-FiniteElementBase<dim>::get_prolongation_matrix (const unsigned int child) const
+FiniteElement<dim>::get_prolongation_matrix (const unsigned int child) const
 {
   Assert (child<GeometryInfo<dim>::children_per_cell,
          ExcIndexRange(child, 0, GeometryInfo<dim>::children_per_cell));
@@ -292,7 +314,7 @@ FiniteElementBase<dim>::get_prolongation_matrix (const unsigned int child) const
 
 template <int dim>
 bool
-FiniteElementBase<dim>::prolongation_is_implemented () const
+FiniteElement<dim>::prolongation_is_implemented () const
 {
   for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
     {
@@ -313,7 +335,7 @@ FiniteElementBase<dim>::prolongation_is_implemented () const
 
 template <int dim>
 bool
-FiniteElementBase<dim>::restriction_is_implemented () const
+FiniteElement<dim>::restriction_is_implemented () const
 {
   for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
     {
@@ -334,7 +356,7 @@ FiniteElementBase<dim>::restriction_is_implemented () const
 
 template <int dim>
 bool
-FiniteElementBase<dim>::constraints_are_implemented () const
+FiniteElement<dim>::constraints_are_implemented () const
 {
   return (this->dofs_per_face  == 0) || (interface_constraints.m() != 0);
 }
@@ -343,7 +365,7 @@ FiniteElementBase<dim>::constraints_are_implemented () const
 
 template <int dim>
 const FullMatrix<double> &
-FiniteElementBase<dim>::constraints () const
+FiniteElement<dim>::constraints () const
 {
   Assert ((this->dofs_per_face  == 0) || (interface_constraints.m() != 0),
           ExcConstraintsVoid());
@@ -360,7 +382,7 @@ FiniteElementBase<dim>::constraints () const
 
 template <int dim>
 TableIndices<2>
-FiniteElementBase<dim>::interface_constraints_size () const 
+FiniteElement<dim>::interface_constraints_size () const 
 {
   switch (dim)
     {
@@ -386,15 +408,15 @@ FiniteElementBase<dim>::interface_constraints_size () const
 
 template <int dim>
 void
-FiniteElementBase<dim>::
-get_interpolation_matrix (const FiniteElementBase<dim> &,
+FiniteElement<dim>::
+get_interpolation_matrix (const FiniteElement<dim> &,
                          FullMatrix<double>           &) const
 {
                                   // by default, no interpolation
                                   // implemented. so throw exception,
                                   // as documentation says
   AssertThrow (false,
-               typename FiniteElementBase<dim>::
+               typename FiniteElement<dim>::
                ExcInterpolationNotImplemented());
 }
 
@@ -403,7 +425,7 @@ get_interpolation_matrix (const FiniteElementBase<dim> &,
 
 template <int dim>
 bool
-FiniteElementBase<dim>::operator == (const FiniteElementBase<dim> &f) const
+FiniteElement<dim>::operator == (const FiniteElement<dim> &f) const
 {
   return ((static_cast<const FiniteElementData<dim>&>(*this) ==
           static_cast<const FiniteElementData<dim>&>(f)) &&
@@ -414,7 +436,7 @@ FiniteElementBase<dim>::operator == (const FiniteElementBase<dim> &f) const
 
 template <int dim>
 const std::vector<Point<dim> > &
-FiniteElementBase<dim>::get_unit_support_points () const
+FiniteElement<dim>::get_unit_support_points () const
 {
                                   // a finite element may define
                                   // support points, but only if
@@ -430,7 +452,7 @@ FiniteElementBase<dim>::get_unit_support_points () const
 
 template <int dim>
 bool
-FiniteElementBase<dim>::has_support_points () const
+FiniteElement<dim>::has_support_points () const
 {
   return (unit_support_points.size() != 0);
 }
@@ -439,7 +461,7 @@ FiniteElementBase<dim>::has_support_points () const
 
 template <int dim>
 const std::vector<Point<dim> > &
-FiniteElementBase<dim>::get_generalized_support_points () const
+FiniteElement<dim>::get_generalized_support_points () const
 {
                                   // a finite element may define
                                   // support points, but only if
@@ -455,7 +477,7 @@ FiniteElementBase<dim>::get_generalized_support_points () const
 
 template <int dim>
 bool
-FiniteElementBase<dim>::has_generalized_support_points () const
+FiniteElement<dim>::has_generalized_support_points () const
 {
   return (generalized_support_points.size() != 0);
 }
@@ -464,7 +486,7 @@ FiniteElementBase<dim>::has_generalized_support_points () const
 
 template <int dim>
 Point<dim>
-FiniteElementBase<dim>::unit_support_point (const unsigned index) const
+FiniteElement<dim>::unit_support_point (const unsigned index) const
 {
   Assert (index < this->dofs_per_cell,
           ExcIndexRange (index, 0, this->dofs_per_cell));
@@ -477,7 +499,7 @@ FiniteElementBase<dim>::unit_support_point (const unsigned index) const
 
 template <int dim>
 const std::vector<Point<dim-1> > &
-FiniteElementBase<dim>::get_unit_face_support_points () const
+FiniteElement<dim>::get_unit_face_support_points () const
 {
                                   // a finite element may define
                                   // support points, but only if
@@ -493,7 +515,7 @@ FiniteElementBase<dim>::get_unit_face_support_points () const
 
 template <int dim>
 bool
-FiniteElementBase<dim>::has_face_support_points () const
+FiniteElement<dim>::has_face_support_points () const
 {
   return (unit_face_support_points.size() != 0);
 }
@@ -502,7 +524,7 @@ FiniteElementBase<dim>::has_face_support_points () const
 
 template <int dim>
 const std::vector<Point<dim-1> > &
-FiniteElementBase<dim>::get_generalized_face_support_points () const
+FiniteElement<dim>::get_generalized_face_support_points () const
 {
                                   // a finite element may define
                                   // support points, but only if
@@ -518,7 +540,7 @@ FiniteElementBase<dim>::get_generalized_face_support_points () const
 
 template <int dim>
 bool
-FiniteElementBase<dim>::has_generalized_face_support_points () const
+FiniteElement<dim>::has_generalized_face_support_points () const
 {
   return (generalized_face_support_points.size() != 0);
 }
@@ -527,7 +549,7 @@ FiniteElementBase<dim>::has_generalized_face_support_points () const
 
 template <int dim>
 Point<dim-1>
-FiniteElementBase<dim>::unit_face_support_point (const unsigned index) const
+FiniteElement<dim>::unit_face_support_point (const unsigned index) const
 {
   Assert (index < this->dofs_per_face,
           ExcIndexRange (index, 0, this->dofs_per_face));
@@ -540,7 +562,7 @@ FiniteElementBase<dim>::unit_face_support_point (const unsigned index) const
 
 template <int dim>
 void
-FiniteElementBase<dim>::interpolate(
+FiniteElement<dim>::interpolate(
   std::vector<double>&       local_dofs,
   const std::vector<double>& values) const
 {
@@ -562,7 +584,7 @@ FiniteElementBase<dim>::interpolate(
 
 template <int dim>
 void
-FiniteElementBase<dim>::interpolate(
+FiniteElement<dim>::interpolate(
   std::vector<double>&    local_dofs,
   const std::vector<Vector<double> >& values,
   unsigned int offset) const
@@ -588,7 +610,7 @@ FiniteElementBase<dim>::interpolate(
 
 template <int dim>
 void
-FiniteElementBase<dim>::interpolate(
+FiniteElement<dim>::interpolate(
   std::vector<double>& local_dofs,
   const VectorSlice<const std::vector<std::vector<double> > >& values) const
 {
@@ -613,7 +635,7 @@ FiniteElementBase<dim>::interpolate(
 
 template <int dim>
 unsigned int
-FiniteElementBase<dim>::memory_consumption () const
+FiniteElement<dim>::memory_consumption () const
 {
   return (sizeof(FiniteElementData<dim>) +
          MemoryConsumption::
@@ -637,7 +659,7 @@ FiniteElementBase<dim>::memory_consumption () const
 
 template <int dim>
 void
-FiniteElementBase<dim>::
+FiniteElement<dim>::
 compute_2nd (const Mapping<dim>                   &mapping,
             const typename Triangulation<dim>::cell_iterator &cell,
             const unsigned int,
@@ -780,7 +802,7 @@ compute_2nd (const Mapping<dim>                   &mapping,
 
 template <int dim>
 std::vector<unsigned int>
-FiniteElementBase<dim>::compute_n_nonzero_components (
+FiniteElement<dim>::compute_n_nonzero_components (
   const std::vector<std::vector<bool> > &nonzero_components)
 {
   std::vector<unsigned int> retval (nonzero_components.size());
@@ -795,40 +817,6 @@ FiniteElementBase<dim>::compute_n_nonzero_components (
 
 /*------------------------------- FiniteElement ----------------------*/
 
-template <int dim>
-FiniteElement<dim>::FiniteElement (const FiniteElementData<dim> &fe_data,
-                                  const std::vector<bool> &restriction_is_additive_flags,
-                                  const std::vector<std::vector<bool> > &nonzero_components)
-                :
-               FiniteElementBase<dim> (fe_data,
-                                       restriction_is_additive_flags,
-                                       nonzero_components)
-{}
-
-
-
-template <int dim>
-FiniteElement<dim>::FiniteElement (const FiniteElement<dim> &)
-                :
-                FiniteElementBase<dim> (FiniteElementData<dim>(),
-                                        std::vector<bool> (),
-                                        std::vector<std::vector<bool> >())
-{
-  Assert (false,
-          ExcMessage ("Finite element objects don't support copying "
-                      "semantics through the copy constructor. If "
-                      "you want to copy a finite element, use the "
-                      "clone() function."));
-}
-
-
-
-template <int dim>
-FiniteElement<dim>::~FiniteElement ()
-{}
-
-
-
 template <int dim>
 typename Mapping<dim>::InternalDataBase *
 FiniteElement<dim>::get_face_data (const UpdateFlags       flags,
@@ -853,16 +841,6 @@ FiniteElement<dim>::get_subface_data (const UpdateFlags        flags,
 
 
 
-template <int dim>
-unsigned int
-FiniteElement<dim>::memory_consumption () const
-{
-  return FiniteElementBase<dim>::memory_consumption ();
-}
-
-
-
-
 template <int dim>
 const FiniteElement<dim>&
 FiniteElement<dim>::base_element(const unsigned index) const
@@ -873,5 +851,4 @@ FiniteElement<dim>::base_element(const unsigned index) const
 
 /*------------------------------- Explicit Instantiations -------------*/
 
-template class FiniteElementBase<deal_II_dimension>;
 template class FiniteElement<deal_II_dimension>;
index 32f0ea11c14b798ff3dc3454303f1692a686af60..36746d0a8e5f01ca80919f2ba086009e95603b18 100644 (file)
@@ -179,7 +179,7 @@ FE_DGPMonomial<dim>::clone() const
 template <int dim>
 void
 FE_DGPMonomial<dim>::
-get_interpolation_matrix (const FiniteElementBase<dim> &source_fe,
+get_interpolation_matrix (const FiniteElement<dim> &source_fe,
                          FullMatrix<double>           &interpolation_matrix) const
 {  
   const FE_DGPMonomial<dim> *source_dgp_monomial
index 81dfb6658a6ba0fef5557f7a4fea707f2fe024d9..02bf00a1c2bb006127c45df4b43e83bdd245d2bc 100644 (file)
@@ -361,7 +361,7 @@ FE_DGQ<dim>::rotate_indices (std::vector<unsigned int> &numbers,
 template <int dim>
 void
 FE_DGQ<dim>::
-get_interpolation_matrix (const FiniteElementBase<dim> &x_source_fe,
+get_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
                          FullMatrix<double>           &interpolation_matrix) const
 {
                                   // this is only implemented, if the
@@ -370,7 +370,7 @@ get_interpolation_matrix (const FiniteElementBase<dim> &x_source_fe,
   AssertThrow ((x_source_fe.get_name().find ("FE_DGQ<") == 0)
                ||
                (dynamic_cast<const FE_DGQ<dim>*>(&x_source_fe) != 0),
-               typename FiniteElementBase<dim>::
+               typename FiniteElement<dim>::
                ExcInterpolationNotImplemented());
   
                                   // ok, source is a Q element, so
index 1849c87baad5ffb4f312fa22bf6f75acb2ae04a9..7c92ab7d751caac8cdcd44631e26d2034242fc4f 100644 (file)
@@ -40,7 +40,7 @@ double
 FE_PolyTensor<POLY,dim>::shape_value (
   const unsigned int, const Point<dim> &) const
 {
-  Assert(false, typename FiniteElementBase<dim>::ExcFENotPrimitive());
+  Assert(false, typename FiniteElement<dim>::ExcFENotPrimitive());
   return 0.;
 }
 
@@ -78,7 +78,7 @@ Tensor<1,dim>
 FE_PolyTensor<POLY,dim>::shape_grad (
   const unsigned int, const Point<dim> &) const
 {
-  Assert(false, typename FiniteElementBase<dim>::ExcFENotPrimitive());
+  Assert(false, typename FiniteElement<dim>::ExcFENotPrimitive());
   return Tensor<1,dim>();
 }
 
@@ -118,7 +118,7 @@ Tensor<2,dim>
 FE_PolyTensor<POLY,dim>::shape_grad_grad (
   const unsigned int, const Point<dim> &) const
 {
-  Assert(false, typename FiniteElementBase<dim>::ExcFENotPrimitive());
+  Assert(false, typename FiniteElement<dim>::ExcFENotPrimitive());
   return Tensor<2,dim>();
 }
 
index d534e4e83f981143731dcdd754e38afe38d43a36..bdddc26a3dea91fbcf5856b013dcce73c5f406d5 100644 (file)
@@ -249,7 +249,7 @@ FE_Q<dim>::clone() const
 template <int dim>
 void
 FE_Q<dim>::
-get_interpolation_matrix (const FiniteElementBase<dim> &x_source_fe,
+get_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
                          FullMatrix<double>           &interpolation_matrix) const
 {
                                   // this is only implemented, if the
@@ -258,7 +258,7 @@ get_interpolation_matrix (const FiniteElementBase<dim> &x_source_fe,
   AssertThrow ((x_source_fe.get_name().find ("FE_Q<") == 0)
                ||
                (dynamic_cast<const FE_Q<dim>*>(&x_source_fe) != 0),
-               typename FiniteElementBase<dim>::
+               typename FiniteElement<dim>::
                ExcInterpolationNotImplemented());
   
                                   // ok, source is a Q element, so
@@ -607,7 +607,7 @@ FE_Q<3>::initialize_constraints ()
                                   // In the following the points x_i
                                   // are constructed in the order as
                                   // described in the documentation
-                                  // of the FiniteElementBase class
+                                  // of the FiniteElement class
                                   // (fe_base.h), i.e.
                                   //   *--13--3--14--*
                                   //   |      |      |
@@ -765,7 +765,7 @@ FE_Q<3>::initialize_constraints ()
                                            //
                                            // For a different explanation of
                                            // the problem, see the discussion
-                                           // in the FiniteElementBase class
+                                           // in the FiniteElement class
                                            // for constraint matrices in 3d.
           mirror[k] = (constraint_point(k) > 0.5);
           if (mirror[k])
index e4b1df845ec3a9bd9c8825c292e33b3c3ff8463d..6becd13f1a0ab8228f76644cd7af99661dd7a8b0 100644 (file)
@@ -287,7 +287,7 @@ FE_RaviartThomas<dim>::shape_grad_grad_component (const unsigned int i,
 template <>
 void
 FE_RaviartThomas<1>::
-get_interpolation_matrix (const FiniteElementBase<1> &,
+get_interpolation_matrix (const FiniteElement<1> &,
                          FullMatrix<double>         &) const
 {
   Assert (false, ExcImpossibleInDim(1));
@@ -299,7 +299,7 @@ get_interpolation_matrix (const FiniteElementBase<1> &,
 template <int dim>
 void
 FE_RaviartThomas<dim>::
-get_interpolation_matrix (const FiniteElementBase<dim> &x_source_fe,
+get_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
                          FullMatrix<double>           &interpolation_matrix) const
 {
                                   // this is only implemented, if the
@@ -310,7 +310,7 @@ get_interpolation_matrix (const FiniteElementBase<dim> &x_source_fe,
   AssertThrow ((x_source_fe.get_name().find ("FE_RaviartThomas<") == 0)
                ||
                (dynamic_cast<const FE_RaviartThomas<dim>*>(&x_source_fe) != 0),
-               typename FiniteElementBase<dim>::
+               typename FiniteElement<dim>::
                ExcInterpolationNotImplemented());
   
                                   // ok, source is a RT element, so
index f120aa76f1f446ea217391a65f3ce04453cb770b..9ed9b15973727c23ef409b47ee6a9a30ff1af0d8 100644 (file)
@@ -62,7 +62,7 @@ FESystem<dim>::InternalData::~InternalData()
 
 
 template <int dim>
-typename FiniteElementBase<dim>::InternalDataBase &
+typename FiniteElement<dim>::InternalDataBase &
 FESystem<dim>::
 InternalData::get_fe_data (const unsigned int base_no) const
 {
@@ -77,7 +77,7 @@ template <int dim>
 void
 FESystem<dim>::
 InternalData::set_fe_data (const unsigned int base_no,
-                          typename FiniteElementBase<dim>::InternalDataBase *ptr)
+                          typename FiniteElement<dim>::InternalDataBase *ptr)
 {
   Assert(base_no<base_fe_datas.size(),
         ExcIndexRange(base_no,0,base_fe_datas.size()));
@@ -132,7 +132,7 @@ FESystem<dim>::InternalData::clear_first_cell ()
 {
                                    // call respective function of base
                                    // class
-  FiniteElementBase<dim>::InternalDataBase::clear_first_cell ();
+  FiniteElement<dim>::InternalDataBase::clear_first_cell ();
                                    // then the functions of all the
                                    // sub-objects
   for (unsigned int i=0; i<base_fe_datas.size(); ++i)
@@ -299,7 +299,7 @@ FESystem<dim>::shape_value (const unsigned int i,
 {
   Assert (i<this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
   Assert (this->is_primitive(i), 
-         typename FiniteElementBase<dim>::ExcShapeFunctionNotPrimitive(i));
+         typename FiniteElement<dim>::ExcShapeFunctionNotPrimitive(i));
 
   return (base_element(this->system_to_base_table[i].first.first)
          .shape_value(this->system_to_base_table[i].second, p));
@@ -351,7 +351,7 @@ FESystem<dim>::shape_grad (const unsigned int i,
 {
   Assert (i<this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
   Assert (this->is_primitive(i),
-         typename FiniteElementBase<dim>::ExcShapeFunctionNotPrimitive(i));
+         typename FiniteElement<dim>::ExcShapeFunctionNotPrimitive(i));
 
   return (base_element(this->system_to_base_table[i].first.first)
          .shape_grad(this->system_to_base_table[i].second, p));
@@ -403,7 +403,7 @@ FESystem<dim>::shape_grad_grad (const unsigned int i,
 {
   Assert (i<this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
   Assert (this->is_primitive(i), 
-         typename FiniteElementBase<dim>::ExcShapeFunctionNotPrimitive(i));
+         typename FiniteElement<dim>::ExcShapeFunctionNotPrimitive(i));
 
   return (base_element(this->system_to_base_table[i].first.first)
          .shape_grad_grad(this->system_to_base_table[i].second, p));
@@ -451,7 +451,7 @@ FESystem<dim>::shape_grad_grad_component (const unsigned int i,
 template <int dim>
 void
 FESystem<dim>::
-get_interpolation_matrix (const FiniteElementBase<dim> &x_source_fe,
+get_interpolation_matrix (const FiniteElement<dim> &x_source_fe,
                          FullMatrix<double>           &interpolation_matrix) const
 {
   Assert (interpolation_matrix.m() == this->dofs_per_cell,
@@ -470,7 +470,7 @@ get_interpolation_matrix (const FiniteElementBase<dim> &x_source_fe,
   AssertThrow ((x_source_fe.get_name().find ("FESystem<") == 0)
                ||
                (dynamic_cast<const FESystem<dim>*>(&x_source_fe) != 0),
-               typename FiniteElementBase<dim>::
+               typename FiniteElement<dim>::
                ExcInterpolationNotImplemented());
   
                                   // ok, source is a system element,
@@ -481,7 +481,7 @@ get_interpolation_matrix (const FiniteElementBase<dim> &x_source_fe,
                                    // condition 2: same number of
                                    // basis elements
   AssertThrow (n_base_elements() == source_fe.n_base_elements(),
-               typename FiniteElementBase<dim>::
+               typename FiniteElement<dim>::
                ExcInterpolationNotImplemented());
 
                                    // condition 3: same number of
@@ -489,7 +489,7 @@ get_interpolation_matrix (const FiniteElementBase<dim> &x_source_fe,
   for (unsigned int i=0; i<n_base_elements(); ++i)
     AssertThrow (element_multiplicity(i) ==
                  source_fe.element_multiplicity(i),
-                 typename FiniteElementBase<dim>::
+                 typename FiniteElement<dim>::
                  ExcInterpolationNotImplemented());
 
                                    // ok, so let's try whether it
@@ -617,8 +617,8 @@ FESystem<dim>::get_data (const UpdateFlags      flags_,
       typename Mapping<dim>::InternalDataBase *base_fe_data_base =
        base_element(base_no).get_data(sub_flags, mapping, quadrature);
 
-      typename FiniteElementBase<dim>::InternalDataBase *base_fe_data =
-       dynamic_cast<typename FiniteElementBase<dim>::InternalDataBase *>
+      typename FiniteElement<dim>::InternalDataBase *base_fe_data =
+       dynamic_cast<typename FiniteElement<dim>::InternalDataBase *>
        (base_fe_data_base);
       
       data->set_fe_data(base_no, base_fe_data);
@@ -823,7 +823,7 @@ compute_fill (const Mapping<dim>                   &mapping,
        {
          const FiniteElement<dim> &
             base_fe      = base_element(base_no);
-         typename FiniteElementBase<dim>::InternalDataBase &
+         typename FiniteElement<dim>::InternalDataBase &
            base_fe_data = fe_data.get_fe_data(base_no);
          FEValuesData<dim> &
             base_data    = fe_data.get_fe_values_data(base_no);
@@ -2389,7 +2389,7 @@ FESystem<dim>::unit_support_point (const unsigned index) const
           ExcIndexRange (index, 0, this->dofs_per_cell));
   Assert ((this->unit_support_points.size() == this->dofs_per_cell) ||
           (this->unit_support_points.size() == 0),
-          typename FiniteElementBase<dim>::ExcFEHasNoSupportPoints ());
+          typename FiniteElement<dim>::ExcFEHasNoSupportPoints ());
 
                                    // let's see whether we have the
                                    // information pre-computed
@@ -2413,7 +2413,7 @@ FESystem<dim>::unit_face_support_point (const unsigned index) const
           ExcIndexRange (index, 0, this->dofs_per_face));
   Assert ((this->unit_face_support_points.size() == this->dofs_per_face) ||
           (this->unit_face_support_points.size() == 0),
-          typename FiniteElementBase<dim>::ExcFEHasNoSupportPoints ());
+          typename FiniteElement<dim>::ExcFEHasNoSupportPoints ());
 
                                    // let's see whether we have the
                                    // information pre-computed
index 92499e2ac0c849ec8e3f3784518a81bb373cbc07..961873332f84aec5ce8c1df05e6f8dce7275f3f8 100644 (file)
@@ -298,7 +298,7 @@ void FETools::get_interpolation_matrix (const FiniteElement<dim> &fe1,
     {
       gim_forwarder (fe1, fe2, interpolation_matrix);
     }
-  catch (typename FiniteElementBase<dim>::ExcInterpolationNotImplemented &)
+  catch (typename FiniteElement<dim>::ExcInterpolationNotImplemented &)
     {
                                        // too bad....
       fe_implements_interpolation = false;
@@ -321,7 +321,7 @@ void FETools::get_interpolation_matrix (const FiniteElement<dim> &fe1,
     fe2_support_points = fe2.get_unit_support_points ();
 
   Assert(fe2_support_points.size()==fe2.dofs_per_cell,
-        typename FiniteElementBase<dim>::ExcFEHasNoSupportPoints());
+        typename FiniteElement<dim>::ExcFEHasNoSupportPoints());
 
   for (unsigned int i=0; i<fe2.dofs_per_cell; ++i)     
     {
index 2ad757b59d810f61657c51805afc222f356fd54b..ab3f43389a14e14b1bb9d1baec0cc09ca2252cbe 100644 (file)
@@ -632,6 +632,13 @@ inconvenience this causes.
 
 <ol>
 
+  <li> <p> Removed: The class <code
+  class="class">FiniteElementBase</code> has been removed and all its
+  functions are now in <code class="class">FiniteElement</code>.
+  <br>
+  (GK, 2005/08/25)
+  </p>
+
   <li> <p> New: class <code class="class">DoFTools</code> now has two
   functions <code class="member">compute_row_length_vector</code>, one
   for equations and one for systems. These give a much fine estimate

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.