From 769931b1f0c01795aea0a3bdf9d7683a6d36b424 Mon Sep 17 00:00:00 2001 From: "angela.klewinghaus" Date: Mon, 25 Feb 2013 16:08:49 +0000 Subject: [PATCH] the new trace element git-svn-id: https://svn.dealii.org/branches/branch_face_elements@28556 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/source/fe/fe_trace.cc | 247 +++++++++++++++++++++++++++++ deal.II/source/fe/fe_trace.inst.in | 21 +++ 2 files changed, 268 insertions(+) create mode 100644 deal.II/source/fe/fe_trace.cc create mode 100644 deal.II/source/fe/fe_trace.inst.in diff --git a/deal.II/source/fe/fe_trace.cc b/deal.II/source/fe/fe_trace.cc new file mode 100644 index 0000000000..390e0e0319 --- /dev/null +++ b/deal.II/source/fe/fe_trace.cc @@ -0,0 +1,247 @@ +#include +#include +#include +#include +#include +#include + +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +//using namespace dealii; + +#include +#include +#include + +DEAL_II_NAMESPACE_OPEN + +/** FE-Class for continuous face functions +* adaption of FE_Q and FE_FaceQ +*/ + +template +class FE_TraceQ : public FE_PolyFace, dim, spacedim> +{ + public: + /** + * Constructor for tensor product + * polynomials of degree + * p. The shape + * functions created using this + * constructor correspond to + * Legendre polynomials in each + * coordinate direction. + */ + FE_TraceQ(unsigned int p); + + virtual FiniteElement* clone() const; + + /** + * Return a string that uniquely + * identifies a finite + * element. This class returns + * FE_DGQ(degree) , with + * dim and degree + * replaced by appropriate + * values. + */ + virtual std::string get_name () const; + + /** + * Check for non-zero values on a face. + * + * This function returns + * @p true, if the shape + * function @p shape_index has + * non-zero values on the face + * @p face_index. + * + * Implementation of the + * interface in + * FiniteElement + */ + virtual bool has_support_on_face (const unsigned int shape_index, + const unsigned int face_index) const; + + private: + /** + * Return vector with dofs per + * vertex, line, quad, hex. + */ + static std::vector get_dpo_vector (const unsigned int deg); +}; + + +template +FE_TraceQ::FE_TraceQ (const unsigned int degree) + : + FE_PolyFace, dim, spacedim> ( + TensorProductPolynomials(Polynomials::LagrangeEquidistant::generate_complete_basis(degree)), + FiniteElementData(get_dpo_vector(degree), 1, degree, FiniteElementData::L2), + std::vector(1,true)) +{ + std::vector renumber (this->dofs_per_face); + FETools::hierarchic_to_lexicographic_numbering (degree, renumber); + this->poly_space.set_numbering(renumber); +} + + +template +FiniteElement* +FE_TraceQ::clone() const +{ + return new FE_TraceQ(this->degree); +} + + +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 + + std::ostringstream namebuf; + namebuf << "FE_TraceQ<" << dim << ">(" << this->degree << ")"; + + return namebuf.str(); +} + +template +bool +FE_TraceQ::has_support_on_face (const unsigned int shape_index, + const unsigned int face_index) const +{ + Assert (shape_index < this->dofs_per_cell, + ExcIndexRange (shape_index, 0, this->dofs_per_cell)); + Assert (face_index < GeometryInfo::faces_per_cell, + ExcIndexRange (face_index, 0, GeometryInfo::faces_per_cell)); + + // let's see whether this is a + // vertex + if (shape_index < this->first_line_index) + { + // for Q elements, there is + // one dof per vertex, so + // shape_index==vertex_number. check + // whether this vertex is + // on the given face. thus, + // for each face, give a + // list of vertices + const unsigned int vertex_no = shape_index; + Assert (vertex_no < GeometryInfo::vertices_per_cell, + ExcInternalError()); + + for (unsigned int v=0; v::vertices_per_face; ++v) + if (GeometryInfo::face_to_cell_vertices(face_index, v) == vertex_no) + return true; + + return false; + } + else if (shape_index < this->first_quad_index) + // ok, dof is on a line + { + const unsigned int line_index + = (shape_index - this->first_line_index) / this->dofs_per_line; + Assert (line_index < GeometryInfo::lines_per_cell, + ExcInternalError()); + + // in 2d, the line is the + // face, so get the line + // index + if (dim == 2) + return (line_index == face_index); + else if (dim == 3) + { + // see whether the + // given line is on the + // given face. + for (unsigned int l=0; l::lines_per_face; ++l) + if (GeometryInfo<3>::face_to_cell_lines(face_index, l) == line_index) + return true; + + return false; + } + else + Assert (false, ExcNotImplemented()); + } + else if (shape_index < this->first_hex_index) + // dof is on a quad + { + const unsigned int quad_index + = (shape_index - this->first_quad_index) / this->dofs_per_quad; + Assert (static_cast(quad_index) < + static_cast(GeometryInfo::quads_per_cell), + ExcInternalError()); + + // in 2d, cell bubble are + // zero on all faces. but + // we have treated this + // case above already + Assert (dim != 2, ExcInternalError()); + + // in 3d, + // quad_index=face_index + if (dim == 3) + return (quad_index == face_index); + else + Assert (false, ExcNotImplemented()); + } + else + // dof on hex + { + // can only happen in 3d, + // but this case has + // already been covered + // above + Assert (false, ExcNotImplemented()); + return false; + } + + // we should not have gotten here + Assert (false, ExcInternalError()); + return false; +} + +template +std::vector +FE_TraceQ::get_dpo_vector (const unsigned int deg) +{ + std::vector dpo(dim+1, 1U); + dpo[dim]=0; + dpo[0]=1; + for (unsigned int i=1;i 1 + template class FE_TraceQ; +#endif + } + -- 2.39.5