]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
add FE_FaceQ and FE_PolyFace for testing
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 25 Nov 2009 05:22:16 +0000 (05:22 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 25 Nov 2009 05:22:16 +0000 (05:22 +0000)
git-svn-id: https://svn.dealii.org/trunk@20164 0785d39b-7218-0410-832d-ea1e28bc413d

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

index a50375e7af0043cef331abac7d903715e99aa84f..e0b44726b1f5109fa10d58c4a1351575f5719370 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
+//    Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2009 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -13,9 +13,6 @@
 #ifndef __deal2__fe_dgq_h
 #define __deal2__fe_dgq_h
 
-
-//TODO: we should probably turn the derivation around, i.e. make FE_DGQ a special case of FE_DGQ_Generic/Arbitrary, and then have other classes that only have the constructor and clone() that use other sets of quadrature nodes.
-
 #include <base/config.h>
 #include <base/tensor_product_polynomials.h>
 #include <fe/fe_poly.h>
diff --git a/deal.II/deal.II/include/fe/fe_face.h b/deal.II/deal.II/include/fe/fe_face.h
new file mode 100644 (file)
index 0000000..97cab1f
--- /dev/null
@@ -0,0 +1,86 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2009 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_face_h
+#define __deal2__fe_face_h
+
+#include <base/config.h>
+#include <base/tensor_product_polynomials.h>
+#include <fe/fe_poly_face.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+
+/**
+ * @warning This class has not been tested
+ *
+ * A finite element, which is a tensor product polynomial on each face
+ * and undefined in the interior of the cells.
+ *
+ * This finite element is the trace space of FE_RavirtThomas on the
+ * faces and serves in hybridized methods.
+ *
+ * @author Guido Kanschat, 2009
+ */
+template <int dim, int spacedim=dim>
+class FE_FaceQ : FE_PolyFace<TensorProductPolynomials<dim-1>, dim, spacedim>
+{
+  public:
+                                    /**
+                                     * Constructor for tensor product
+                                     * polynomials of degree
+                                     * <tt>p</tt>. The shape
+                                     * functions created using this
+                                     * constructor correspond to
+                                     * Legendre polynomials in each
+                                     * coordinate direction.
+                                     */
+    FE_FaceQ(unsigned int p);
+    
+                                    /**
+                                     * Return a string that uniquely
+                                     * identifies a finite
+                                     * element. This class returns
+                                     * <tt>FE_DGQ<dim>(degree)</tt> , with
+                                     * <tt>dim</tt> and <tt>degree</tt>
+                                     * replaced by appropriate
+                                     * values.
+                                     */
+    virtual std::string get_name () 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.
+                                     *
+                                     * Implementation of the
+                                     * interface in
+                                     * FiniteElement
+                                     */
+    virtual bool has_support_on_face (const unsigned int shape_index,
+                                     const unsigned int face_index) const;
+
+  private:
+                                    /**
+                                     * Return vector with dofs per
+                                     * vertex, line, quad, hex.
+                                     */
+    static std::vector<unsigned int> get_dpo_vector (const unsigned int deg);
+};
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index ea3d772b965ab918d340c0ddcb3075f3afd2a8c7..cc0875e653e409c822874d1832561bc8dc52332f 100644 (file)
@@ -26,10 +26,12 @@ DEAL_II_NAMESPACE_OPEN
  * FiniteElement classes based on a polynomial spaces like the
  * TensorProductPolynomials or a PolynomialSpace classes.
  *
- * Every class that implements following functions can be used as
+ * Every class conforming to the following interface can be used as
  * template parameter POLY.
  *
  * @code
+ * static const unsigned int dimension;
+ *
  * double compute_value (const unsigned int i,
  *                       const Point<dim> &p) const;
  *
@@ -55,7 +57,10 @@ DEAL_II_NAMESPACE_OPEN
  * functions depend on the cells in real space, the update_once() and
  * update_each() functions must be overloaded.
  *
- * @author Ralf Hartmann 2004
+ * @todo Since nearly all functions for spacedim != dim are
+ * specialized, this class needs cleaning up.
+ *
+ * @author Ralf Hartmann 2004, Guido Kanschat, 2009
  **/
 template <class POLY, int dim=POLY::dimension, int spacedim=dim>
 class FE_Poly : public FiniteElement<dim,spacedim>
diff --git a/deal.II/deal.II/include/fe/fe_poly_face.h b/deal.II/deal.II/include/fe/fe_poly_face.h
new file mode 100644 (file)
index 0000000..5f5e96e
--- /dev/null
@@ -0,0 +1,387 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2009 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_face_h
+#define __deal2__fe_poly_face_h
+
+
+#include <fe/fe.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+/*!@addtogroup febase */
+/*@{*/
+
+/**
+ * @warning This class is not sufficiently tested yet!
+ *
+ * This class gives a unified framework for the implementation of
+ * FiniteElement classes only located on faces of the mesh. Theu are
+ * based on polynomial spaces like the TensorProductPolynomials or a
+ * PolynomialSpace classes.
+ *
+ * Every class that implements following functions can be used as
+ * template parameter POLY.
+ *
+ * @code
+ * double compute_value (const unsigned int i,
+ *                       const Point<dim> &p) const;
+ * @endcode
+ * Example classes are TensorProductPolynomials, PolynomialSpace or
+ * PolynomialsP.
+ *
+ * This class is not a fully implemented FiniteElement class. Instead
+ * there are several pure virtual functions declared in the
+ * FiniteElement and FiniteElement classes which cannot
+ * implemented by this class but are left for implementation in
+ * derived classes.
+ *
+ * Furthermore, this class assumes that shape functions of the
+ * FiniteElement under consideration do <em>not</em> depend on the
+ * actual shape of the cells in real space, i.e. update_once()
+ * includes <tt>update_values</tt>. For FiniteElements whose shape
+ * functions depend on the cells in real space, the update_once() and
+ * update_each() functions must be overloaded.
+ *
+ * @author Guido Kanschat, 2009
+ **/
+template <class POLY, int dim=POLY::dimension+1, int spacedim=dim>
+class FE_PolyFace : public FiniteElement<dim,spacedim>
+{
+  public:
+                                    /**
+                                     * Constructor.
+                                     */
+    FE_PolyFace (const POLY& poly_space,
+                const FiniteElementData<dim> &fe_data,
+                const std::vector<bool> &restriction_is_additive_flags);
+
+                                    /**
+                                     * Return the polynomial degree
+                                     * of this finite element,
+                                     * i.e. the value passed to the
+                                     * constructor.
+                                     */
+    unsigned int get_degree () const;
+
+                                    /**
+                                     * Return the value of the
+                                     * <tt>i</tt>th shape function at
+                                     * the point <tt>p</tt>. See the
+                                     * FiniteElement 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
+                                     * <tt>component</tt>th vector
+                                     * component of the <tt>i</tt>th
+                                     * shape function at the point
+                                     * <tt>p</tt>. See the
+                                     * FiniteElement 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
+                                     * <tt>_component</tt> 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
+                                     * <tt>i</tt>th shape function at
+                                     * the point <tt>p</tt>. See the
+                                     * FiniteElement 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
+                                     * <tt>component</tt>th vector
+                                     * component of the <tt>i</tt>th
+                                     * shape function at the point
+                                     * <tt>p</tt>. See the
+                                     * FiniteElement 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
+                                     * <tt>_component</tt> 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
+                                     * <tt>i</tt>th shape function at
+                                     * point <tt>p</tt> on the unit
+                                     * cell. See the
+                                     * FiniteElement 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 <tt>component</tt>th
+                                     * vector component of the
+                                     * <tt>i</tt>th shape function at
+                                     * the point <tt>p</tt>. See the
+                                     * FiniteElement 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
+                                     * <tt>_component</tt> 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,
+                                     * <tt>base_element(0)</tt> is
+                                     * <tt>this</tt>, and all other
+                                     * indices throw an error.
+                                     */
+    virtual const FiniteElement<dim,spacedim> &
+    base_element (const unsigned int index) const;
+
+                                     /**
+                                      * Multiplicity of base element
+                                      * <tt>index</tt>. Since this is
+                                      * a scalar element,
+                                      * <tt>element_multiplicity(0)</tt>
+                                      * returns one, and all other
+                                      * indices will throw an error.
+                                      */
+    virtual unsigned int element_multiplicity (const unsigned int index) const;
+
+    
+  protected:
+      
+    virtual
+    typename Mapping<dim,spacedim>::InternalDataBase *
+    get_data (const UpdateFlags,
+             const Mapping<dim,spacedim>& mapping,
+             const Quadrature<dim>& quadrature) const ;
+
+    typename Mapping<dim,spacedim>::InternalDataBase *
+    get_face_data (const UpdateFlags,
+                  const Mapping<dim,spacedim>& mapping,
+                  const Quadrature<dim-1>& quadrature) const ;
+    
+    typename Mapping<dim,spacedim>::InternalDataBase *
+    get_subface_data (const UpdateFlags,
+                     const Mapping<dim,spacedim>& mapping,
+                     const Quadrature<dim-1>& quadrature) const ;
+    
+    virtual void
+    fill_fe_values (const Mapping<dim,spacedim>                           &mapping,
+                   const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+                   const Quadrature<dim>                                 &quadrature,
+                   typename Mapping<dim,spacedim>::InternalDataBase      &mapping_internal,
+                   typename Mapping<dim,spacedim>::InternalDataBase      &fe_internal,
+                   FEValuesData<dim,spacedim>                            &data,
+                   CellSimilarity::Similarity                       &cell_similarity) const;
+    
+    virtual void
+    fill_fe_face_values (const Mapping<dim,spacedim> &mapping,
+                        const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+                        const unsigned int                    face_no,
+                        const Quadrature<dim-1>                &quadrature,
+                        typename Mapping<dim,spacedim>::InternalDataBase      &mapping_internal,
+                        typename Mapping<dim,spacedim>::InternalDataBase      &fe_internal,
+                        FEValuesData<dim,spacedim>& data) const ;
+    
+    virtual void
+    fill_fe_subface_values (const Mapping<dim,spacedim> &mapping,
+                           const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+                           const unsigned int                    face_no,
+                           const unsigned int                    sub_no,
+                           const Quadrature<dim-1>                &quadrature,
+                           typename Mapping<dim,spacedim>::InternalDataBase      &mapping_internal,
+                           typename Mapping<dim,spacedim>::InternalDataBase      &fe_internal,
+                           FEValuesData<dim,spacedim>& data) const ;
+
+    
+                                    /**
+                                     * Determine the values that need
+                                     * to be computed on the unit
+                                     * cell to be able to compute all
+                                     * values required by
+                                     * <tt>flags</tt>.
+                                     *
+                                     * For the purpuse of this
+                                     * function, refer to the
+                                     * documentation in
+                                     * FiniteElement.
+                                     *
+                                     * This class assumes that shape
+                                     * functions of this
+                                     * FiniteElement do <em>not</em>
+                                     * depend on the actual shape of
+                                     * the cells in real
+                                     * space. Therefore, the effect
+                                     * in this element is as follows:
+                                     * if <tt>update_values</tt> is
+                                     * set in <tt>flags</tt>, 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
+                                     * <tt>flags</tt>.
+                                     *
+                                     * For the purpuse of this
+                                     * function, refer to the
+                                     * documentation in
+                                     * FiniteElement.
+                                     *
+                                     * This class assumes that shape
+                                     * functions of this
+                                     * FiniteElement do <em>not</em>
+                                     * depend on the actual shape of
+                                     * the cells in real
+                                     * space.
+                                     *
+                                     * The effect in this element is
+                                     * as follows:
+                                     * <ul>
+
+                                     * <li> if
+                                     * <tt>update_gradients</tt> is
+                                     * set, the result will contain
+                                     * <tt>update_gradients</tt> and
+                                     * <tt>update_covariant_transformation</tt>.
+                                     * The latter is required to
+                                     * transform the gradient on the
+                                     * unit cell to the real
+                                     * cell. Remark, that the action
+                                     * required by
+                                     * <tt>update_covariant_transformation</tt>
+                                     * is actually performed by the
+                                     * Mapping object used in
+                                     * conjunction with this finite
+                                     * element.  <li> if
+                                     * <tt>update_hessians</tt>
+                                     * is set, the result will
+                                     * contain
+                                     * <tt>update_hessians</tt>
+                                     * and
+                                     * <tt>update_covariant_transformation</tt>.
+                                     * 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.
+                                     *
+                                     * </ul>
+                                     */
+    virtual UpdateFlags update_each (const UpdateFlags flags) const;
+
+
+                                    /**
+                                     * Fields of cell-independent data.
+                                     *
+                                     * For information about the
+                                     * general purpose of this class,
+                                     * see the documentation of the
+                                     * base class.
+                                     */
+    class InternalData : public FiniteElement<dim,spacedim>::InternalDataBase
+    {
+      public:
+                                        /**
+                                         * Array with shape function
+                                         * values in quadrature
+                                         * points on one face. 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 one face of 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.
+                                         *
+                                         * In particular, we can
+                                         * simply copy the same set
+                                         * of values to each of the
+                                         * faces.
+                                         */
+       std::vector<std::vector<double> > shape_values;
+    };
+    
+                                     /**
+                                      * The polynomial space. Its type
+                                      * is given by the template
+                                      * parameter POLY.
+                                      */    
+    POLY poly_space;
+};
+
+/*@}*/
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
diff --git a/deal.II/deal.II/include/fe/fe_poly_face.templates.h b/deal.II/deal.II/include/fe/fe_poly_face.templates.h
new file mode 100644 (file)
index 0000000..16aaee2
--- /dev/null
@@ -0,0 +1,281 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2009 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/qprojector.h>
+#include <base/polynomial_space.h>
+#include <fe/fe_values.h>
+#include <fe/fe_poly_face.h>
+
+
+DEAL_II_NAMESPACE_OPEN
+
+template <class POLY, int dim, int spacedim>
+FE_PolyFace<POLY,dim,spacedim>::FE_PolyFace (
+  const POLY& poly_space,
+  const FiniteElementData<dim> &fe_data,
+  const std::vector<bool> &restriction_is_additive_flags):
+               FiniteElement<dim,spacedim> (fe_data,
+                                            restriction_is_additive_flags,
+                                            std::vector<std::vector<bool> > (1, std::vector<bool>(1,true))),
+               poly_space(poly_space)
+{
+  AssertDimension(dim, POLY::dimension+1);
+}
+
+
+template <class POLY, int dim, int spacedim>
+unsigned int
+FE_PolyFace<POLY,dim,spacedim>::get_degree () const
+{
+  return this->degree;
+}
+
+
+//---------------------------------------------------------------------------
+// Auxiliary functions
+//---------------------------------------------------------------------------
+
+
+
+
+template <class POLY, int dim, int spacedim>
+UpdateFlags
+FE_PolyFace<POLY,dim,spacedim>::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, int spacedim>
+UpdateFlags
+FE_PolyFace<POLY,dim,spacedim>::update_each (const UpdateFlags flags) const
+{
+  UpdateFlags out = update_default;
+
+  if (flags & update_gradients)
+    out |= update_gradients | update_covariant_transformation;
+  if (flags & update_hessians)
+    out |= update_hessians | update_covariant_transformation;
+  if (flags & update_cell_normal_vectors)
+    out |= update_cell_normal_vectors | update_JxW_values;
+  
+  return out;
+}
+
+
+
+//---------------------------------------------------------------------------
+// Data field initialization
+//---------------------------------------------------------------------------
+
+template <class POLY, int dim, int spacedim>
+typename Mapping<dim,spacedim>::InternalDataBase *
+FE_PolyFace<POLY,dim,spacedim>::get_data (
+  const UpdateFlags,
+  const Mapping<dim,spacedim>&,
+  const Quadrature<dim>&) const
+{
+  Assert(false, ExcNotImplemented());
+  return 0;
+}
+
+
+template <class POLY, int dim, int spacedim>
+typename Mapping<dim,spacedim>::InternalDataBase *
+FE_PolyFace<POLY,dim,spacedim>::get_face_data (
+  const UpdateFlags update_flags,
+  const Mapping<dim,spacedim>&,
+  const Quadrature<dim-1>& 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.size();
+
+                                  // some scratch arrays
+  std::vector<double> values(0);
+  std::vector<Tensor<1,dim-1> > grads(0);
+  std::vector<Tensor<2,dim-1> > 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.resize (this->dofs_per_cell,
+                                std::vector<double> (n_q_points));
+      for (unsigned int i=0; i<n_q_points; ++i)
+       {
+         poly_space.compute(quadrature.point(i),
+                            values, grads, grad_grads);
+         
+         for (unsigned int k=0; k<this->dofs_per_cell; ++k)
+           data->shape_values[k][i] = values[k];
+       }
+    }
+                                  // No derivatives of this element
+                                  // are implemented.
+  if (flags & update_gradients || flags & update_hessians)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  
+  return data;
+}
+
+
+
+template <class POLY, int dim, int spacedim>
+typename Mapping<dim,spacedim>::InternalDataBase *
+FE_PolyFace<POLY,dim,spacedim>::get_subface_data (
+  const UpdateFlags flags,
+  const Mapping<dim,spacedim>& mapping,
+  const Quadrature<dim-1>& quadrature) const
+{
+  return get_face_data (flags, mapping,
+                       QProjector<dim-1>::project_to_all_children(quadrature));
+}
+
+
+
+
+
+//---------------------------------------------------------------------------
+// Fill data of FEValues
+//---------------------------------------------------------------------------
+template <class POLY, int dim, int spacedim>
+void
+FE_PolyFace<POLY,dim,spacedim>::fill_fe_values 
+  (const Mapping<dim,spacedim>&,
+   const typename Triangulation<dim,spacedim>::cell_iterator&,
+   const Quadrature<dim>&,
+   typename Mapping<dim,spacedim>::InternalDataBase&,
+   typename Mapping<dim,spacedim>::InternalDataBase&,
+   FEValuesData<dim,spacedim>&,
+   CellSimilarity::Similarity&) const
+{
+  Assert(false, ExcNotImplemented());
+}
+
+
+
+template <class POLY, int dim, int spacedim>
+void
+FE_PolyFace<POLY,dim,spacedim>::fill_fe_face_values (
+  const Mapping<dim,spacedim>&,
+  const typename Triangulation<dim,spacedim>::cell_iterator&,
+  const unsigned int,
+  const Quadrature<dim-1>& quadrature,
+  typename Mapping<dim,spacedim>::InternalDataBase&,
+  typename Mapping<dim,spacedim>::InternalDataBase& fedata,
+  FEValuesData<dim,spacedim>& data) const
+{
+                                  // convert data object to internal
+                                  // data for this class. fails with
+                                  // an exception if that is not
+                                  // possible
+  Assert (dynamic_cast<InternalData *> (&fedata) != 0, ExcInternalError());
+  InternalData &fe_data = static_cast<InternalData &> (fedata);
+  
+  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.size(); ++i)
+         data.shape_values(k,i) = fe_data.shape_values[k][i];
+    }
+}
+
+
+template <class POLY, int dim, int spacedim>
+void
+FE_PolyFace<POLY,dim,spacedim>::fill_fe_subface_values (
+  const Mapping<dim,spacedim>&,
+  const typename Triangulation<dim,spacedim>::cell_iterator&,
+  const unsigned int,
+  const unsigned int subface,
+  const Quadrature<dim-1>& quadrature,
+  typename Mapping<dim,spacedim>::InternalDataBase&,
+  typename Mapping<dim,spacedim>::InternalDataBase& fedata,
+  FEValuesData<dim,spacedim>& data) const
+{
+                                  // convert data object to internal
+                                  // data for this class. fails with
+                                  // an exception if that is not
+                                  // possible
+  Assert (dynamic_cast<InternalData *> (&fedata) != 0, ExcInternalError());
+  InternalData &fe_data = static_cast<InternalData &> (fedata);
+
+  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
+
+  unsigned int offset = subface*quadrature.size();
+  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
+    {
+      if (flags & update_values)
+        for (unsigned int i=0; i<quadrature.size(); ++i)
+         data.shape_values(k,i) = fe_data.shape_values[k][i+offset];
+      
+    }
+  Assert (!(flags & update_gradients), ExcNotImplemented());
+  Assert (!(flags & update_hessians), ExcNotImplemented());
+}
+
+
+
+template <class POLY, int dim, int spacedim>
+unsigned int
+FE_PolyFace<POLY,dim,spacedim>::n_base_elements () const
+{
+  return 1;
+}
+
+
+
+template <class POLY, int dim, int spacedim>
+const FiniteElement<dim,spacedim> &
+FE_PolyFace<POLY,dim,spacedim>::base_element (const unsigned int index) const
+{
+  Assert (index==0, ExcIndexRange(index, 0, 1));
+  return *this;
+}
+
+
+
+template <class POLY, int dim, int spacedim>
+unsigned int
+FE_PolyFace<POLY,dim,spacedim>::element_multiplicity (const unsigned int index) const
+{
+  Assert (index==0, ExcIndexRange(index, 0, 1));
+  return 1;
+}
+
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/deal.II/deal.II/source/fe/fe_face.cc b/deal.II/deal.II/source/fe/fe_face.cc
new file mode 100644 (file)
index 0000000..43d966b
--- /dev/null
@@ -0,0 +1,78 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2002, 2003, 2004, 2005, 2006, 2007, 2008 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 <fe/fe_face.h>
+#include <fe/fe_poly_face.templates.h>
+#include <sstream>
+
+DEAL_II_NAMESPACE_OPEN
+
+template <int dim, int spacedim>
+FE_FaceQ<dim,spacedim>::FE_FaceQ (const unsigned int degree)
+               :
+                FE_PolyFace<TensorProductPolynomials<dim-1>, dim, spacedim> (
+                 TensorProductPolynomials<dim-1>(Polynomials::Legendre::generate_complete_basis(degree)),
+                 FiniteElementData<dim>(get_dpo_vector(degree), 1, degree, FiniteElementData<dim>::L2),
+                 std::vector<bool>(1,true))
+{}
+
+
+template <int dim, int spacedim>
+std::string
+FE_FaceQ<dim,spacedim>::get_name () const
+{
+                                  // note that the
+                                  // FETools::get_fe_from_name
+                                  // function depends on the
+                                  // particular format of the string
+                                  // this function returns, so they
+                                  // have to be kept in synch
+
+  std::ostringstream namebuf;  
+  namebuf << "FE_FaceQ<" << dim << ">(" << this->degree << ")";
+
+  return namebuf.str();
+}
+
+template <int dim, int spacedim>
+bool
+FE_FaceQ<dim,spacedim>::has_support_on_face (
+  const unsigned int shape_index,
+  const unsigned int face_index) const
+{
+  return (face_index == (shape_index/this->dofs_per_face));
+}
+                 
+                 
+template <int dim, int spacedim>
+std::vector<unsigned int>
+FE_FaceQ<dim,spacedim>::get_dpo_vector (const unsigned int deg)
+{
+  std::vector<unsigned int> dpo(dim+1, 0U);
+  dpo[dim-1] = deg+1;
+  for (unsigned int i=1;i<dim-1;++i)
+    dpo[dim-1] *= deg+1;
+  return dpo;
+}
+
+
+                 
+
+#if deal_II_dimension > 1
+                 template class FE_PolyFace<TensorProductPolynomials<deal_II_dimension-1> >;
+//template class FE_PolyFace<PolynomialSpace<deal_II_dimension>, deal_II_dimension>;
+//template class FE_PolyFace<PolynomialsP<deal_II_dimension>, deal_II_dimension>;
+                 template class FE_FaceQ<deal_II_dimension>;
+#endif
+                 
+DEAL_II_NAMESPACE_CLOSE

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.