get_interpolation_matrix (const FiniteElement<dim> &source,
FullMatrix<double> &matrix) const;
//@}
+
+ /**
+ * @name Functions to support hp
+ * @{
+ */
+
+ /**
+ * If, on a vertex, several
+ * finite elements are active,
+ * the hp code first assigns the
+ * degrees of freedom of each of
+ * these FEs different global
+ * indices. It then calls this
+ * function to find out which of
+ * them should get identical
+ * values, and consequently can
+ * receive the same global DoF
+ * index. This function therefore
+ * returns a list of identities
+ * between DoFs of the present
+ * finite element object with the
+ * DoFs of @p fe_other, which is
+ * a reference to a finite
+ * element object representing
+ * one of the other finite
+ * elements active on this
+ * particular vertex. The
+ * function computes which of the
+ * degrees of freedom of the two
+ * finite element objects are
+ * equivalent, and returns a list
+ * of pairs of global dof indices
+ * in @p identities. The first
+ * index of each pair denotes one
+ * of the vertex dofs of the
+ * present element, whereas the
+ * second is the corresponding
+ * index of the other finite
+ * element.
+ */
+ virtual
+ std::vector<std::pair<unsigned int, unsigned int> >
+ hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+ /**
+ * Same as
+ * hp_vertex_dof_indices(),
+ * except that the function
+ * treats degrees of freedom on
+ * lines.
+ */
+ virtual
+ std::vector<std::pair<unsigned int, unsigned int> >
+ hp_line_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+ /**
+ * Same as
+ * hp_vertex_dof_indices(),
+ * except that the function
+ * treats degrees of freedom on
+ * quads.
+ */
+ virtual
+ std::vector<std::pair<unsigned int, unsigned int> >
+ hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+ //@}
/**
* Comparison operator. We also
virtual bool has_support_on_face (const unsigned int shape_index,
const unsigned int face_index) const;
+ /**
+ * @name Functions to support hp
+ * @{
+ */
+
+ /**
+ * If, on a vertex, several
+ * finite elements are active,
+ * the hp code first assigns the
+ * degrees of freedom of each of
+ * these FEs different global
+ * indices. It then calls this
+ * function to find out which of
+ * them should get identical
+ * values, and consequently can
+ * receive the same global DoF
+ * index. This function therefore
+ * returns a list of identities
+ * between DoFs of the present
+ * finite element object with the
+ * DoFs of @p fe_other, which is
+ * a reference to a finite
+ * element object representing
+ * one of the other finite
+ * elements active on this
+ * particular vertex. The
+ * function computes which of the
+ * degrees of freedom of the two
+ * finite element objects are
+ * equivalent, and returns a list
+ * of pairs of global dof indices
+ * in @p identities. The first
+ * index of each pair denotes one
+ * of the vertex dofs of the
+ * present element, whereas the
+ * second is the corresponding
+ * index of the other finite
+ * element.
+ */
+ virtual
+ std::vector<std::pair<unsigned int, unsigned int> >
+ hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+ /**
+ * Same as
+ * hp_vertex_dof_indices(),
+ * except that the function
+ * treats degrees of freedom on
+ * lines.
+ */
+ virtual
+ std::vector<std::pair<unsigned int, unsigned int> >
+ hp_line_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+ /**
+ * Same as
+ * hp_vertex_dof_indices(),
+ * except that the function
+ * treats degrees of freedom on
+ * quads.
+ */
+ virtual
+ std::vector<std::pair<unsigned int, unsigned int> >
+ hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
+
+ //@}
+
/**
* Determine an estimate for the
* memory consumption (in bytes)
ExcInterpolationNotImplemented());
}
-
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FiniteElement<dim>::
+hp_vertex_dof_identities (const FiniteElement<dim> &) const
+{
+ Assert (false, ExcNotImplemented());
+ return std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FiniteElement<dim>::
+hp_line_dof_identities (const FiniteElement<dim> &) const
+{
+ Assert (false, ExcNotImplemented());
+ return std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FiniteElement<dim>::
+hp_quad_dof_identities (const FiniteElement<dim> &) const
+{
+ Assert (false, ExcNotImplemented());
+ return std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
template <int dim>
}
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_Q<dim>::
+hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const
+{
+ // we can presently only compute
+ // these identities if both FEs are
+ // FE_Qs. in that case, there
+ // should be exactly one single DoF
+ // of each FE at a vertex, and they
+ // should have identical value
+ const FE_Q<dim> *fe_q_other = dynamic_cast<const FE_Q<dim>*>(&fe_other);
+ if (fe_q_other != 0)
+ return
+ std::vector<std::pair<unsigned int, unsigned int> > (1, std::make_pair (0U, 0U));
+ else
+ {
+ Assert (false, ExcNotImplemented());
+ return std::vector<std::pair<unsigned int, unsigned int> > ();
+ }
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_Q<dim>::
+hp_line_dof_identities (const FiniteElement<dim> &/*fe_other*/) const
+{
+ Assert (false, ExcNotImplemented());
+ return std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
+template <int dim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_Q<dim>::
+hp_quad_dof_identities (const FiniteElement<dim> &/*fe_other*/) const
+{
+ Assert (false, ExcNotImplemented());
+ return std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
//---------------------------------------------------------------------------
// Auxiliary functions
//---------------------------------------------------------------------------