]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
New (intermediate) FiniteElement class which allows implementation of FiniteElements...
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 16 Mar 2004 16:35:30 +0000 (16:35 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 16 Mar 2004 16:35:30 +0000 (16:35 +0000)
git-svn-id: https://svn.dealii.org/trunk@8775 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_poly.h [new file with mode: 0644]
deal.II/deal.II/source/fe/fe_poly.cc [new file with mode: 0644]

diff --git a/deal.II/deal.II/include/fe/fe_poly.h b/deal.II/deal.II/include/fe/fe_poly.h
new file mode 100644 (file)
index 0000000..1bdb386
--- /dev/null
@@ -0,0 +1,384 @@
+//---------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2004 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------
+#ifndef __deal2__fe_poly_h
+#define __deal2__fe_poly_h
+
+
+#include <fe/fe.h>
+
+
+/**
+ * This class gives a unified framework for implementing a
+ * FiniteElement class based on a polynomial space like a
+ * @p TensorProductSpace or a @p PolynomialSpace.
+ *
+ * Every class that implements following functions can be used as
+ * template parameter @p POLY. Example classes are @p
+ * TensorProductSpace and @p PolynomialSpace.
+ *
+ * @code
+ * double compute_value (const unsigned int i,
+ *                       const Point<dim> &p) const;
+ *
+ * Tensor<1,dim> compute_grad (const unsigned int i,
+ *                             const Point<dim> &p) const;
+ *
+ * Tensor<2,dim> compute_grad_grad (const unsigned int i,
+ *                                  const Point<dim> &p) const;
+ * @endcode
+ *
+ * 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
+ * implemented by this class but are left for implementation in
+ * derived classes.
+ *
+ * Todos:
+ * - checke dim of POLY
+ * - templatisiere nur auf POLY, dim ergibt sich durch POLY::dim
+ *
+ * @author Ralf Hartmann 2004
+ **/
+template <class POLY, int dim/*=POLY::dim*/>
+class FE_Poly : public FiniteElement<dim>
+{
+  public:
+                                    /**
+                                     * Constructor.
+                                     */
+    FE_Poly (const POLY& poly_space,
+            const FiniteElementData<dim> &fe_data,
+            const std::vector<bool> &restriction_is_additive_flags,
+            const std::vector<std::vector<bool> > &nonzero_components);
+
+                                    /**
+                                     * Return the value of the
+                                     * @p{i}th shape function at the
+                                     * point @p{p}. See the
+                                     * @ref{FiniteElementBase} base
+                                     * class for more information
+                                     * about the semantics of this
+                                     * function.
+                                     */
+    virtual double shape_value (const unsigned int i,
+                               const Point<dim> &p) const;
+    
+                                    /**
+                                     * Return the value of the
+                                     * @p{component}th vector
+                                     * component of the @p{i}th shape
+                                     * function at the point
+                                     * @p{p}. See the
+                                     * @ref{FiniteElementBase} base
+                                     * class for more information
+                                     * about the semantics of this
+                                     * function.
+                                     *
+                                     * Since this element is scalar,
+                                     * the returned value is the same
+                                     * as if the function without the
+                                     * @p{_component} suffix were
+                                     * called, provided that the
+                                     * specified component is zero.
+                                     */
+    virtual double shape_value_component (const unsigned int i,
+                                         const Point<dim> &p,
+                                         const unsigned int component) const;
+
+                                    /**
+                                     * Return the gradient of the
+                                     * @p{i}th shape function at the
+                                     * point @p{p}. See the
+                                     * @ref{FiniteElementBase} base
+                                     * class for more information
+                                     * about the semantics of this
+                                     * function.
+                                     */
+    virtual Tensor<1,dim> shape_grad (const unsigned int  i,
+                                     const Point<dim>   &p) const;
+
+                                    /**
+                                     * Return the gradient of the
+                                     * @p{component}th vector
+                                     * component of the @p{i}th shape
+                                     * function at the point
+                                     * @p{p}. See the
+                                     * @ref{FiniteElementBase} base
+                                     * class for more information
+                                     * about the semantics of this
+                                     * function.
+                                     *
+                                     * Since this element is scalar,
+                                     * the returned value is the same
+                                     * as if the function without the
+                                     * @p{_component} suffix were
+                                     * called, provided that the
+                                     * specified component is zero.
+                                     */
+    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{i}th
+                                     * shape function at point @p{p}
+                                     * on the unit cell. See the
+                                     * @ref{FiniteElementBase} base
+                                     * class for more information
+                                     * about the semantics of this
+                                     * function.
+                                     */
+    virtual Tensor<2,dim> shape_grad_grad (const unsigned int  i,
+                                          const Point<dim> &p) const;
+
+                                    /**
+                                     * Return the second derivative
+                                     * of the @p{component}th vector
+                                     * component of the @p{i}th shape
+                                     * function at the point
+                                     * @p{p}. See the
+                                     * @ref{FiniteElementBase} base
+                                     * class for more information
+                                     * about the semantics of this
+                                     * function.
+                                     *
+                                     * Since this element is scalar,
+                                     * the returned value is the same
+                                     * as if the function without the
+                                     * @p{_component} suffix were
+                                     * called, provided that the
+                                     * specified component is zero.
+                                     */
+    virtual Tensor<2,dim> shape_grad_grad_component (const unsigned int i,
+                                                    const Point<dim> &p,
+                                                    const unsigned int component) const;
+
+                                     /**
+                                     * Number of base elements in a
+                                     * mixed discretization. Since
+                                     * this is a scalar element,
+                                     * return one.
+                                     */
+    virtual unsigned int n_base_elements () const;
+    
+                                    /**
+                                     * Access to base element
+                                     * objects. Since this element is
+                                     * scalar, @p{base_element(0)} is
+                                     * @p{this}, and all other
+                                     * indices throw an error.
+                                     */
+    virtual const FiniteElement<dim> &
+    base_element (const unsigned int index) const;
+
+                                     /**
+                                      * Multiplicity of base element
+                                      * @p{index}. Since this is a
+                                      * scalar element,
+                                      * @p{element_multiplicity(0)}
+                                      * returns one, and all other
+                                      * indices will throw an error.
+                                      */
+    virtual unsigned int element_multiplicity (const unsigned int index) const;
+
+    
+  protected:
+      
+                                    /**
+                                     * Prepare internal data
+                                     * structures and fill in values
+                                     * independent of the cell.
+                                     */
+    virtual
+    typename Mapping<dim>::InternalDataBase *
+    get_data (const UpdateFlags,
+             const Mapping<dim>& mapping,
+             const Quadrature<dim>& quadrature) const ;
+
+                                    /**
+                                     * Implementation of the same
+                                     * function in
+                                     * @ref{FiniteElement}.
+                                     */
+    virtual void
+    fill_fe_values (const Mapping<dim> &mapping,
+                   const typename DoFHandler<dim>::cell_iterator &cell,
+                   const Quadrature<dim>                &quadrature,
+                   typename Mapping<dim>::InternalDataBase      &mapping_internal,
+                   typename Mapping<dim>::InternalDataBase      &fe_internal,
+                   FEValuesData<dim>& data) const;
+    
+                                    /**
+                                     * Implementation of the same
+                                     * function in
+                                     * @ref{FiniteElement}.
+                                     */
+    virtual void
+    fill_fe_face_values (const Mapping<dim> &mapping,
+                        const typename DoFHandler<dim>::cell_iterator &cell,
+                        const unsigned int                    face_no,
+                        const Quadrature<dim-1>                &quadrature,
+                        typename Mapping<dim>::InternalDataBase      &mapping_internal,
+                        typename Mapping<dim>::InternalDataBase      &fe_internal,
+                        FEValuesData<dim>& data) const ;
+    
+                                    /**
+                                     * Implementation of the same
+                                     * function in
+                                     * @ref{FiniteElement}.
+                                     */
+    virtual void
+    fill_fe_subface_values (const Mapping<dim> &mapping,
+                           const typename DoFHandler<dim>::cell_iterator &cell,
+                           const unsigned int                    face_no,
+                           const unsigned int                    sub_no,
+                           const Quadrature<dim-1>                &quadrature,
+                           typename Mapping<dim>::InternalDataBase      &mapping_internal,
+                           typename Mapping<dim>::InternalDataBase      &fe_internal,
+                           FEValuesData<dim>& data) const ;
+
+    
+                                    /**
+                                     * Determine the values that need
+                                     * to be computed on the unit
+                                     * cell to be able to compute all
+                                     * values required by @p{flags}.
+                                     *
+                                     * For the purpuse of this
+                                     * function, refer to the
+                                     * documentation in
+                                     * @p{FiniteElement}.
+                                     *
+                                     * The effect in this element is
+                                     * as follows: if
+                                     * @p{update_values} is set in
+                                     * @p{flags}, copy it to the
+                                     * result. All other flags of the
+                                     * result are cleared, since
+                                     * everything else must be
+                                     * computed for each cell.
+                                     */
+    virtual UpdateFlags update_once (const UpdateFlags flags) const;
+  
+                                    /**
+                                     * Determine the values that need
+                                     * to be computed on every
+                                     * cell to be able to compute all
+                                     * values required by @p{flags}.
+                                     *
+                                     * For the purpuse of this
+                                     * function, refer to the
+                                     * documentation in
+                                     * @p{FiniteElement}.
+                                     *
+                                     * The effect in this element is
+                                     * as follows:
+                                     * @begin{itemize}
+                                     * @item if @p{update_gradients}
+                                     * is set, the result will
+                                     * contain @p{update_gradients}
+                                     * and
+                                     * @p{update_covariant_transformation}.
+                                     * The latter is required to
+                                     * transform the gradient on the
+                                     * unit cell to the real
+                                     * cell. Remark, that the action
+                                     * required by
+                                     * @p{update_covariant_transformation}
+                                     * is actually performed by the
+                                     * @p{Mapping} object used in
+                                     * conjunction with this finite
+                                     * element.
+                                     * @item if
+                                     * @p{update_second_derivatives}
+                                     * is set, the result will
+                                     * contain
+                                     * @p{update_second_derivatives}
+                                     * and
+                                     * @p{update_covariant_transformation}.
+                                     * The rationale is the same as
+                                     * above and no higher
+                                     * derivatives of the
+                                     * transformation are required,
+                                     * since we use difference
+                                     * quotients for the actual
+                                     * computation.
+                                     * @end{itemize}
+                                     */
+    virtual UpdateFlags update_each (const UpdateFlags flags) const;
+
+
+                                    /**
+                                     * Fields of cell-independent data.
+                                     *
+                                     * For information about the
+                                     * general purpose of this class,
+                                     * see the documentation of the
+                                     * base class.
+                                     */
+    class InternalData : public FiniteElementBase<dim>::InternalDataBase
+    {
+      public:
+                                        /**
+                                         * Array with shape function
+                                         * values in quadrature
+                                         * points. There is one
+                                         * row for each shape
+                                         * function, containing
+                                         * values for each quadrature
+                                         * point.
+                                         *
+                                         * In this array, we store
+                                         * the values of the shape
+                                         * function in the quadrature
+                                         * points on the unit
+                                         * cell. Since these values
+                                         * do not change under
+                                         * transformation to the real
+                                         * cell, we only need to copy
+                                         * them over when visiting a
+                                         * concrete cell.
+                                         */
+       Table<2,double> shape_values;
+
+                                        /**
+                                         * Array with shape function
+                                         * gradients in quadrature
+                                         * points. There is one
+                                         * row for each shape
+                                         * function, containing
+                                         * values for each quadrature
+                                         * point.
+                                         *
+                                         * We store the gradients in
+                                         * the quadrature points on
+                                         * the unit cell. We then
+                                         * only have to apply the
+                                         * transformation (which is a
+                                         * matrix-vector
+                                         * multiplication) when
+                                         * visiting an actual cell.
+                                         */      
+       Table<2,Tensor<1,dim> > shape_gradients;
+    };
+
+                                     /**
+                                      * The polynomial space.
+                                      */    
+    POLY poly_space;
+};
+
+
+
+#endif
diff --git a/deal.II/deal.II/source/fe/fe_poly.cc b/deal.II/deal.II/source/fe/fe_poly.cc
new file mode 100644 (file)
index 0000000..8bde074
--- /dev/null
@@ -0,0 +1,419 @@
+//----------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2004 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------------------------------------------
+
+//  #include <base/quadrature.h>
+//  #include <base/polynomial.h>
+//  #include <base/tensor_product_polynomials.h>
+//  #include <grid/tria.h>
+//  #include <grid/tria_iterator.h>
+//  #include <dofs/dof_accessor.h>
+#include <base/quadrature.h>
+#include <base/tensor_product_polynomials.h>
+#include <base/polynomials_p.h>
+#include <fe/fe_poly.h>
+//#include <fe/fe.h>
+//  #include <fe/mapping.h>
+//  #include <fe/fe_q.h>
+#include <fe/fe_values.h>
+
+//  #ifdef HAVE_STD_STRINGSTREAM
+//  #  include <sstream>
+//  #else
+//  #  include <strstream>
+//  #endif
+
+
+template <class POLY, int dim>
+FE_Poly<POLY,dim>::FE_Poly (const POLY& poly_space,
+                           const FiniteElementData<dim> &fe_data,
+                           const std::vector<bool> &restriction_is_additive_flags,
+                           const std::vector<std::vector<bool> > &nonzero_components):
+               FiniteElement<dim> (fe_data,
+                                   restriction_is_additive_flags,
+                                   nonzero_components),
+               poly_space(poly_space)
+{}
+
+
+
+template <class POLY, int dim>
+double
+FE_Poly<POLY,dim>::shape_value (const unsigned int i,
+                               const Point<dim> &p) const
+{
+  Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
+  return poly_space.compute_value(i, p);
+}
+
+
+template <class POLY, int dim>
+double
+FE_Poly<POLY,dim>::shape_value_component (const unsigned int i,
+                                          const Point<dim> &p,
+                                          const unsigned int component) const
+{
+  Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
+  Assert (component == 0, ExcIndexRange (component, 0, 1));
+  return poly_space.compute_value(i, p);
+}
+
+
+
+template <class POLY, int dim>
+Tensor<1,dim>
+FE_Poly<POLY,dim>::shape_grad (const unsigned int i,
+                               const Point<dim> &p) const
+{
+  Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
+  return poly_space.compute_grad(i, p);
+}
+
+
+
+template <class POLY, int dim>
+Tensor<1,dim>
+FE_Poly<POLY,dim>::shape_grad_component (const unsigned int i,
+                                         const Point<dim> &p,
+                                         const unsigned int component) const
+{
+  Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
+  Assert (component == 0, ExcIndexRange (component, 0, 1));
+  return poly_space.compute_grad(i, p);
+}
+
+
+
+template <class POLY, int dim>
+Tensor<2,dim>
+FE_Poly<POLY,dim>::shape_grad_grad (const unsigned int i,
+                                    const Point<dim> &p) const
+{
+  Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
+  return poly_space.compute_grad_grad(i, p);
+}
+
+
+
+template <class POLY, int dim>
+Tensor<2,dim>
+FE_Poly<POLY,dim>::shape_grad_grad_component (const unsigned int i,
+                                              const Point<dim> &p,
+                                              const unsigned int component) const
+{
+  Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
+  Assert (component == 0, ExcIndexRange (component, 0, 1));
+  return poly_space.compute_grad_grad(i, p);
+}
+
+
+
+//----------------------------------------------------------------------
+// Auxiliary functions
+//----------------------------------------------------------------------
+
+
+
+
+template <class POLY, int dim>
+UpdateFlags
+FE_Poly<POLY,dim>::update_once (const UpdateFlags flags) const
+{
+                                  // for this kind of elements, only
+                                  // the values can be precomputed
+                                  // once and for all. set this flag
+                                  // if the values are requested at
+                                  // all
+  return (update_default | (flags & update_values));
+}
+
+
+
+template <class POLY, int dim>
+UpdateFlags
+FE_Poly<POLY,dim>::update_each (const UpdateFlags flags) const
+{
+  UpdateFlags out = update_default;
+
+  if (flags & update_gradients)
+    out |= update_gradients | update_covariant_transformation;
+  if (flags & update_second_derivatives)
+    out |= update_second_derivatives | update_covariant_transformation;
+
+  return out;
+}
+
+
+
+//----------------------------------------------------------------------
+// Data field initialization
+//----------------------------------------------------------------------
+
+template <class POLY, int dim>
+typename Mapping<dim>::InternalDataBase *
+FE_Poly<POLY,dim>::get_data (const UpdateFlags      update_flags,
+                            const Mapping<dim>    &mapping,
+                            const Quadrature<dim> &quadrature) const
+{
+                                  // generate a new data object and
+                                  // initialize some fields
+  InternalData* data = new InternalData;
+
+                                  // check what needs to be
+                                  // initialized only once and what
+                                  // on every cell/face/subface we
+                                  // visit
+  data->update_once = update_once(update_flags);
+  data->update_each = update_each(update_flags);
+  data->update_flags = data->update_once | data->update_each;
+
+  const UpdateFlags flags(data->update_flags);
+  const unsigned int n_q_points = quadrature.n_quadrature_points;
+
+                                  // some scratch arrays
+  std::vector<double> values(0);
+  std::vector<Tensor<1,dim> > grads(0);
+  std::vector<Tensor<2,dim> > grad_grads(0);
+
+                                  // initialize fields only if really
+                                  // necessary. otherwise, don't
+                                  // allocate memory
+  if (flags & update_values)
+    {
+      values.resize (this->dofs_per_cell);
+      data->shape_values.reinit (this->dofs_per_cell,
+                                n_q_points);
+    }
+
+  if (flags & update_gradients)
+    {
+      grads.resize (this->dofs_per_cell);
+      data->shape_gradients.reinit (this->dofs_per_cell,
+                                   n_q_points);
+    }
+
+                                  // if second derivatives through
+                                  // finite differencing is required,
+                                  // then initialize some objects for
+                                  // that
+  if (flags & update_second_derivatives)
+    data->initialize_2nd (this, mapping, quadrature);
+
+                                  // next already fill those fields
+                                  // of which we have information by
+                                  // now. note that the shape
+                                  // gradients are only those on the
+                                  // unit cell, and need to be
+                                  // transformed when visiting an
+                                  // actual cell
+  if (flags & (update_values | update_gradients))
+    for (unsigned int i=0; i<n_q_points; ++i)
+      {
+       poly_space.compute(quadrature.point(i),
+                          values, grads, grad_grads);
+       
+       if (flags & update_values)
+         for (unsigned int k=0; k<this->dofs_per_cell; ++k)
+           data->shape_values[k][i] = values[k];
+       
+       if (flags & update_gradients)
+         for (unsigned int k=0; k<this->dofs_per_cell; ++k)
+           data->shape_gradients[k][i] = grads[k];
+      }
+  return data;
+}
+
+
+
+
+//----------------------------------------------------------------------
+// Fill data of FEValues
+//----------------------------------------------------------------------
+
+template <class POLY, int dim>
+void
+FE_Poly<POLY,dim>::fill_fe_values (const Mapping<dim>                   &mapping,
+                          const typename DoFHandler<dim>::cell_iterator &cell,
+                          const Quadrature<dim>                &quadrature,
+                          typename Mapping<dim>::InternalDataBase &mapping_data,
+                          typename Mapping<dim>::InternalDataBase &fedata,
+                          FEValuesData<dim>                    &data) const
+{
+                                  // convert data object to internal
+                                  // data for this class. fails with
+                                  // an exception if that is not
+                                  // possible
+  InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+  
+  const UpdateFlags flags(fe_data.current_update_flags());
+
+  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
+    {
+      if (flags & update_values)
+       for (unsigned int i=0; i<quadrature.n_quadrature_points; ++i)
+         data.shape_values(k,i) = fe_data.shape_values[k][i];
+      
+      if (flags & update_gradients)
+       {
+         Assert (data.shape_gradients[k].size() <=
+                 fe_data.shape_gradients[k].size(),
+                 ExcInternalError());    
+         mapping.transform_covariant(data.shape_gradients[k].begin(),
+                                     data.shape_gradients[k].end(),
+                                     fe_data.shape_gradients[k].begin(),
+                                     mapping_data);
+       };
+    }
+
+  const typename QProjector<dim>::DataSetDescriptor dsd;
+  if (flags & update_second_derivatives)
+    this->compute_2nd (mapping, cell, dsd.cell(),
+                      mapping_data, fe_data, data);
+}
+
+
+
+template <class POLY, int dim>
+void
+FE_Poly<POLY,dim>::fill_fe_face_values (const Mapping<dim>                   &mapping,
+                               const typename DoFHandler<dim>::cell_iterator &cell,
+                               const unsigned int                    face,
+                               const Quadrature<dim-1>              &quadrature,
+                               typename Mapping<dim>::InternalDataBase       &mapping_data,
+                               typename Mapping<dim>::InternalDataBase       &fedata,
+                               FEValuesData<dim>                    &data) const
+{
+                                  // convert data object to internal
+                                  // data for this class. fails with
+                                  // an exception if that is not
+                                  // possible
+  InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+
+                                  // offset determines which data set
+                                  // to take (all data sets for all
+                                  // faces are stored contiguously)
+  
+  const typename QProjector<dim>::DataSetDescriptor dsd;
+  const typename QProjector<dim>::DataSetDescriptor offset
+    = dsd.face (face, cell->face_orientation(face),
+               quadrature.n_quadrature_points);
+  
+  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
+
+  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
+    {
+      if (flags & update_values)
+        for (unsigned int i=0; i<quadrature.n_quadrature_points; ++i)
+         data.shape_values(k,i) = fe_data.shape_values[k][i+offset];
+      
+      if (flags & update_gradients)
+       {
+         Assert (data.shape_gradients[k].size() + offset <=
+                 fe_data.shape_gradients[k].size(),
+                 ExcInternalError());
+         mapping.transform_covariant(data.shape_gradients[k].begin(),
+                                     data.shape_gradients[k].end(),
+                                     fe_data.shape_gradients[k].begin()+offset,
+                                     mapping_data);
+       };
+    }
+
+  if (flags & update_second_derivatives)
+    this->compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
+}
+
+
+
+template <class POLY, int dim>
+void
+FE_Poly<POLY,dim>::fill_fe_subface_values (const Mapping<dim>                   &mapping,
+                                  const typename DoFHandler<dim>::cell_iterator &cell,
+                                  const unsigned int                    face,
+                                  const unsigned int                    subface,
+                                  const Quadrature<dim-1>              &quadrature,
+                                  typename Mapping<dim>::InternalDataBase       &mapping_data,
+                                  typename Mapping<dim>::InternalDataBase       &fedata,
+                                  FEValuesData<dim>                    &data) const
+{
+                                  // convert data object to internal
+                                  // data for this class. fails with
+                                  // an exception if that is not
+                                  // possible
+  InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+
+                                  // offset determines which data set
+                                  // to take (all data sets for all
+                                // sub-faces are stored contiguously)
+
+  const typename QProjector<dim>::DataSetDescriptor dsd;
+  const typename QProjector<dim>::DataSetDescriptor offset
+    = dsd.sub_face (face, subface, cell->face_orientation(face),
+                   quadrature.n_quadrature_points);
+
+  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
+
+  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
+    {
+      if (flags & update_values)
+        for (unsigned int i=0; i<quadrature.n_quadrature_points; ++i)
+         data.shape_values(k,i) = fe_data.shape_values[k][i+offset];
+      
+      if (flags & update_gradients)
+       {
+         Assert (data.shape_gradients[k].size() + offset <=
+                 fe_data.shape_gradients[k].size(),
+                 ExcInternalError());
+         mapping.transform_covariant(data.shape_gradients[k].begin(),
+                                     data.shape_gradients[k].end(),
+                                     fe_data.shape_gradients[k].begin()+offset,
+                                     mapping_data);
+       };
+    }
+  
+  if (flags & update_second_derivatives)
+    this->compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
+}
+
+
+
+template <class POLY, int dim>
+unsigned int
+FE_Poly<POLY,dim>::n_base_elements () const
+{
+  return 1;
+}
+
+
+
+template <class POLY, int dim>
+const FiniteElement<dim> &
+FE_Poly<POLY,dim>::base_element (const unsigned int index) const
+{
+  Assert (index==0, ExcIndexRange(index, 0, 1));
+  return *this;
+}
+
+
+
+template <class POLY, int dim>
+unsigned int
+FE_Poly<POLY,dim>::element_multiplicity (const unsigned int index) const
+{
+  Assert (index==0, ExcIndexRange(index, 0, 1));
+  return 1;
+}
+
+
+
+template class FE_Poly<TensorProductPolynomials<deal_II_dimension>, deal_II_dimension>;
+template class FE_Poly<PolynomialSpace<deal_II_dimension>, deal_II_dimension>;
+template class FE_Poly<PolynomialsP<deal_II_dimension>, deal_II_dimension>;

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.