--- /dev/null
+//---------------------------------------------------------------------------
+// $Id$
+// Version: $Name$
+//
+// Copyright (C) 2005 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_tensor_h
+#define __deal2__fe_poly_tensor_h
+
+
+#include <lac/full_matrix.h>
+#include <fe/fe.h>
+
+/*!@addtogroup fe */
+/*@{*/
+
+/**
+ * This class gives a unified framework for the implementation of
+ * FiniteElement classes based on Tensor valued polynomial spaces like
+ * PolynomialsBDM and PolynomialsRaviartThomas.
+ *
+ * Every class that implements following function can be used as
+ * template parameter POLY.
+ *
+ * @code
+ * void compute (const Point<dim> &unit_point,
+ * std::vector<Tensor<1,dim> > &values,
+ * std::vector<Tensor<2,dim> > &grads,
+ * std::vector<Tensor<3,dim> > &grad_grads) const;
+ * @endcode
+ *
+ * The polynomial spaces are usually described as direct sum of
+ * simpler spaces. In such a case, the usual basis of node functionals
+ * is not dual to the basis of the polynomial space. Therefore, the
+ * matrix #inverse_node_matrix can be filled by the constructor of a
+ * derived class such that the usual interpolation condition
+ * <i>N<sub>i</sub>(v<sub>j</sub>)</i> holds on the reference cell.
+ *
+ * In many cases, the node functionals depend on the shape of the mesh
+ * cell, since they evaluate normal or tangential components on the
+ * faces. In order to allow for a set of transformations, the variable
+ * #mapping_type has been introduced. It should also be set in the
+ * constructor of a derived class.
+ *
+ * This class is not a fully implemented FiniteElement class, but
+ * implements some common features of vector valued elements based on
+ * vector valued polynomial classes. What's missing here in particular
+ * is information on the topological location of the node values.
+ *
+ * @see PolynomialsBDM, PolynomialsRaviartThomas
+ *
+ * @author Guido Kanschat, 2005
+ **/
+template <class POLY, int dim>
+class FE_PolyTensor : public FiniteElement<dim>
+{
+ public:
+ /**
+ * Constructor.
+ */
+ FE_PolyTensor (unsigned int degree,
+ const FiniteElementData<dim> &fe_data,
+ const std::vector<bool> &restriction_is_additive_flags,
+ const std::vector<std::vector<bool> > &nonzero_components);
+
+ /**
+ * Since these elements are
+ * vector valued, an exception is
+ * thrown.
+ */
+ virtual double shape_value (const unsigned int i,
+ const Point<dim> &p) const;
+
+ virtual double shape_value_component (const unsigned int i,
+ const Point<dim> &p,
+ const unsigned int component) const;
+
+ /**
+ * Since these elements are
+ * vector valued, an exception is
+ * thrown.
+ */
+ virtual Tensor<1,dim> shape_grad (const unsigned int i,
+ const Point<dim> &p) const;
+
+ virtual Tensor<1,dim> shape_grad_component (const unsigned int i,
+ const Point<dim> &p,
+ const unsigned int component) const;
+
+ /**
+ * Since these elements are
+ * vector valued, an exception is
+ * thrown.
+ */
+ virtual Tensor<2,dim> shape_grad_grad (const unsigned int i,
+ const Point<dim> &p) const;
+
+ 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 not a composed element,
+ * return one.
+ */
+ virtual unsigned int n_base_elements () const;
+
+ /**
+ * Access to base element
+ * objects. Since this element is
+ * not composed of several elements,
+ * <tt>base_element(0)</tt> is
+ * <tt>this</tt>, and all other
+ * indices throw an error.
+ */
+ virtual const FiniteElement<dim> &
+ base_element (const unsigned int index) const;
+
+ /**
+ * Multiplicity of base element
+ * <tt>index</tt>. Since this is
+ * not a composed 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:
+ /**
+ * Different options for
+ * transforming the basis
+ * functions from the reference
+ * cell to the actual mesh cell.
+ *
+ * Most vector valued elements
+ * either transform shape
+ * functions to keep node values
+ * on edges meaningful. Still, in
+ * special cases, it may be
+ * possible to avoid the mapping.
+ */
+ enum MappingType {
+ /// Shape functions do not depend on actual mesh cell
+ independent,
+ /// Shape functions are transformed covariant.
+ covariant,
+ /// Shape functions are transformed contravariant.
+ contravariant
+ };
+
+ /**
+ * The mapping type to be used to
+ * map shape functions from the
+ * reference cell to the msh
+ * cell.
+ */
+ MappingType mapping_type;
+
+ virtual
+ typename Mapping<dim>::InternalDataBase *
+ get_data (const UpdateFlags,
+ const Mapping<dim>& mapping,
+ const Quadrature<dim>& quadrature) const ;
+
+ virtual void
+ fill_fe_values (const Mapping<dim> &mapping,
+ const typename Triangulation<dim>::cell_iterator &cell,
+ const Quadrature<dim> &quadrature,
+ typename Mapping<dim>::InternalDataBase &mapping_internal,
+ typename Mapping<dim>::InternalDataBase &fe_internal,
+ FEValuesData<dim>& data) const;
+
+ virtual void
+ fill_fe_face_values (const Mapping<dim> &mapping,
+ const typename Triangulation<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 ;
+
+ virtual void
+ fill_fe_subface_values (const Mapping<dim> &mapping,
+ const typename Triangulation<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 ;
+
+ /**
+ * Fields of cell-independent
+ * data for FE_PolyTensor. Stores
+ * the values of the shape
+ * functions and their
+ * derivatives on the reference
+ * cell for later use.
+ *
+ * All tables are organized in a
+ * way, that the value for shape
+ * function <i>i</i> at
+ * quadrature point <i>k</i> is
+ * accessed by indices
+ * <i>(i,k)</i>.
+ */
+ class InternalData : public FiniteElementBase<dim>::InternalDataBase
+ {
+ public:
+ /**
+ * Array with shape function
+ * values in quadrature
+ * points. There is one
+ * row for each shape
+ * function, containing
+ * values for each quadrature
+ * point.
+ */
+ std::vector<std::vector<Tensor<1,dim> > > shape_values;
+
+ /**
+ * Array with shape function
+ * gradients in quadrature
+ * points. There is one
+ * row for each shape
+ * function, containing
+ * values for each quadrature
+ * point.
+ */
+ std::vector<std::vector<Tensor<2,dim> > > shape_grads;
+ };
+
+ /**
+ * Degree of the polynomials.
+ */
+ unsigned int degree;
+
+ /**
+ * The polynomial space. Its type
+ * is given by the template
+ * parameter POLY.
+ */
+ POLY poly_space;
+ /**
+ * The inverse of the matrix
+ * <i>a<sub>ij</sub></i> of node
+ * values <i>N<sub>i</sub></i>
+ * applied to polynomial
+ * <i>p<sub>j</sub></i>. This
+ * matrix is used to convert
+ * polynomials in the "raw" basis
+ * provided in #poly_space to the
+ * basis dual to the node
+ * functionals on the reference cell.
+ *
+ * This object is not filled by
+ * FE_PolyTensor, but is a chance
+ * for a derived class to allow
+ * for reorganization of the
+ * basis functions. If it is left
+ * empty, the basis in
+ * #poly_space is used.
+ */
+ FullMatrix<double> inverse_node_matrix;
+
+ /**
+ * If a shape function is
+ * computed at a single point, we
+ * must compute all of them to
+ * apply #inverse_node_matrix. In
+ * order to avoid too much
+ * overhead, we cache the point
+ * and the function values for
+ * the next evaluation.
+ */
+ mutable Point<dim> cached_point;
+
+ /**
+ * Cached shape function values after
+ * call to
+ * shape_value_component().
+ */
+ mutable std::vector<Tensor<1,dim> > cached_values;
+
+ /**
+ * Cached shape function gradients after
+ * call to
+ * shape_grad_component().
+ */
+ mutable std::vector<Tensor<2,dim> > cached_grads;
+
+ /**
+ * Cached second derivatives of
+ * shape functions after call to
+ * shape_grad_grad_component().
+ */
+ mutable std::vector<Tensor<3,dim> > cached_grad_grads;
+};
+
+/*@}*/
+
+#endif
--- /dev/null
+//---------------------------------------------------------------------------
+// $Id$
+// Version: $Name$
+//
+// Copyright (C) 2005 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/polynomials_bdm.h>
+#include <base/polynomials_raviart_thomas.h>
+#include <fe/fe_poly_tensor.h>
+#include <fe/fe_values.h>
+
+
+template <class POLY, int dim>
+FE_PolyTensor<POLY,dim>::FE_PolyTensor (
+ unsigned int degree,
+ const FiniteElementData<dim> &fe_data,
+ const std::vector<bool> &restriction_is_additive_flags,
+ const std::vector<std::vector<bool> > &nonzero_components):
+ FiniteElement<dim> (fe_data,
+ restriction_is_additive_flags,
+ nonzero_components),
+ degree(degree),
+ poly_space(POLY(degree))
+{
+ cached_point(0) = -1;
+}
+
+
+template <class POLY, int dim>
+double
+FE_PolyTensor<POLY,dim>::shape_value (
+ const unsigned int, const Point<dim> &) const
+{
+ Assert(false, typename FiniteElementBase<dim>::ExcFENotPrimitive());
+ return 0.;
+}
+
+
+template <class POLY, int dim>
+double
+FE_PolyTensor<POLY,dim>::shape_value_component (
+ const unsigned int i,
+ const Point<dim> &p,
+ const unsigned int component) const
+{
+ Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
+ Assert (component < dim, ExcIndexRange (component, 0, dim));
+
+ if (cached_point != p || cached_values.size() == 0)
+ {
+ cached_point = p;
+ cached_values.resize(poly_space.n());
+ poly_space.compute(p, cached_values, cached_grads, cached_grad_grads);
+ }
+
+ double s = 0;
+ if (inverse_node_matrix.n_cols() == 0)
+ return cached_values[i][component];
+ else
+ for (unsigned int j=0;j<inverse_node_matrix.n_cols();++j)
+ s += inverse_node_matrix(i,j) * cached_values[j][component];
+ return s;
+}
+
+
+
+template <class POLY, int dim>
+Tensor<1,dim>
+FE_PolyTensor<POLY,dim>::shape_grad (
+ const unsigned int, const Point<dim> &) const
+{
+ Assert(false, typename FiniteElementBase<dim>::ExcFENotPrimitive());
+ return Tensor<1,dim>();
+}
+
+
+
+template <class POLY, int dim>
+Tensor<1,dim>
+FE_PolyTensor<POLY,dim>::shape_grad_component (
+ const unsigned int i,
+ const Point<dim> &p,
+ const unsigned int component) const
+{
+ Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
+ Assert (component < dim, ExcIndexRange (component, 0, dim));
+
+ if (cached_point != p || cached_grads.size() == 0)
+ {
+ cached_point = p;
+ cached_grads.resize(poly_space.n());
+ poly_space.compute(p, cached_values, cached_grads, cached_grad_grads);
+ }
+
+ Tensor<1,dim> s;
+ if (inverse_node_matrix.n_cols() == 0)
+ return cached_grads[i][component];
+ else
+ for (unsigned int j=0;j<inverse_node_matrix.n_cols();++j)
+ s += inverse_node_matrix(i,j) * cached_grads[j][component];
+
+ return s;
+}
+
+
+
+template <class POLY, int dim>
+Tensor<2,dim>
+FE_PolyTensor<POLY,dim>::shape_grad_grad (
+ const unsigned int, const Point<dim> &) const
+{
+ Assert(false, typename FiniteElementBase<dim>::ExcFENotPrimitive());
+ return Tensor<2,dim>();
+}
+
+
+
+template <class POLY, int dim>
+Tensor<2,dim>
+FE_PolyTensor<POLY,dim>::shape_grad_grad_component (
+ const unsigned int i,
+ const Point<dim> &p,
+ const unsigned int component) const
+{
+ Assert (i<this->dofs_per_cell, ExcIndexRange(i,0,this->dofs_per_cell));
+ Assert (component < dim, ExcIndexRange (component, 0, dim));
+
+ if (cached_point != p || cached_grad_grads.size() == 0)
+ {
+ cached_point = p;
+ cached_grad_grads.resize(poly_space.n());
+ poly_space.compute(p, cached_values, cached_grads, cached_grad_grads);
+ }
+
+ Tensor<2,dim> s;
+ if (inverse_node_matrix.n_cols() == 0)
+ return cached_grad_grads[i][component];
+ else
+ for (unsigned int j=0;j<inverse_node_matrix.n_cols();++j)
+ s += inverse_node_matrix(i,j) * cached_grad_grads[j][component];
+
+ return s;
+}
+
+
+
+//---------------------------------------------------------------------------
+// Data field initialization
+//---------------------------------------------------------------------------
+
+template <class POLY, int dim>
+typename Mapping<dim>::InternalDataBase *
+FE_PolyTensor<POLY,dim>::get_data (const UpdateFlags update_flags,
+ const Mapping<dim> &mapping,
+ const Quadrature<dim> &quadrature) const
+{
+ // generate a new data object and
+ // initialize some fields
+ InternalData* data = new InternalData;
+
+ // check what needs to be
+ // initialized only once and what
+ // on every cell/face/subface we
+ // visit
+ data->update_once = update_once(update_flags);
+ data->update_each = update_each(update_flags);
+ data->update_flags = data->update_once | data->update_each;
+
+ const UpdateFlags flags(data->update_flags);
+ const unsigned int n_q_points = quadrature.n_quadrature_points;
+
+ // some scratch arrays
+ std::vector<Tensor<1,dim> > values(0);
+ std::vector<Tensor<2,dim> > grads(0);
+ std::vector<Tensor<3,dim> > grad_grads(0);
+
+ // initialize fields only if really
+ // necessary. otherwise, don't
+ // allocate memory
+ if (flags & update_values)
+ {
+ values.resize (this->dofs_per_cell);
+ data->shape_values.resize (this->dofs_per_cell);
+ for (unsigned int i=0;i<this->dofs_per_cell;++i)
+ data->shape_values[i].resize (n_q_points);
+ }
+
+ if (flags & update_gradients)
+ {
+ grads.resize (this->dofs_per_cell);
+ data->shape_grads.resize (this->dofs_per_cell);
+ for (unsigned int i=0;i<this->dofs_per_cell;++i)
+ data->shape_grads[i].resize (n_q_points);
+ }
+
+ // if second derivatives through
+ // finite differencing is required,
+ // then initialize some objects for
+ // that
+ if (flags & update_second_derivatives)
+ {
+// grad_grads.resize (this->dofs_per_cell);
+ data->initialize_2nd (this, mapping, quadrature);
+ }
+
+ // Compute shape function values
+ // and derivatives on the reference
+ // cell. Make sure, that for the
+ // node values N_i holds
+ // N_i(v_j)=\delta_ij for all basis
+ // functions v_j
+ if (flags & (update_values | update_gradients))
+ for (unsigned int k=0; k<n_q_points; ++k)
+ {
+ poly_space.compute(quadrature.point(k),
+ values, grads, grad_grads);
+
+ if (flags & update_values)
+ if (inverse_node_matrix.n_cols() == 0)
+ for (unsigned int i=0; i<this->dofs_per_cell; ++i)
+ data->shape_values[i][k] = values[i];
+ else
+ for (unsigned int i=0; i<this->dofs_per_cell; ++i)
+ for (unsigned int j=0; j<this->dofs_per_cell; ++j)
+ data->shape_values[i][k] = inverse_node_matrix(i,j) * values[j];
+
+ if (flags & update_gradients)
+ if (inverse_node_matrix.n_cols() == 0)
+ for (unsigned int i=0; i<this->dofs_per_cell; ++i)
+ data->shape_grads[i][k] = grads[i];
+ else
+ for (unsigned int i=0; i<this->dofs_per_cell; ++i)
+ for (unsigned int j=0; j<this->dofs_per_cell; ++j)
+ data->shape_grads[i][k] = inverse_node_matrix(i,j) * grads[j];
+ }
+ return data;
+}
+
+
+
+
+//---------------------------------------------------------------------------
+// Fill data of FEValues
+//---------------------------------------------------------------------------
+
+template <class POLY, int dim>
+void
+FE_PolyTensor<POLY,dim>::fill_fe_values (const Mapping<dim> &mapping,
+ const typename Triangulation<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 unsigned int n_dofs = this->dofs_per_cell;
+ const unsigned int n_quad = quadrature.n_quadrature_points;
+ const UpdateFlags flags(fe_data.current_update_flags());
+
+ Assert(mapping_type == independent, ExcNotImplemented());
+
+ Assert(!(flags & update_values) || fe_data.shape_values.size() == n_dofs,
+ ExcDimensionMismatch(fe_data.shape_values.size(), n_dofs));
+ Assert(!(flags & update_values) || fe_data.shape_values[0].size() == n_quad,
+ ExcDimensionMismatch(fe_data.shape_values[0].size(), n_quad));
+
+ for (unsigned int i=0; i<n_dofs; ++i)
+ {
+ const unsigned int first = data.shape_function_to_row_table[i];
+
+ if (flags & update_values)
+ for (unsigned int k=0; k<n_quad; ++k)
+ for (unsigned int d=0;d<dim;++d)
+ data.shape_values(first+d,k) = fe_data.shape_values[i][k][d];
+
+ if (flags & update_gradients)
+ {
+ std::vector<Tensor<2,dim> > shape_grads1 (n_quad);
+ mapping.transform_covariant(fe_data.shape_grads[i], 0,
+ shape_grads1,
+ mapping_data);
+ for (unsigned int k=0; k<n_quad; ++k)
+ for (unsigned int d=0;d<dim;++d)
+ data.shape_gradients[first+d][k] = shape_grads1[k][d];
+ }
+ }
+
+ const typename QProjector<dim>::DataSetDescriptor dsd;
+ if (flags & update_second_derivatives)
+ this->compute_2nd (mapping, cell, dsd.cell(),
+ mapping_data, fe_data, data);
+}
+
+
+
+template <class POLY, int dim>
+void
+FE_PolyTensor<POLY,dim>::fill_fe_face_values (const Mapping<dim> &mapping,
+ const typename Triangulation<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 unsigned int n_dofs = this->dofs_per_cell;
+ const unsigned int n_quad = quadrature.n_quadrature_points;
+ // offset determines which data set
+ // to take (all data sets for all
+ // faces are stored contiguously)
+
+ const typename QProjector<dim>::DataSetDescriptor dsd;
+ const typename QProjector<dim>::DataSetDescriptor offset
+ = dsd.face (face, cell->face_orientation(face),
+ n_quad);
+
+ const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
+
+ Assert(mapping_type == independent, ExcNotImplemented());
+//TODO: Size assertions
+
+ for (unsigned int i=0; i<n_dofs; ++i)
+ {
+ const unsigned int first = data.shape_function_to_row_table[i];
+
+ if (flags & update_values)
+ for (unsigned int k=0; k<n_quad; ++k)
+ for (unsigned int d=0;d<dim;++d)
+ data.shape_values(first+d,k) = fe_data.shape_values[i][k+offset][d];
+
+ if (flags & update_gradients)
+ {
+ std::vector<Tensor<2,dim> > shape_grads1 (n_quad);
+ mapping.transform_covariant(fe_data.shape_grads[i], offset,
+ shape_grads1, mapping_data);
+ for (unsigned int k=0; k<n_quad; ++k)
+ for (unsigned int d=0;d<dim;++d)
+ data.shape_gradients[first+d][k] = shape_grads1[k][d];
+ }
+ }
+
+ if (flags & update_second_derivatives)
+ this->compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
+}
+
+
+
+template <class POLY, int dim>
+void
+FE_PolyTensor<POLY,dim>::fill_fe_subface_values (const Mapping<dim> &mapping,
+ const typename Triangulation<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 unsigned int n_dofs = this->dofs_per_cell;
+ const unsigned int n_quad = quadrature.n_quadrature_points;
+ // offset determines which data set
+ // to take (all data sets for all
+ // sub-faces are stored contiguously)
+
+ const typename QProjector<dim>::DataSetDescriptor dsd;
+ const typename QProjector<dim>::DataSetDescriptor offset
+ = dsd.sub_face (face, subface, cell->face_orientation(face),
+ n_quad);
+
+ const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
+
+ Assert(mapping_type == independent, ExcNotImplemented());
+//TODO: Size assertions
+
+ for (unsigned int i=0; i<n_dofs; ++i)
+ {
+ const unsigned int first = data.shape_function_to_row_table[i];
+
+ if (flags & update_values)
+ for (unsigned int k=0; k<n_quad; ++k)
+ for (unsigned int d=0;d<dim;++d)
+ data.shape_values(first+d,k) = fe_data.shape_values[i][k+offset][d];
+
+ if (flags & update_gradients)
+ {
+ std::vector<Tensor<2,dim> > shape_grads1 (n_quad);
+ mapping.transform_covariant(fe_data.shape_grads[i], offset,
+ shape_grads1, mapping_data);
+ for (unsigned int k=0; k<n_quad; ++k)
+ for (unsigned int d=0;d<dim;++d)
+ data.shape_gradients[first+d][k] = shape_grads1[k][d];
+ }
+ }
+
+ if (flags & update_second_derivatives)
+ this->compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
+}
+
+
+
+template <class POLY, int dim>
+unsigned int
+FE_PolyTensor<POLY,dim>::n_base_elements () const
+{
+ return 1;
+}
+
+
+
+template <class POLY, int dim>
+const FiniteElement<dim> &
+FE_PolyTensor<POLY,dim>::base_element (const unsigned int index) const
+{
+ Assert (index==0, ExcIndexRange(index, 0, 1));
+ return *this;
+}
+
+
+
+template <class POLY, int dim>
+unsigned int
+FE_PolyTensor<POLY,dim>::element_multiplicity (const unsigned int index) const
+{
+ Assert (index==0, ExcIndexRange(index, 0, 1));
+ return 1;
+}
+
+
+
+template class FE_PolyTensor<PolynomialsRaviartThomas<deal_II_dimension>,deal_II_dimension>;