]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
the new trace element
authorangela.klewinghaus <angela.klewinghaus@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 25 Feb 2013 16:08:44 +0000 (16:08 +0000)
committerangela.klewinghaus <angela.klewinghaus@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 25 Feb 2013 16:08:44 +0000 (16:08 +0000)
git-svn-id: https://svn.dealii.org/branches/branch_face_elements@28555 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/fe/fe_trace.h [new file with mode: 0644]

diff --git a/deal.II/include/deal.II/fe/fe_trace.h b/deal.II/include/deal.II/fe/fe_trace.h
new file mode 100644 (file)
index 0000000..e98b6c8
--- /dev/null
@@ -0,0 +1,88 @@
+
+#ifndef __deal2__fe_trace_h
+#define __deal2__fe_trace_h
+
+#include <deal.II/base/config.h>
+#include <deal.II/base/tensor_product_polynomials.h>
+#include <deal.II/fe/fe_poly_face.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+/**
+ * A finite element, which is the trace of Fe_Q elements, that is
+ * a tensor product of polynomials on the faces,
+ * undefined in the interior of the cells and continuous. The basis functions on
+ * the faces are from Polynomials::LagrangeEquidistant
+ *
+ * This finite element is the trace space of FE_Q on the
+ * faces.
+ *
+ * @note Since these are only finite elements on faces, only
+ * FEFaceValues and FESubfaceValues will be able to extract reasonable
+ * values from any face polynomial. In order to make the use of
+ * FESystem simpler, FEValues objects will not fail using this finite
+ * element space, but all shape function values extracted will equal
+ * to zero.
+ *
+ * @todo Polynomials::LagrangeEquidistant should be and will be
+ * replaced by Polynomials::LagrangeGaussLobatto as soon as such a
+ * polynomial set exists.
+ * @todo so far only implemented for 2D
+ *
+ */
+
+template <int dim, int spacedim=dim>
+class FE_TraceQ : public FE_PolyFace<TensorProductPolynomials<dim-1>, dim, spacedim>
+{
+public:
+  /**
+   * Constructor for tensor product
+   * polynomials of degree
+   * <tt>p</tt>. The shape
+   * functions created using this
+   * constructor correspond to
+   * Legendre polynomials in each
+   * coordinate direction.
+   */
+  FE_TraceQ(unsigned int p);
+
+  virtual FiniteElement<dim,spacedim> *clone() const;
+
+  /**
+   * Return a string that uniquely
+   * identifies a finite
+   * element. This class returns
+   * <tt>FE_DGQ<dim>(degree)</tt>, with
+   * <tt>dim</tt> and <tt>degree</tt>
+   * 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<unsigned int> get_dpo_vector (const unsigned int deg);
+};
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif

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.