]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
new class FE_DGPNonparametric
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 19 Sep 2002 22:21:54 +0000 (22:21 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 19 Sep 2002 22:21:54 +0000 (22:21 +0000)
git-svn-id: https://svn.dealii.org/trunk@6488 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_dgp_nonparametric.h [new file with mode: 0644]
deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc [new file with mode: 0644]
deal.II/doc/news/2002/c-3-4.html

diff --git a/deal.II/deal.II/include/fe/fe_dgp_nonparametric.h b/deal.II/deal.II/include/fe/fe_dgp_nonparametric.h
new file mode 100644 (file)
index 0000000..b755739
--- /dev/null
@@ -0,0 +1,457 @@
+//----------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002 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_dgp_nonparametric_h
+#define __deal2__fe_dgp_nonparametric_h
+
+#include <base/config.h>
+#include <base/polynomial.h>
+#include <base/polynomial_space.h>
+#include <fe/fe.h>
+
+template <int dim> class PolynomialSpace;
+template <int dim> class MappingQ;
+
+
+/**
+ * Discontinuous finite elements evaluated at the mapped quadrature points.
+ *
+ * This finite element implements complete polynomial spaces, that is,
+ * $d$-dimensional polynomials of order $k$.
+ *
+ * The polynomials are not mapped. Therefore, they are constant,
+ * linear, quadratic, etc. on any grid cell.
+ *
+ * Since the polynomials are evaluated at the quadrature points of the
+ * actual grid cell, no grid transfer and interpolation matrices are
+ * available.
+ *
+ * The purpose of this class is experimental, therefore the
+ * implementation will remain incomplete.
+ *
+ * @author Guido Kanschat, 2002
+ */
+template <int dim>
+class FE_DGPNonparametric : public FiniteElement<dim>
+{
+  public:
+                                    /**
+                                     * Constructor for tensor product
+                                     * polynomials of degree @p{k}.
+                                     */
+    FE_DGPNonparametric (const unsigned int k);
+    
+                                    /**
+                                     * 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;
+
+                                    /**
+                                     * Return the polynomial degree
+                                     * of this finite element,
+                                     * i.e. the value passed to the
+                                     * constructor.
+                                     */
+    unsigned int get_degree () 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;
+    
+                                    /**
+                                     * 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
+                                     * @ref{FiniteElement}
+                                     */
+    virtual bool has_support_on_face (const unsigned int shape_index,
+                                     const unsigned int face_index) const;
+
+                                    /**
+                                     * Determine an estimate for the
+                                     * memory consumption (in bytes)
+                                     * of this object.
+                                     *
+                                     * This function is made virtual,
+                                     * since finite element objects
+                                     * are usually accessed through
+                                     * pointers to their base class,
+                                     * rather than the class itself.
+                                     */
+    virtual unsigned int memory_consumption () const;
+
+
+                                    /**
+                                     * Declare a nested class which
+                                     * will hold static definitions of
+                                     * various matrices such as
+                                     * constraint and embedding
+                                     * matrices. The definition of
+                                     * the various static fields are
+                                     * in the files @p{fe_dgp_[123]d.cc}
+                                     * in the source directory.
+                                     */
+    struct Matrices
+    {
+                                        /**
+                                         * Pointers to the embedding
+                                         * matrices, one for each
+                                         * polynomial degree starting
+                                         * from constant elements
+                                         */
+       static const double * const embedding[][GeometryInfo<dim>::children_per_cell];
+
+                                        /**
+                                         * Number of elements (first
+                                         * index) the above field
+                                         * has. Equals the highest
+                                         * polynomial degree plus one
+                                         * for which the embedding
+                                         * matrices have been
+                                         * computed.
+                                         */
+       static const unsigned int n_embedding_matrices;
+
+                                        /**
+                                         * As @p{embedding} but for
+                                         * projection matrices.
+                                         */
+       static const double * const projection_matrices[][GeometryInfo<dim>::children_per_cell];
+
+                                        /**
+                                         * As
+                                         * @p{n_embedding_matrices}
+                                         * but for projection
+                                         * matrices.
+                                         */
+       static const unsigned int n_projection_matrices;
+    };
+
+  protected:
+
+                                    /**
+                                     * @p{clone} function instead of
+                                     * a copy constructor.
+                                     *
+                                     * This function is needed by the
+                                     * constructors of @p{FESystem}.
+                                     */
+    virtual FiniteElement<dim> *clone() const;
+  
+                                    /**
+                                     * 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 ;
+
+  private:
+    
+                                    /**
+                                     * Only for internal use. Its
+                                     * full name is
+                                     * @p{get_dofs_per_object_vector}
+                                     * function and it creates the
+                                     * @p{dofs_per_object} vector that is
+                                     * needed within the constructor to
+                                     * be passed to the constructor of
+                                     * @p{FiniteElementData}.
+                                     */
+    static std::vector<unsigned int> get_dpo_vector(unsigned int degree);
+    
+                                    /**
+                                     * Given a set of flags indicating
+                                     * what quantities are requested
+                                     * from a @p{FEValues} object,
+                                     * return which of these can be
+                                     * precomputed once and for
+                                     * all. Often, the values of
+                                     * shape function at quadrature
+                                     * points can be precomputed, for
+                                     * example, in which case the
+                                     * return value of this function
+                                     * would be the logical and of
+                                     * the input @p{flags} and
+                                     * @p{update_values}.
+                                     *
+                                     * For the present kind of finite
+                                     * element, this is exactly the
+                                     * case.
+                                     */
+    virtual UpdateFlags update_once (const UpdateFlags flags) const;
+  
+                                    /**
+                                     * This is the opposite to the
+                                     * above function: given a set of
+                                     * flags indicating what we want
+                                     * to know, return which of these
+                                     * need to be computed each time
+                                     * we visit a new cell.
+                                     *
+                                     * If for the computation of one
+                                     * quantity something else is
+                                     * also required (for example, we
+                                     * often need the covariant
+                                     * transformation when gradients
+                                     * need to be computed), include
+                                     * this in the result as well.
+                                     */
+    virtual UpdateFlags update_each (const UpdateFlags flags) const;
+  
+                                    /**
+                                     * Degree of the polynomials.
+                                     */  
+    const unsigned int degree;
+
+                                    /**
+                                     * Pointer to an object
+                                     * representing the polynomial
+                                     * space used here.
+                                     */
+    const PolynomialSpace<dim> polynomial_space;
+
+                                    /**
+                                     * 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:
+                                      // have some scratch arrays
+      std::vector<double> values;
+      std::vector<Tensor<1,dim> > grads;
+      std::vector<Tensor<2,dim> > grad_grads;
+    };
+    
+                                    /**
+                                     * Allow access from other dimensions.
+                                     */
+    template <int dim1> friend class FE_DGPNonparametric;
+
+                                    /**
+                                     * Allows @p{MappingQ} class to
+                                     * access to build_renumbering
+                                     * function.
+                                     */
+    friend class MappingQ<dim>;
+};
+
+
+// declaration of explicit specializations of member variables, if the
+// compiler allows us to do that (the standard says we must)
+#ifndef DEAL_II_MEMBER_VAR_SPECIALIZATION_BUG
+template <> 
+const double * const FE_DGPNonparametric<1>::Matrices::embedding[][GeometryInfo<1>::children_per_cell];
+
+template <>
+const unsigned int FE_DGPNonparametric<1>::Matrices::n_embedding_matrices;
+
+template <>
+const double * const FE_DGPNonparametric<1>::Matrices::projection_matrices[][GeometryInfo<1>::children_per_cell];
+
+template <>
+const unsigned int FE_DGPNonparametric<1>::Matrices::n_projection_matrices;
+
+template <> 
+const double * const FE_DGPNonparametric<2>::Matrices::embedding[][GeometryInfo<2>::children_per_cell];
+
+template <>
+const unsigned int FE_DGPNonparametric<2>::Matrices::n_embedding_matrices;
+
+template <>
+const double * const FE_DGPNonparametric<2>::Matrices::projection_matrices[][GeometryInfo<2>::children_per_cell];
+
+template <>
+const unsigned int FE_DGPNonparametric<2>::Matrices::n_projection_matrices;
+
+template <> 
+const double * const FE_DGPNonparametric<3>::Matrices::embedding[][GeometryInfo<3>::children_per_cell];
+
+template <>
+const unsigned int FE_DGPNonparametric<3>::Matrices::n_embedding_matrices;
+
+template <>
+const double * const FE_DGPNonparametric<3>::Matrices::projection_matrices[][GeometryInfo<3>::children_per_cell];
+
+template <>
+const unsigned int FE_DGPNonparametric<3>::Matrices::n_projection_matrices;
+#endif
+
+#endif
diff --git a/deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc b/deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc
new file mode 100644 (file)
index 0000000..0d18b01
--- /dev/null
@@ -0,0 +1,439 @@
+//----------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002 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 <grid/tria.h>
+#include <grid/tria_iterator.h>
+#include <dofs/dof_accessor.h>
+#include <fe/fe.h>
+#include <fe/mapping.h>
+#include <fe/fe_dgp_nonparametric.h>
+#include <fe/fe_values.h>
+
+
+
+template <int dim>
+FE_DGPNonparametric<dim>::FE_DGPNonparametric (const unsigned int degree)
+               :
+               FiniteElement<dim> (FiniteElementData<dim>(get_dpo_vector(degree),1),
+                                   std::vector<bool>(FiniteElementData<dim>(get_dpo_vector(degree),1).dofs_per_cell,true),
+                                   std::vector<std::vector<bool> >(FiniteElementData<dim>(get_dpo_vector(degree),1).dofs_per_cell,
+                                                                   std::vector<bool>(1,true))),
+                degree(degree),
+                polynomial_space (Legendre<double>::generate_complete_basis(degree))
+{
+                                  // if defined, copy over matrices
+                                  // from precomputed arrays
+//    if ((degree < Matrices::n_embedding_matrices) &&
+//        (Matrices::embedding[degree][0] != 0))
+//      for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
+//        this->prolongation[c].fill (Matrices::embedding[degree][c]);
+//    else
+//      for (unsigned int i=0; i<GeometryInfo<dim>::children_per_cell;++i)
+//        this->prolongation[i].reinit(0,0);
+
+                                   // restriction can be defined
+                                   // through projection for
+                                   // discontinuous elements, but is
+                                   // presently not implemented for DGPNonparametric
+                                   // elements.
+                                   //
+                                   // if it were, then the following
+                                   // snippet would be the right code
+//    if ((degree < Matrices::n_projection_matrices) &&
+//        (Matrices::projection_matrices[degree] != 0))
+//      {
+//        restriction[0].fill (Matrices::projection_matrices[degree]);
+//      }
+//    else
+//                                  // matrix undefined, set size to zero
+//      for (unsigned int i=0;i<GeometryInfo<dim>::children_per_cell;++i)
+//        restriction[i].reinit(0, 0);
+                                   // since not implemented, set to
+                                   // "empty"
+  for (unsigned int i=0;i<GeometryInfo<dim>::children_per_cell;++i)
+    restriction[i].reinit(0, 0);
+
+                                   // note further, that these
+                                   // elements have neither support
+                                   // nor face-support points, so
+                                   // leave these fields empty
+};
+
+
+
+template <int dim>
+FiniteElement<dim> *
+FE_DGPNonparametric<dim>::clone() const
+{
+  return new FE_DGPNonparametric<dim>(degree);
+}
+
+
+
+template <int dim>
+double
+FE_DGPNonparametric<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 polynomial_space.compute_value(i, p);
+}
+
+
+
+template <int dim>
+double
+FE_DGPNonparametric<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 polynomial_space.compute_value(i, p);
+}
+
+
+
+template <int dim>
+Tensor<1,dim>
+FE_DGPNonparametric<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 polynomial_space.compute_grad(i, p);
+}
+
+
+template <int dim>
+Tensor<1,dim>
+FE_DGPNonparametric<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 polynomial_space.compute_grad(i, p);
+}
+
+
+
+template <int dim>
+Tensor<2,dim>
+FE_DGPNonparametric<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 polynomial_space.compute_grad_grad(i, p);
+}
+
+
+
+template <int dim>
+Tensor<2,dim>
+FE_DGPNonparametric<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 polynomial_space.compute_grad_grad(i, p);
+}
+
+
+//----------------------------------------------------------------------
+// Auxiliary functions
+//----------------------------------------------------------------------
+
+
+template <int dim>
+std::vector<unsigned int>
+FE_DGPNonparametric<dim>::get_dpo_vector(unsigned int deg)
+{
+  std::vector<unsigned int> dpo(dim+1, static_cast<unsigned int>(0));
+  dpo[dim] = ++deg;
+  for (unsigned int i=1;i<dim;++i)
+    {
+      dpo[dim] *= deg+i;
+      dpo[dim] /= i+1;
+    }
+  return dpo;
+}
+
+
+template <int dim>
+UpdateFlags
+FE_DGPNonparametric<dim>::update_once (const UpdateFlags) 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;
+}
+
+
+template <int dim>
+UpdateFlags
+FE_DGPNonparametric<dim>::update_each (const UpdateFlags flags) const
+{
+  UpdateFlags out = flags & update_values;
+
+  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 <int dim>
+typename Mapping<dim>::InternalDataBase *
+FE_DGPNonparametric<dim>::get_data (const UpdateFlags      update_flags,
+                      const Mapping<dim>    &mapping,
+                      const Quadrature<dim> &quadrature) const
+{
+                                  // generate a new data object
+  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);
+
+                                  // initialize fields only if really
+                                  // necessary. otherwise, don't
+                                  // allocate memory
+  if (flags & update_values)
+    {
+      data->values.resize (this->dofs_per_cell);
+    }
+
+  if (flags & update_gradients)
+    {
+      data->grads.resize (this->dofs_per_cell);
+    }
+
+                                  // 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);
+  return data;
+}
+
+
+
+//----------------------------------------------------------------------
+// Fill data of FEValues
+//----------------------------------------------------------------------
+
+template <int dim>
+void
+FE_DGPNonparametric<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());
+
+  const unsigned int n_q_points = data.quadrature_points.size();
+  
+  if (flags & (update_values | update_gradients))
+    for (unsigned int i=0; i<n_q_points; ++i)
+      {
+       polynomial_space.compute(data.quadrature_points[i],
+                                fe_data.values, fe_data.grads, fe_data.grad_grads);
+       for (unsigned int k=0; k<this->dofs_per_cell; ++k)
+         {
+           if (flags & update_values)
+             data.shape_values[k][i] = fe_data.values[k];
+           if (flags & update_gradients)
+             data.shape_gradients[k][i] = fe_data.grads[k];
+         }
+      }
+  
+  if (flags & update_second_derivatives)
+    compute_2nd (mapping, cell, 0, mapping_data, fe_data, data);
+  
+  fe_data.first_cell = false;
+}
+
+
+
+template <int dim>
+void
+FE_DGPNonparametric<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);
+
+  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
+
+  const unsigned int n_q_points = data.quadrature_points.size();
+  
+  if (flags & (update_values | update_gradients))
+    for (unsigned int i=0; i<n_q_points; ++i)
+      {
+       polynomial_space.compute(data.quadrature_points[i],
+                                fe_data.values, fe_data.grads, fe_data.grad_grads);
+       for (unsigned int k=0; k<this->dofs_per_cell; ++k)
+         {
+           if (flags & update_values)
+             data.shape_values[k][i] = fe_data.values[k];
+           if (flags & update_gradients)
+             data.shape_gradients[k][i] = fe_data.grads[k];
+         }
+      }
+
+                                  // offset determines which data set
+                                  // to take (all data sets for all
+                                  // faces are stored contiguously)
+  const unsigned int offset = face * quadrature.n_quadrature_points;
+  
+  if (flags & update_second_derivatives)
+    compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
+  
+  fe_data.first_cell = false;
+}
+
+
+
+template <int dim>
+void
+FE_DGPNonparametric<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);
+  
+  const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
+
+  const unsigned int n_q_points = data.quadrature_points.size();
+  
+  if (flags & (update_values | update_gradients))
+    for (unsigned int i=0; i<n_q_points; ++i)
+      {
+       polynomial_space.compute(data.quadrature_points[i],
+                                fe_data.values, fe_data.grads, fe_data.grad_grads);
+       for (unsigned int k=0; k<this->dofs_per_cell; ++k)
+         {
+           if (flags & update_values)
+             data.shape_values[k][i] = fe_data.values[k];
+           if (flags & update_gradients)
+             data.shape_gradients[k][i] = fe_data.grads[k];
+         }
+      }
+  
+                                  // offset determines which data set
+                                  // to take (all data sets for all
+                                  // sub-faces are stored contiguously)
+  const unsigned int offset = (face * GeometryInfo<dim>::subfaces_per_face + subface)
+                             * quadrature.n_quadrature_points;
+
+  if (flags & update_second_derivatives)
+    compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
+  
+  fe_data.first_cell = false;
+}
+
+
+
+template <int dim>
+unsigned int
+FE_DGPNonparametric<dim>::n_base_elements () const
+{
+  return 1;
+};
+
+
+
+template <int dim>
+const FiniteElement<dim> &
+FE_DGPNonparametric<dim>::base_element (const unsigned int index) const
+{
+  Assert (index==0, ExcIndexRange(index, 0, 1));
+  return *this;
+};
+
+
+
+template <int dim>
+bool
+FE_DGPNonparametric<dim>::has_support_on_face (const unsigned int,
+                                 const unsigned int) const
+{
+  return true;
+}
+
+
+
+template <int dim>
+unsigned int
+FE_DGPNonparametric<dim>::memory_consumption () const
+{
+  Assert (false, ExcNotImplemented ());
+  return 0;
+}
+
+
+
+template <int dim>
+unsigned int
+FE_DGPNonparametric<dim>::get_degree () const
+{
+  return degree;
+};
+
+
+
+
+template class FE_DGPNonparametric<deal_II_dimension>;
index 4c0927a9e04e901b669552b498ca55f187399e89..b6b2595ab4c9e33257d52d144862cc4db58a8827 100644 (file)
@@ -299,6 +299,16 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
 <h3>deal.II</h3>
 
 <ol>
+  <li> <p> New: the class <code
+  class="class">FE_DGPNonparametric</code> implements finite elements
+  where shape functions are polynomials of order <i>k</i> on the
+  actual grid cell. This is achieved by evaluating the polynomials at
+  the mapped quadrature points. No grid transfer matrices are
+  available for this class.
+  <br>
+  (GK 09/19/2002)
+  </p>
+
   <li> <p> 
        Fixed: Some of the various instances of the <code
        class="member">VectorTools::interpolate_boundary_values</code> functions

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.