From: simfeld Date: Thu, 6 Jul 2017 12:24:24 +0000 (+0200) Subject: Implemented the functions hp_vertex_dof_identities, X-Git-Tag: v9.0.0-rc1~716^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=2de8297ea5b2dbf8afd42170c23ac129801e61cf;p=dealii.git Implemented the functions hp_vertex_dof_identities, 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. --- diff --git a/include/deal.II/fe/fe_face.h b/include/deal.II/fe/fe_face.h index 82943decfe..3946e4a322 100644 --- a/include/deal.II/fe/fe_face.h +++ b/include/deal.II/fe/fe_face.h @@ -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 > + hp_vertex_dof_identities (const FiniteElement &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 > + hp_line_dof_identities (const FiniteElement &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 > + hp_quad_dof_identities (const FiniteElement &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 > + hp_vertex_dof_identities (const FiniteElement &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 > + hp_line_dof_identities (const FiniteElement &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 > + hp_quad_dof_identities (const FiniteElement &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 diff --git a/source/fe/fe_face.cc b/source/fe/fe_face.cc index dad2e319f2..eb1d5efbe6 100644 --- a/source/fe/fe_face.cc +++ b/source/fe/fe_face.cc @@ -263,6 +263,51 @@ FE_FaceQ::hp_constraints_are_implemented () const +template +std::vector > +FE_FaceQ:: +hp_vertex_dof_identities (const FiniteElement &/*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 > (); +} + + + +template +std::vector > +FE_FaceQ:: +hp_line_dof_identities (const FiniteElement &/*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 > (); +} + + + +template +std::vector > +FE_FaceQ:: +hp_quad_dof_identities (const FiniteElement &/*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 > (); +} + + + template FiniteElementDomination::Domination FE_FaceQ:: @@ -591,6 +636,51 @@ FE_FaceP::hp_constraints_are_implemented () const +template +std::vector > +FE_FaceP:: +hp_vertex_dof_identities (const FiniteElement &/*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 > (); +} + + + +template +std::vector > +FE_FaceP:: +hp_line_dof_identities (const FiniteElement &/*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 > (); +} + + + +template +std::vector > +FE_FaceP:: +hp_quad_dof_identities (const FiniteElement &/*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 > (); +} + + + template FiniteElementDomination::Domination FE_FaceP::