--- /dev/null
+//----------------------------------------------------------------
+// $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
--- /dev/null
+//----------------------------------------------------------------
+// $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>;