]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implemented the functions hp_vertex_dof_identities,
authorsimfeld <sml.imfeld@gmail.com>
Thu, 6 Jul 2017 12:24:24 +0000 (14:24 +0200)
committersimfeld <sml.imfeld@gmail.com>
Thu, 6 Jul 2017 12:24:24 +0000 (14:24 +0200)
hp_line_dof_identities and hp_quad_dof_identities
for FE_FaceQ and FE_FaceP. These functions are required
when using FE_Face elements in an hp finite element
context. Since these elements are discontinuous, there
is nothing to be returned.

include/deal.II/fe/fe_face.h
source/fe/fe_face.cc

index 82943decfee39bd68fbecdeb9604e83448727873..3946e4a322acd4c47ce5539eb2b4dfe6f3bbbe44 100644 (file)
@@ -102,6 +102,55 @@ public:
   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,
+   * both numbered between zero and the corresponding value of dofs_per_vertex
+   * of the two finite elements. 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.
+   *
+   * This being a discontinuous element, the set of such constraints is of
+   * course empty.
+   */
+  virtual
+  std::vector<std::pair<unsigned int, unsigned int> >
+  hp_vertex_dof_identities (const FiniteElement<dim, spacedim> &fe_other) const;
+
+  /**
+   * Same as hp_vertex_dof_indices(), except that the function treats degrees
+   * of freedom on lines.
+   *
+   * This being a discontinuous element, the set of such constraints is of
+   * course empty.
+   */
+  virtual
+  std::vector<std::pair<unsigned int, unsigned int> >
+  hp_line_dof_identities (const FiniteElement<dim, spacedim> &fe_other) const;
+
+  /**
+   * Same as hp_vertex_dof_indices(), except that the function treats degrees
+   * of freedom on quads.
+   *
+   * This being a discontinuous element, the set of such constraints is of
+   * course empty.
+   */
+  virtual
+  std::vector<std::pair<unsigned int, unsigned int> >
+  hp_quad_dof_identities (const FiniteElement<dim, spacedim> &fe_other) const;
+
   /**
    * Return whether this element implements its hanging node constraints in
    * the new way, which has to be used to make elements "hp compatible".
@@ -408,6 +457,50 @@ public:
    */
   virtual bool hp_constraints_are_implemented () const;
 
+  /**
+   * 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,
+   * both numbered between zero and the corresponding value of dofs_per_vertex
+   * of the two finite elements. 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.
+   *
+   * This being a discontinuous element, the set of such constraints is of
+   * course empty.
+   */
+  virtual
+  std::vector<std::pair<unsigned int, unsigned int> >
+  hp_vertex_dof_identities (const FiniteElement<dim, spacedim> &fe_other) const;
+
+  /**
+   * Same as hp_vertex_dof_indices(), except that the function treats degrees
+   * of freedom on lines.
+   *
+   * This being a discontinuous element, the set of such constraints is of
+   * course empty.
+   */
+  virtual
+  std::vector<std::pair<unsigned int, unsigned int> >
+  hp_line_dof_identities (const FiniteElement<dim, spacedim> &fe_other) const;
+
+  /**
+   * Same as hp_vertex_dof_indices(), except that the function treats degrees
+   * of freedom on quads.
+   *
+   * This being a discontinuous element, the set of such constraints is of
+   * course empty.
+   */
+  virtual
+  std::vector<std::pair<unsigned int, unsigned int> >
+  hp_quad_dof_identities (const FiniteElement<dim, spacedim> &fe_other) const;
+
   /**
    * Return whether this element dominates the one given as argument when they
    * meet at a common face, whether it is the other way around, whether
index dad2e319f28ab5d5cf0578094e83a12c68c2007a..eb1d5efbe666e77df2b4c16ec6896b88d82fa7a2 100644 (file)
@@ -263,6 +263,51 @@ FE_FaceQ<dim,spacedim>::hp_constraints_are_implemented () const
 
 
 
+template <int dim, int spacedim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_FaceQ<dim, spacedim>::
+hp_vertex_dof_identities (const FiniteElement<dim, spacedim> &/*fe_other*/) const
+{
+  // this element is discontinuous, so by definition there can
+  // be no identities between its dofs and those of any neighbor
+  // (of whichever type the neighbor may be -- after all, we have
+  // no face dofs on this side to begin with)
+  return
+    std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
+template <int dim, int spacedim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_FaceQ<dim, spacedim>::
+hp_line_dof_identities (const FiniteElement<dim, spacedim> &/*fe_other*/) const
+{
+  // this element is discontinuous, so by definition there can
+  // be no identities between its dofs and those of any neighbor
+  // (of whichever type the neighbor may be -- after all, we have
+  // no face dofs on this side to begin with)
+  return
+    std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
+template <int dim, int spacedim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_FaceQ<dim, spacedim>::
+hp_quad_dof_identities (const FiniteElement<dim, spacedim> &/*fe_other*/) const
+{
+  // this element is discontinuous, so by definition there can
+  // be no identities between its dofs and those of any neighbor
+  // (of whichever type the neighbor may be -- after all, we have
+  // no face dofs on this side to begin with)
+  return
+    std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
 template <int dim, int spacedim>
 FiniteElementDomination::Domination
 FE_FaceQ<dim,spacedim>::
@@ -591,6 +636,51 @@ FE_FaceP<dim,spacedim>::hp_constraints_are_implemented () const
 
 
 
+template <int dim, int spacedim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_FaceP<dim, spacedim>::
+hp_vertex_dof_identities (const FiniteElement<dim, spacedim> &/*fe_other*/) const
+{
+  // this element is discontinuous, so by definition there can
+  // be no identities between its dofs and those of any neighbor
+  // (of whichever type the neighbor may be -- after all, we have
+  // no face dofs on this side to begin with)
+  return
+    std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
+template <int dim, int spacedim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_FaceP<dim, spacedim>::
+hp_line_dof_identities (const FiniteElement<dim, spacedim> &/*fe_other*/) const
+{
+  // this element is discontinuous, so by definition there can
+  // be no identities between its dofs and those of any neighbor
+  // (of whichever type the neighbor may be -- after all, we have
+  // no face dofs on this side to begin with)
+  return
+    std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
+template <int dim, int spacedim>
+std::vector<std::pair<unsigned int, unsigned int> >
+FE_FaceP<dim, spacedim>::
+hp_quad_dof_identities (const FiniteElement<dim, spacedim> &/*fe_other*/) const
+{
+  // this element is discontinuous, so by definition there can
+  // be no identities between its dofs and those of any neighbor
+  // (of whichever type the neighbor may be -- after all, we have
+  // no face dofs on this side to begin with)
+  return
+    std::vector<std::pair<unsigned int, unsigned int> > ();
+}
+
+
+
 template <int dim, int spacedim>
 FiniteElementDomination::Domination
 FE_FaceP<dim,spacedim>::

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.