From 10b70fa0a1364b665ff4fd04c3825458bec087a4 Mon Sep 17 00:00:00 2001 From: bangerth Date: Thu, 27 Jul 2006 00:39:44 +0000 Subject: [PATCH] Add functions to compute vertex dof equivalences git-svn-id: https://svn.dealii.org/trunk@13439 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe.h | 67 +++++++++++++++++++++++++++++++ deal.II/deal.II/include/fe/fe_q.h | 67 +++++++++++++++++++++++++++++++ deal.II/deal.II/source/fe/fe.cc | 34 +++++++++++++++- deal.II/deal.II/source/fe/fe_q.cc | 46 +++++++++++++++++++++ 4 files changed, 213 insertions(+), 1 deletion(-) diff --git a/deal.II/deal.II/include/fe/fe.h b/deal.II/deal.II/include/fe/fe.h index e5a792a703..7dd6811cdf 100644 --- a/deal.II/deal.II/include/fe/fe.h +++ b/deal.II/deal.II/include/fe/fe.h @@ -898,6 +898,73 @@ class FiniteElement : public Subscriptor, get_interpolation_matrix (const FiniteElement &source, FullMatrix &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 > + 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. + */ + 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. + */ + virtual + std::vector > + hp_quad_dof_identities (const FiniteElement &fe_other) const; + + //@} /** * Comparison operator. We also diff --git a/deal.II/deal.II/include/fe/fe_q.h b/deal.II/deal.II/include/fe/fe_q.h index fa44e97be7..9ebe9dab1d 100644 --- a/deal.II/deal.II/include/fe/fe_q.h +++ b/deal.II/deal.II/include/fe/fe_q.h @@ -289,6 +289,73 @@ class FE_Q : public FE_Poly,dim> 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 > + 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. + */ + 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. + */ + virtual + std::vector > + hp_quad_dof_identities (const FiniteElement &fe_other) const; + + //@} + /** * Determine an estimate for the * memory consumption (in bytes) diff --git a/deal.II/deal.II/source/fe/fe.cc b/deal.II/deal.II/source/fe/fe.cc index 59081535ef..13b1c8ea87 100644 --- a/deal.II/deal.II/source/fe/fe.cc +++ b/deal.II/deal.II/source/fe/fe.cc @@ -440,7 +440,39 @@ get_interpolation_matrix (const FiniteElement &, ExcInterpolationNotImplemented()); } - + + +template +std::vector > +FiniteElement:: +hp_vertex_dof_identities (const FiniteElement &) const +{ + Assert (false, ExcNotImplemented()); + return std::vector > (); +} + + + +template +std::vector > +FiniteElement:: +hp_line_dof_identities (const FiniteElement &) const +{ + Assert (false, ExcNotImplemented()); + return std::vector > (); +} + + + +template +std::vector > +FiniteElement:: +hp_quad_dof_identities (const FiniteElement &) const +{ + Assert (false, ExcNotImplemented()); + return std::vector > (); +} + template diff --git a/deal.II/deal.II/source/fe/fe_q.cc b/deal.II/deal.II/source/fe/fe_q.cc index d21246db90..3b2fbf8f0a 100644 --- a/deal.II/deal.II/source/fe/fe_q.cc +++ b/deal.II/deal.II/source/fe/fe_q.cc @@ -324,6 +324,52 @@ get_interpolation_matrix (const FiniteElement &x_source_fe, } + +template +std::vector > +FE_Q:: +hp_vertex_dof_identities (const FiniteElement &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 *fe_q_other = dynamic_cast*>(&fe_other); + if (fe_q_other != 0) + return + std::vector > (1, std::make_pair (0U, 0U)); + else + { + Assert (false, ExcNotImplemented()); + return std::vector > (); + } +} + + + +template +std::vector > +FE_Q:: +hp_line_dof_identities (const FiniteElement &/*fe_other*/) const +{ + Assert (false, ExcNotImplemented()); + return std::vector > (); +} + + + +template +std::vector > +FE_Q:: +hp_quad_dof_identities (const FiniteElement &/*fe_other*/) const +{ + Assert (false, ExcNotImplemented()); + return std::vector > (); +} + + //--------------------------------------------------------------------------- // Auxiliary functions //--------------------------------------------------------------------------- -- 2.39.5