From: Martin Kronbichler Date: Mon, 23 Feb 2015 14:29:36 +0000 (+0100) Subject: Comment on the behavior of FE_TraceQ as compared to FE_Q. X-Git-Tag: v8.3.0-rc1~430^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F590%2Fhead;p=dealii.git Comment on the behavior of FE_TraceQ as compared to FE_Q. --- diff --git a/include/deal.II/fe/fe_trace.h b/include/deal.II/fe/fe_trace.h index 48c3fd8ffe..757d30a6e8 100644 --- a/include/deal.II/fe/fe_trace.h +++ b/include/deal.II/fe/fe_trace.h @@ -40,8 +40,6 @@ DEAL_II_NAMESPACE_OPEN * * @todo Polynomials::LagrangeEquidistant should be and will be replaced by * Polynomials::LagrangeGaussLobatto as soon as such a polynomial set exists. - * @todo so far, hanging nodes are not implemented - * */ template @@ -55,6 +53,11 @@ public: */ FE_TraceQ(unsigned int p); + /** + * @p clone function instead of a copy constructor. + * + * This function is needed by the constructors of @p FESystem. + */ virtual FiniteElement *clone() const; /** @@ -78,9 +81,6 @@ public: virtual std::pair, std::vector > get_constant_modes () const; - - - /** * Return whether this element implements its hanging node constraints in * the new way, which has to be used to make elements "hp compatible". @@ -130,6 +130,7 @@ private: * Store a copy of FE_Q for delegating the hp-constraints functionality. */ FE_Q fe_q; + /** * Return vector with dofs per vertex, line, quad, hex. */ diff --git a/source/fe/fe_trace.cc b/source/fe/fe_trace.cc index 4d5306fbc2..0acd179ef7 100644 --- a/source/fe/fe_trace.cc +++ b/source/fe/fe_trace.cc @@ -68,12 +68,9 @@ template std::string FE_TraceQ::get_name () const { - // note that the - // FETools::get_fe_from_name - // function depends on the - // particular format of the string - // this function returns, so they - // have to be kept in synch + // note that the FETools::get_fe_from_name function depends on the + // particular format of the string this function returns, so they have to be + // kept in synch std::ostringstream namebuf; namebuf << "FE_TraceQ<" @@ -95,6 +92,11 @@ FE_TraceQ::has_support_on_face (const unsigned int shape_index, Assert (face_index < GeometryInfo::faces_per_cell, ExcIndexRange (face_index, 0, GeometryInfo::faces_per_cell)); + // FE_TraceQ shares the numbering of elemental degrees of freedom with FE_Q + // except for the missing interior ones (quad dofs in 2D and hex dofs in + // 3D). Therefore, it is safe to ask fe_q for the corresponding + // information. The assertion 'shape_index < this->dofs_per_cell' will make + // sure that we only access the trace dofs. return fe_q.has_support_on_face (shape_index, face_index); } @@ -117,6 +119,9 @@ template std::vector FE_TraceQ::get_dpo_vector (const unsigned int deg) { + // This constructs FE_TraceQ in exactly the same way as FE_Q except for the + // interior degrees of freedom that are not present here (lnine in 1D, quad + // in 2D, hex in 3D). AssertThrow(deg>0,ExcMessage("FE_TraceQ needs to be of degree > 0.")); std::vector dpo(dim+1, 1U); dpo[dim]=0;